m37init.F File Reference

#include "implicit_f.inc"
#include "com04_c.inc"
#include "param_c.inc"
#include "vect01_c.inc"
#include "mvsiz_p.inc"

Go to the source code of this file.

Functions/Subroutines
subroutine	m37init (ipm, pm, ngl, nuvar, uvar, uparam, x, mat, iparg, iform, ix, nix, iloadp, facload, gbuf, nel)

Function/Subroutine Documentation

m37init()

subroutine m37init	(	integer, dimension(npropmi,nummat)	ipm,
			pm,
		integer, dimension(mvsiz)	ngl,
		integer	nuvar,
			uvar,
			uparam,
			x,
		integer, dimension(mvsiz)	mat,
		integer, dimension(nparg)	iparg,
		integer	iform,
		integer, dimension(nix,*)	ix,
		integer	nix,
		integer, dimension(sizloadp,*), intent(in)	iloadp,
		dimension(lfacload,*), intent(in)	facload,
		type(g_bufel_), intent(inout), target	gbuf,
		integer, intent(in)	nel )

Definition at line 29 of file m37init.F.

C-----------------------------------------------
C   D e s c r i p t i o n
C-----------------------------------------------
C This subroutine is initializing cell based on law37
C In case of /LOAD/HYDRO is attached to it then
C   Pressure is computed so that P(t=0) = P0 + DP 
C     where DP=rho0.g.h
C     and h is signed distance from water/air free surface
C     which is defined by B and n in input
C-----------------------------------------------
C   M o d u l e s
C-----------------------------------------------
      USE elbufdef_mod 
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      "com04_c.inc"
#include      "param_c.inc"
#include      "vect01_c.inc"
#include      "mvsiz_p.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER :: IPM(NPROPMI,NUMMAT),NGL(MVSIZ)      ,MAT(MVSIZ), IPARG(NPARG)   ,IFORM,NIX,IX(NIX,*)
      my_real :: pm(npropm,nummat),x(3,numnod), uvar(llt,nuvar)
      INTEGER,INTENT(IN) :: ILOADP(SIZLOADP,*),NEL
      my_real,INTENT(IN) :: facload(lfacload,*)
      TYPE(G_BUFEL_), INTENT(INOUT),TARGET :: GBUF  
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER  :: I,J,NUVAR, GG1, GG2, GG3, ISFLUID,NL
      INTEGER  :: NPH,IFLG,IE,IADBUF,NUPARAM, M_UID,M_IID,IGRAV
      my_real  :: mfrac, a, mu1p1, mu1p2, rho1, rho2, rho10, rho20, c1, gam, pmin, p0, mu2p1, rho,p,a1,psh,depth,z(3)
      my_real  :: uparam(*),xbas,ybas,zbas,nx,ny,nz,grav, bz(3),dotprod, delta_p, rho0_liq,rho0_gas, rho_liq,rho_gas
      my_real  :: pold,height
      my_real  :: n(3,8)
C-----------------------------------------------
 
      !---------------------------------!
      ! CELL INIT.   (UVAR ARRAY)       !
      !---------------------------------!
      !UVAR(I,1) : massic percentage of liquid * global density  (rho1*V1/V : it needs to give liquid mass multiplying by element volume in aleconve.F)
      !UVAR(I,2) : density of gas
      !UVAR(I,3) : density of liquid
      !UVAR(I,4) : volumetric fraction of liquid
      !UVAR(I,5) : volumetric fraction of gas
      DO i=lft,llt
        ie        =  i+nft
        c1        =  uparam(4)
        gam       =  uparam(5)
        pmin      =  uparam(8)
        p0        =  uparam(9)
        rho10     =  uparam(11)
        rho20     =  uparam(12)
        a1        =  uparam(10)         
        p         = -p0  
        rho       =  rho10 * a1 + (one-a1)*rho20     
        IF(gam*c1>=em30)THEN !if Liquid and gas correctly defined
          mu1p1     = one       !(P/C1+ONE)
          mu2p1     = one       !( ONE+P/P0)**(ONE/GAM)
          rho1      = rho10*mu1p1
          rho2      = rho20*mu2p1
          a         = (rho-rho2)/(rho1-rho2)
          uvar(i,1) = a*rho1
          uvar(i,2) = rho2
          uvar(i,3) = rho1
          uvar(i,4) = a
          IF(uvar(i,4)<em20)uvar(i,4)=zero
          uvar(i,5) = one-uvar(i,4)
          !SOUNDSP(I)  = A*SQRT(C1/RHO1) + (ONE-A) * SQRT(GAM*P0/RHO2)  -> check sigeps37.F (starter)
        ELSE !boundary element
          uvar(i,3) = rho
        ENDIF
      ENDDO
 
C-----------------------------------------------
      RETURN