s4derit3.F File Reference

#include "implicit_f.inc"
#include "comlock.inc"
#include "mvsiz_p.inc"
#include "param_c.inc"
#include "scr17_c.inc"

Go to the source code of this file.

Functions/Subroutines
subroutine	s4derit3 (off, det, ngl, deltax, mxt, x1, x2, x3, x4, y1, y2, y3, y4, z1, z2, z3, z4, px1, px2, px3, px4, py1, py2, py3, py4, pz1, pz2, pz3, pz4, rx, ry, rz, sx, sy, sz, tx, ty, tz, pm, voldp, nel, iformdt)

Function/Subroutine Documentation

s4derit3()

subroutine s4derit3	(		off,
			det,
		integer, dimension(*)	ngl,
			deltax,
		integer, dimension(*)	mxt,
		double precision, dimension(*)	x1,
		double precision, dimension(*)	x2,
		double precision, dimension(*)	x3,
		double precision, dimension(*)	x4,
		double precision, dimension(*)	y1,
		double precision, dimension(*)	y2,
		double precision, dimension(*)	y3,
		double precision, dimension(*)	y4,
		double precision, dimension(*)	z1,
		double precision, dimension(*)	z2,
		double precision, dimension(*)	z3,
		double precision, dimension(*)	z4,
			px1,
			px2,
			px3,
			px4,
			py1,
			py2,
			py3,
			py4,
			pz1,
			pz2,
			pz3,
			pz4,
			rx,
			ry,
			rz,
			sx,
			sy,
			sz,
			tx,
			ty,
			tz,
			pm,
		double precision, dimension(*)	voldp,
		integer, intent(in)	nel,
		integer, intent(in)	iformdt )

Definition at line 31 of file s4derit3.F.

C-----------------------------------------------
C   I m p l i c i t   T y p e s
C-----------------------------------------------
#include      "implicit_f.inc"
#include      "comlock.inc"
C-----------------------------------------------
C   G l o b a l   P a r a m e t e r s
C-----------------------------------------------
#include      "mvsiz_p.inc"
#include      "param_c.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER, INTENT(IN) :: NEL
      INTEGER, INTENT(IN) :: IFORMDT
      INTEGER MXT(*)
      double precision
     .   x1(*), x2(*), x3(*), x4(*),
     .   y1(*), y2(*), y3(*), y4(*),
     .   z1(*), z2(*), z3(*), z4(*),voldp(*)
      my_real
     .   off(*),det(*),deltax(*),
     .   px1(*), px2(*), px3(*), px4(*),  
     .   py1(*), py2(*), py3(*), py4(*),  
     .   pz1(*), pz2(*), pz3(*), pz4(*), 
     .   rx(*), ry(*), rz(*), sx(*), sy(*), sz(*),tx(*), ty(*), tz(*),
     .   pm(npropm,*)
C-----------------------------------------------
C   C o m m o n   B l o c k s
C-----------------------------------------------
#include      "scr17_c.inc"
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER NGL(*), I
      double precision
     .   x41, y41, z41, x42, y42, z42, x43, y43, z43,b1dp,c1dp,d1dp
      my_real
     .   b1(mvsiz), b2(mvsiz), b3(mvsiz), b4(mvsiz), 
     .   c1(mvsiz), c2(mvsiz), c3(mvsiz), c4(mvsiz),
     .   d1(mvsiz), d2(mvsiz), d3(mvsiz), d4(mvsiz),
     .       d,
     .   pxx, pyy, pzz, pxy, pyz, pxz, gfac, aa, bb, p, ld
C-----------------------------------------------
C
      DO i=1,nel
       x43 = x4(i) - x3(i)
       y43 = y4(i) - y3(i)
       z43 = z4(i) - z3(i)
       x41 = x4(i) - x1(i)
       y41 = y4(i) - y1(i)
       z41 = z4(i) - z1(i)
       x42 = x4(i) - x2(i)
       y42 = y4(i) - y2(i)
       z42 = z4(i) - z2(i)
C
       rx(i) =  -x41
       ry(i) =  -y41
       rz(i) =  -z41
       sx(i) =  -x42
       sy(i) =  -y42
       sz(i) =  -z42
C
       tx(i) =  -x43
       ty(i) =  -y43
       tz(i) =  -z43
C
       b1dp  =  y43*z42 - y42*z43
       b1(i) =  b1dp
       b2(i) =  y41*z43 - y43*z41
       b3(i) =  y42*z41 - y41*z42
       b4(i) =  -(b1(i) + b2(i) + b3(i))
C
       c1dp  =  z43*x42 - z42*x43
       c1(i) =  c1dp
       c2(i) =  z41*x43 - z43*x41
       c3(i) =  z42*x41 - z41*x42
       c4(i) =  -(c1(i) + c2(i) + c3(i))
C
       d1dp  =  x43*y42 - x42*y43
       d1(i) =  d1dp
       d2(i) =  x41*y43 - x43*y41
       d3(i) =  x42*y41 - x41*y42
       d4(i) =  -(d1(i) + d2(i) + d3(i))
C
       voldp(i) = (x41*b1dp + y41*c1dp + z41*d1dp)*one_over_6
       det(i) = voldp(i)
C
      ENDDO
C
      CALL schkjab3(
     1   off,     det,     ngl,     nel)
C
      DO i=1,nel
        d = one/det(i)/six
        px1(i)=-b1(i)*d
        py1(i)=-c1(i)*d
        pz1(i)=-d1(i)*d
        px2(i)=-b2(i)*d
        py2(i)=-c2(i)*d
        pz2(i)=-d2(i)*d
        px3(i)=-b3(i)*d
        py3(i)=-c3(i)*d
        pz3(i)=-d3(i)*d
        px4(i)=-b4(i)*d
        py4(i)=-c4(i)*d
        pz4(i)=-d4(i)*d
      END DO
 
      IF(idt1sol==0)THEN
 
        DO i=1,nel
          d = max(px1(i)*px1(i)+py1(i)*py1(i)+pz1(i)*pz1(i),
     .      px2(i)*px2(i)+py2(i)*py2(i)+pz2(i)*pz2(i),
     .      px3(i)*px3(i)+py3(i)*py3(i)+pz3(i)*pz3(i),
     .      px4(i)*px4(i)+py4(i)*py4(i)+pz4(i)*pz4(i))
          deltax(i) = one / sqrt(d)
        END DO
 
 
       ELSEIF(iformdt==0)THEN
         DO i=1,nel
             d = px1(i)*px1(i)+py1(i)*py1(i)+pz1(i)*pz1(i)
     .       + px2(i)*px2(i)+py2(i)*py2(i)+pz2(i)*pz2(i)
     .       + px3(i)*px3(i)+py3(i)*py3(i)+pz3(i)*pz3(i)
     .       + px4(i)*px4(i)+py4(i)*py4(i)+pz4(i)*pz4(i)
           deltax(i) = one / sqrt(d)
         END DO
 
      ELSEIF(iformdt==1)THEN
 
         gfac=pm(105,mxt(1))
         ld  =two*sqrt(max(one-gfac,zero))+one
         DO i=1,nel
           pxx=px1(i)*px1(i)+px2(i)*px2(i)+px3(i)*px3(i)+px4(i)*px4(i)
           pyy=py1(i)*py1(i)+py2(i)*py2(i)+py3(i)*py3(i)+py4(i)*py4(i)
           pzz=pz1(i)*pz1(i)+pz2(i)*pz2(i)+pz3(i)*pz3(i)+pz4(i)*pz4(i)
           pxy=px1(i)*py1(i)+px2(i)*py2(i)+px3(i)*py3(i)+px4(i)*py4(i)
           pxz=px1(i)*pz1(i)+px2(i)*pz2(i)+px3(i)*pz3(i)+px4(i)*pz4(i)
           pyz=py1(i)*pz1(i)+py2(i)*pz2(i)+py3(i)*pz3(i)+py4(i)*pz4(i)
C
           aa = -(pxx+pyy+pzz)
           bb =  (pxx*pyy+pxx*pzz+pyy*pzz-pxy**2-pxz**2-pyz**2) 
           p  = bb-third*aa*aa
           d=two*sqrt(third*max(-p,zero))-third*aa
C
             d=ld*d
C
           deltax(i) = one / sqrt(d)
         END DO
 
       ELSEIF(iformdt==2)THEN
 
         gfac=pm(105,mxt(1))
         DO i=1,nel
           pxx=px1(i)*px1(i)+px2(i)*px2(i)+px3(i)*px3(i)+px4(i)*px4(i)
           pyy=py1(i)*py1(i)+py2(i)*py2(i)+py3(i)*py3(i)+py4(i)*py4(i)
           pzz=pz1(i)*pz1(i)+pz2(i)*pz2(i)+pz3(i)*pz3(i)+pz4(i)*pz4(i)
           pxy=px1(i)*py1(i)+px2(i)*py2(i)+px3(i)*py3(i)+px4(i)*py4(i)
           pxz=px1(i)*pz1(i)+px2(i)*pz2(i)+px3(i)*pz3(i)+px4(i)*pz4(i)
           pyz=py1(i)*pz1(i)+py2(i)*pz2(i)+py3(i)*pz3(i)+py4(i)*pz4(i)
C
           aa = -(pxx+pyy+pzz)
           bb =  gfac*(pxx*pyy+pxx*pzz+pyy*pzz-pxy**2-pxz**2-pyz**2) 
           p  = bb-third*aa*aa
           d  = two*sqrt(third*max(-p,zero))-third*aa
C
           deltax(i) = one / sqrt(d)
         END DO
 
      END IF
C
      RETURN
C