cbadefpinch.F File Reference

#include "implicit_f.inc"
#include "mvsiz_p.inc"
#include "vectorize.inc"

Go to the source code of this file.

Functions/Subroutines
subroutine	cbadefpinch (jft, jlt, ng, vqg, vdef, veta, vksi, tc, nplat, iplat, bcp, bp, vpinchxyz, vdefpinch, tnpg, dbetadxy, vpincht1, vpincht2, bpinchdamp)

Function/Subroutine Documentation

cbadefpinch()

subroutine cbadefpinch	(	integer	jft,
		integer	jlt,
		integer	ng,
			vqg,
			vdef,
			veta,
			vksi,
			tc,
		integer	nplat,
		integer, dimension(*)	iplat,
			bcp,
			bp,
			vpinchxyz,
			vdefpinch,
			tnpg,
			dbetadxy,
			vpincht1,
			vpincht2,
			bpinchdamp )

Definition at line 28 of file cbadefpinch.F.

C-----------------------------------------------
C   I M P L I C I T   T Y P E S
C-----------------------------------------------
#include      "implicit_f.inc"
C-----------------------------------------------
C   G l o b a l   P a r a m e t e r s
C-----------------------------------------------
#include      "mvsiz_p.inc"
C-----------------------------------------------
C   D U M M Y   A R G U M E N T S
C-----------------------------------------------
      INTEGER NG,JFT,JLT,NPLAT,IPLAT(*)
      my_real 
     .   vksi(4,4),veta(4,4),
     .   bcp(mvsiz,8),bp(mvsiz,4),
     .   vpinchxyz(mvsiz,4),vdef(mvsiz,8),vdefpinch(mvsiz,3),
     .   vqg(mvsiz,3,3,4),tnpg(mvsiz,4,4),tc(mvsiz,2,2),
     .   dbetadxy(mvsiz,3),vpincht1(mvsiz,4),vpincht2(mvsiz,4),
     .   bpinchdamp(mvsiz,8)
C-----------------------------------------------
C   L O C A L   V A R I A B L E S
C-----------------------------------------------
      INTEGER I,EP
      my_real bcx,bcy,d1,d2,
     .        dn1dx,dn2dx,dn3dx,dn4dx,
     .        dn1dy,dn2dy,dn3dy,dn4dy
C---------------------------------
C
C----------------------------------
C     MATRIX [B] DUE TO PINCHING
C----------------------------------
#include "vectorize.inc"
       DO i=jft,jlt 
         ep=iplat(i)
C----------------------------------
C     MATRIX [BCP] FOR TRNASVERSE SHEAR DUE TO PINCHING
C----------------------------------
 
C-------Node 1 contribution
 
         d1=vksi(1,ng)*tc(ep,1,1)+veta(1,ng)*tc(ep,2,1)
         d2=vksi(1,ng)*tc(ep,1,2)+veta(1,ng)*tc(ep,2,2)
C
C        dN1/dx, dN1/dy
C
         bcp(ep,1) = vqg(ep,1,1,ng)*d1
         bcp(ep,2) = vqg(ep,1,2,ng)*d2
 
C-------Node 2 contribution
 
         d1=vksi(2,ng)*tc(ep,1,1)+veta(2,ng)*tc(ep,2,1)
         d2=vksi(2,ng)*tc(ep,1,2)+veta(2,ng)*tc(ep,2,2)
 
C        dN2/dx, dN2/dy
         bcp(ep,3) = vqg(ep,1,1,ng)*d1
         bcp(ep,4) = vqg(ep,1,2,ng)*d2
 
C-------Node 3 contribution
 
         d1=vksi(3,ng)*tc(ep,1,1)+veta(3,ng)*tc(ep,2,1)
         d2=vksi(3,ng)*tc(ep,1,2)+veta(3,ng)*tc(ep,2,2)
 
C        dN3/dx, dN3/dy
         bcp(ep,5) = vqg(ep,1,1,ng)*d1
         bcp(ep,6) = vqg(ep,1,2,ng)*d2
 
C-------Node 4 contribution
 
         d1=vksi(4,ng)*tc(ep,1,1)+veta(4,ng)*tc(ep,2,1)
         d2=vksi(4,ng)*tc(ep,1,2)+veta(4,ng)*tc(ep,2,2)
 
C        dN4/dx, dN4/dy
         bcp(ep,7) = vqg(ep,1,1,ng)*d1
         bcp(ep,8) = vqg(ep,1,2,ng)*d2
 
C--------Node 1 contribution to BCX/BCY due to pinching
         bcx=      bcp(ep,1)*vpinchxyz(ep,1)       
         bcy=      bcp(ep,2)*vpinchxyz(ep,1)
 
C--------Node 2 contribution to BCX/BCY due to pinching
         bcx=bcx + bcp(ep,3)*vpinchxyz(ep,2)       
         bcy=bcy + bcp(ep,4)*vpinchxyz(ep,2)
 
C--------Node 3 contribution to BCX/BCY due to pinching
         bcx=bcx + bcp(ep,5)*vpinchxyz(ep,3)       
         bcy=bcy + bcp(ep,6)*vpinchxyz(ep,3)
 
C--------Node 4 contribution to BCX/BCY due to pinching
         bcx=bcx + bcp(ep,7)*vpinchxyz(ep,4)       
         bcy=bcy + bcp(ep,8)*vpinchxyz(ep,4)
 
C        calculate CT rate because of updates in BCX and BCY     
 
         vdefpinch(ep,1)= tc(ep,1,1)*bcx+tc(ep,2,1)*bcx 
         vdefpinch(ep,2)= tc(ep,1,2)*bcx+tc(ep,2,2)*bcy
 
C----------------------------------
C     MATRIX [BP] FOR PINCHING
C----------------------------------
 
C--------even though constant for all elements, 
C        if B_gp = Ni beta_i does not work then will need it here
C        like MITC4 calculation of CT at the edge centre.
         bp(ep,1) = tnpg(ep,1,ng)
         bp(ep,2) = tnpg(ep,2,ng)
         bp(ep,3) = tnpg(ep,3,ng)
         bp(ep,4) = tnpg(ep,4,ng)
 
         vdefpinch(ep,3) = bp(ep,1)*vpinchxyz(ep,1) 
     .                   + bp(ep,2)*vpinchxyz(ep,2)
     .                   + bp(ep,3)*vpinchxyz(ep,3)
     .                   + bp(ep,4)*vpinchxyz(ep,4)
 
C
C----------------------------------
C     CALCULATE DBeta_x/Dx, Dbeta_y/Dx, Dbeta_x/dy and Dbeta_y/dy at this Gp
C----------------------------------
C-------Calculating dN/dx and dN/dy at this gp
        dn1dx = vksi(1,ng)*tc(ep,1,1)+veta(1,ng)*tc(ep,2,1)
        dn1dy = vksi(1,ng)*tc(ep,1,2)+veta(1,ng)*tc(ep,2,2)
        dn2dx = vksi(2,ng)*tc(ep,1,1)+veta(2,ng)*tc(ep,2,1)
        dn2dy = vksi(2,ng)*tc(ep,1,2)+veta(2,ng)*tc(ep,2,2)
        dn3dx = vksi(3,ng)*tc(ep,1,1)+veta(3,ng)*tc(ep,2,1)
        dn3dy = vksi(3,ng)*tc(ep,1,2)+veta(3,ng)*tc(ep,2,2)
        dn4dx = vksi(4,ng)*tc(ep,1,1)+veta(4,ng)*tc(ep,2,1)
        dn4dy = vksi(4,ng)*tc(ep,1,2)+veta(4,ng)*tc(ep,2,2)
C
C       dbeta_x/dx at this gp
        dbetadxy(ep,1) =  dn1dx*vpincht1(ep,1)+dn2dx*vpincht1(ep,2)
     -                     +dn3dx*vpincht1(ep,3)+dn4dx*vpincht1(ep,4)
C       dbeta_y/dy at this gp
        dbetadxy(ep,2) =  dn1dy*vpincht2(ep,1)+dn2dy*vpincht2(ep,2)
     -                     +dn3dy*vpincht2(ep,3)+dn4dy*vpincht2(ep,4)
C       0.5*(dbeta_x/dy +dbeta_y/dx) at this gp
        dbetadxy(ep,3) =  half*(dn1dy*vpincht1(ep,1)+dn2dy*vpincht1(ep,2)
     -                     +dn3dy*vpincht1(ep,3)+dn4dy*vpincht1(ep,4)
     -                     +dn1dx*vpincht2(ep,1)+dn2dx*vpincht2(ep,2)
     -                     +dn3dx*vpincht2(ep,3)+dn4dx*vpincht2(ep,4))
C
        bpinchdamp(ep,1) = dn1dx
        bpinchdamp(ep,2) = dn1dy
        bpinchdamp(ep,3) = dn2dx
        bpinchdamp(ep,4) = dn2dy
        bpinchdamp(ep,5) = dn3dx
        bpinchdamp(ep,6) = dn3dy
        bpinchdamp(ep,7) = dn4dx
        bpinchdamp(ep,8) = dn4dy
C
       ENDDO
       RETURN