init_qd.F File Reference

#include "implicit_f.inc"
#include "units_c.inc"

Go to the source code of this file.

Functions/Subroutines
subroutine	init_qd (iflow, ibuf, elem, x, xs, ys, zs, xd, yd, zd, rflow, normal, ta, af, cosg, dcp)

Function/Subroutine Documentation

init_qd()

subroutine init_qd	(	integer, dimension(*)	iflow,
		integer, dimension(*)	ibuf,
		integer, dimension(5,*)	elem,
			x,
			xs,
			ys,
			zs,
			xd,
			yd,
			zd,
			rflow,
			normal,
			ta,
			af,
			cosg,
			dcp )

Definition at line 28 of file init_qd.F.

C-----------------------------------------------
C   I m p l i c i t   T y p e s
C-----------------------------------------------
#include      "implicit_f.inc"
C-----------------------------------------------
C   C o m m o n   B l o c k s
C-----------------------------------------------
#include      "units_c.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER IFLOW(*), IBUF(*), ELEM(5,*)
      my_real x(3,*), af(*), rflow(*), normal(3,*), ta(*), cosg(*), dcp(*)
      my_real xs, ys, zs, xd, yd, zd
C--------------------------------------------------------------
C   L o c a l   V a r i a b l e s
C--------------------------------------------------------------
      INTEGER ILVOUT, NEL, IWAVE, FREESURF
      INTEGER N1, N2, N3, N4, N5, NN1, NN2, NN3, NN4, IEL, JEL
      my_real x1, y1, z1, x2, y2, z2, x3, y3, z3, x4, y4, z4,
     .        xp, yp, zp, x13, y13, z13, x24, y24, z24,    
     .        nrx, nry, nrz, area2, dcs, ssp, wi(4,2),
     .        xc, yc, zc, dirx, diry, dirz
C
      ilvout  = iflow(17)
      nel     = iflow(6)
      iwave   = iflow(22)
      freesurf= iflow(25)
C-------------------------------------------------------------------
C Compute Area, Normal and Arrival Time
C-------------------------------------------------------------------
      wi(1,1)=fourth
      wi(2,1)=fourth
      wi(3,1)=fourth
      wi(4,1)=fourth
      wi(1,2)=third
      wi(2,2)=third
      wi(3,2)=one_over_6
      wi(4,2)=one_over_6
      ssp=rflow(2)    
      IF(iwave==1) THEN
        xc =rflow(9)    
        yc =rflow(10)    
        zc =rflow(11)
        dcs=rflow(12)
      ELSEIF(iwave==2) THEN
        dirx=rflow(9)
        diry=rflow(10)
        dirz=rflow(11)
      ENDIF
            
      DO iel=1,nel
         n1=elem(1,iel)
         n2=elem(2,iel)
         n3=elem(3,iel)
         n4=elem(4,iel)
         n5=elem(5,iel)
         nn1=ibuf(n1)
         nn2=ibuf(n2)
         nn3=ibuf(n3)
         nn4=ibuf(n4)
         x1=x(1,nn1)
         x2=x(1,nn2)
         x3=x(1,nn3)
         x4=x(1,nn4)
         y1=x(2,nn1)
         y2=x(2,nn2)
         y3=x(2,nn3)
         y4=x(2,nn4)
         z1=x(3,nn1)
         z2=x(3,nn2)
         z3=x(3,nn3)
         z4=x(3,nn4)
C  Normale
         x13=x3-x1
         y13=y3-y1
         z13=z3-z1
         x24=x4-x2
         y24=y4-y2
         z24=z4-z2
         nrx=y13*z24-z13*y24
         nry=z13*x24-x13*z24
         nrz=x13*y24-y13*x24
         area2=sqrt(nrx**2+nry**2+nrz**2)
         normal(1,iel)=nrx/area2
         normal(2,iel)=nry/area2
         normal(3,iel)=nrz/area2
         af(iel)=half*area2
C  Centroid
         xp=wi(1,n5)*x1+wi(2,n5)*x2+wi(3,n5)*x3+wi(4,n5)*x4
         yp=wi(1,n5)*y1+wi(2,n5)*y2+wi(3,n5)*y3+wi(4,n5)*y4
         zp=wi(1,n5)*z1+wi(2,n5)*z2+wi(3,n5)*z3+wi(4,n5)*z4
C  Arrival time
         IF(iwave==1) THEN
            dcp(iel) = sqrt((xp-xc)**2+(yp-yc)**2+(zp-zc)**2)
            ta(iel)  = (dcp(iel)-dcs)/ssp
            cosg(iel)= (nrx*(xp-xc)+nry*(yp-yc)+nrz*(zp-zc))/(area2*dcp(iel))
            IF(freesurf == 2) THEN
               jel=iel+nel
               dcp(jel) = sqrt((xp-xd)**2+(yp-yd)**2+(zp-zd)**2)
               ta(jel)  = (dcp(jel)-dcs)/ssp
               cosg(jel)= (nrx*(xp-xd)+nry*(yp-yd)+nrz*(zp-zd))/(area2*dcp(jel))
            ENDIF
         ELSEIF(iwave==2) THEN
            dcp(iel) = (xp-xs)*dirx+(yp-xs)*diry+(zp-zs)*dirz
            ta(iel)  = dcp(iel)/ssp
            cosg(iel)= (nrx*dirx+nry*diry+nrz*dirz)/area2
         ENDIF
      ENDDO 
 
      IF(ilvout>=2) THEN
C---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----1----+----2----+----3--
         WRITE (iout,'(//7X,2A)') 'ELEMENT   ARRIVAL TIME        AREA     DISTANCE   DIR.COSINE     ',
     .                            'NORMAL-X     NORMAL-Y     NORMAL-Z'
         DO iel = 1,nel
            WRITE (iout,'(5X,I10,7E13.5)')iel,ta(iel),af(iel),dcp(iel),cosg(iel),normal(1,iel),normal(2,iel),normal(3,iel)
         ENDDO
         IF(freesurf == 2) THEN
            WRITE (iout,'(//7X,2A)') 'ELEMENT   ARRIVAL TIME        AREA     DISTANCE   DIR.COSINE     '
            DO iel = 1,nel
               jel=iel+nel
               WRITE (iout,'(5X,I10,7E13.5)')iel,ta(jel),af(iel),dcp(jel),cosg(jel)
            ENDDO
         ENDIF
      ENDIF
            
      RETURN