OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
gruneisen_mod Module Reference

Functions/Subroutines

subroutine gruneisen (iflag, nel, pmin, off, eint, mu, mu2, espe, dvol, df, vnew, rho0, pnew, dpdm, dpde, psh, eos_struct)

Function/Subroutine Documentation

◆ gruneisen()

subroutine gruneisen_mod::gruneisen ( integer, intent(in) iflag,
integer, intent(in) nel,
intent(in) pmin,
dimension(nel), intent(in) off,
dimension(nel), intent(inout) eint,
dimension(nel), intent(in) mu,
dimension(nel), intent(in) mu2,
dimension(nel), intent(in) espe,
dimension(nel), intent(in) dvol,
dimension(nel), intent(in) df,
dimension(nel), intent(in) vnew,
dimension(nel), intent(in) rho0,
dimension(nel), intent(inout) pnew,
dimension(nel), intent(inout) dpdm,
dimension(nel), intent(inout) dpde,
dimension(nel), intent(inout) psh,
type(eos_param_), intent(in) eos_struct )
Parameters
[in]iflagflag for the operation (staggered scheme IFLAG=0,1 ; collocated scheme IFLAG=2)
[in]nelnumber of elements in the current group
[in]eos_structdata structure for EoS parameters

Definition at line 41 of file gruneisen.F.

46C-----------------------------------------------
47C M o d u l e s
48C-----------------------------------------------
49 USE eos_param_mod , ONLY : eos_param_
50 USE precision_mod, ONLY : wp
51C-----------------------------------------------
52C D e s c r i p t i o n
53C-----------------------------------------------
54!----------------------------------------------------------------------------
55!! \details STAGGERED SCHEME IS EXECUTED IN TWO PASSES IN EOSMAIN : IFLG=0 THEN IFLG=1
56!! \details COLLOCATED SCHEME IS DOING A SINGLE PASS : IFLG=2
57!! \details
58!! \details STAGGERED SCHEME
59!! \details EOSMAIN / IFLG = 0 : DERIVATIVE CALCULATION FOR SOUND SPEED ESTIMATION c[n+1] REQUIRED FOR PSEUDO-VISCOSITY (DPDE:partial derivative, DPDM:total derivative)
60!! \details MQVISCB : PSEUDO-VISCOSITY Q[n+1]
61!! \details MEINT : INTERNAL ENERGY INTEGRATION FOR E[n+1] : FIRST PART USING P[n], Q[n], and Q[n+1] CONTRIBUTIONS
62!! \details EOSMAIN / IFLG = 1 : UPDATE P[n+1], T[N+1]
63!! \details INTERNAL ENERGY INTEGRATION FOR E[n+1] : LAST PART USING P[n+1] CONTRIBUTION
64!! \details (second order integration dE = -P.dV where P = 0.5(P[n+1] + P[n]) )
65!! \details COLLOCATED SCHEME
66!! \details EOSMAIN / IFLG = 2 : SINGLE PASS FOR P[n+1] AND DERIVATIVES
67!----------------------------------------------------------------------------
68C-----------------------------------------------
69C I m p l i c i t T y p e s
70C-----------------------------------------------
71#include "implicit_f.inc"
72C-----------------------------------------------
73C C o m m o n B l o c k s
74C-----------------------------------------------
75C
76C-----------------------------------------------
77C D u m m y A r g u m e n t s
78C-----------------------------------------------
79 INTEGER,INTENT(IN) :: IFLAG !< flag for the operation (staggered scheme IFLAG=0,1 ; collocated scheme IFLAG=2)
80 INTEGER,INTENT(IN) :: NEL !< number of elements in the current group
81 my_real, INTENT(IN) :: pmin !< minimum pressure allowed
82 my_real,INTENT(IN) :: off(nel) !< OFF=1 if element is active, OFF=0 if element is inactive
83 my_real,INTENT(IN) :: mu(nel) !< volumetric strain : rho/rho0-1
84 my_real,INTENT(IN) :: mu2(nel) !< MU**2 IF MU > 0 , 0 otherwise
85 my_real,INTENT(IN) :: espe(nel) !< specific internal energy
86 my_real,INTENT(IN) :: dvol(nel) !< volume change
87 my_real,INTENT(IN) :: vnew(nel) !< current volume
88 my_real,INTENT(IN) :: rho0(nel) !< reference density
89 my_real,INTENT(IN) :: df(nel) !< relative volume (v/v0 = rho0/rho)
90 my_real,INTENT(INOUT) :: eint(nel) !< internal energy
91 my_real,INTENT(INOUT) :: psh(nel) !< pressure shift (for relative pressure modeling)
92 my_real,INTENT(INOUT) :: pnew(nel) !< current pressure
93 my_real,INTENT(INOUT) :: dpde(nel) !< partial derivative dP/dE where E=Eint/V0
94 my_real,INTENT(INOUT) :: dpdm(nel) !< total derivative : DP/Dmu = dPdmu + dPdE* P/(1+mu)**2
95 TYPE(EOS_PARAM_),INTENT(IN) :: EOS_STRUCT !< data structure for EoS parameters
96C-----------------------------------------------
97C L o c a l V a r i a b l e s
98C-----------------------------------------------
99 INTEGER I
100 real(kind=wp) :: aa, bb, dvv, ff, fg, fac, xx, dff, dfg, fac1, pp
101 real(kind=wp) :: cc,s1,s2,s3,g0,ga
102C-----------------------------------------------
103C B o d y
104C-----------------------------------------------
105 cc = eos_struct%UPARAM(1)
106 s1 = eos_struct%UPARAM(2)
107 s2 = eos_struct%UPARAM(3)
108 s3 = eos_struct%UPARAM(4)
109 g0 = eos_struct%UPARAM(5)
110 ga = eos_struct%UPARAM(6)
111 psh(:) = eos_struct%PSH
112
113 IF(iflag == 0) THEN !--- SOUND SPEED DERIVATIVE (DPDE is partial derivative, DPDM is total derivative)
114 ! Derivative
115 DO i=1,nel
116 fac=one
117 fac1=one
118 IF(mu(i) > zero) THEN
119 xx= mu(i)/(one+mu(i))
120 ff=one+(one-half*g0)*mu(i)-half*ga*mu2(i)
121 fg=one-(s1-one+s2*xx+s3*xx*xx)*mu(i)
122 fac=ff/(fg*fg)
123 dff=one-half*g0-ga*mu(i)
124 dfg=one-s1+xx*(-two*s2+xx*(s2-three*s3)+two*s3*xx*xx)
125 fac1=fac*(one+mu(i)*(dff/ff-two*dfg/fg))
126 ENDIF
127 aa=fac*rho0(i)*cc*cc*mu(i)
128 bb=g0+ga*mu(i)
129 pp=max(aa+bb*espe(i),pmin)*off(i) ! total pressure need for sound speed calculation (=> no pressure shift)
130 dpdm(i)=fac1*rho0(i)*cc*cc+pp*df(i)*df(i)*bb+ga*espe(i)
131 dpde(i)=bb
132 pnew(i) = max(pp,pmin)*off(i)! P(mu[n+1],E[n])
133 pnew(i) = pnew(i) - psh(i)
134 ENDDO
135
136 ELSEIF(iflag == 1) THEN !--- P[n+1] and EINT[n+1]
137 ! Pressure P[n+1]
138 ! Internal Energy Eint[n+1] (2nd order integration requires P[n+1])
139 DO i=1,nel
140 fac=one
141 IF(mu(i) > zero) THEN
142 xx= mu(i)/(one+mu(i))
143 ff=one+(one-half*g0)*mu(i)-half*ga*mu2(i)
144 fg=one-(s1-one+s2*xx+s3*xx*xx)*mu(i)
145 fac=ff/(fg*fg)
146 ENDIF
147 aa=fac*rho0(i)*cc*cc*mu(i)
148 bb=g0+ga*mu(i)
149 dpde(i) = bb
150 dvv=half*dvol(i)*df(i) / max(em15,vnew(i))
151 pnew(i)=(aa+bb*espe(i))/(one+bb*dvv)
152 pnew(i)= max(pnew(i),pmin)*off(i) ! P(mu[n+1],E[n+1])
153 eint(i)=eint(i) - half*dvol(i)*pnew(i) !total pressure for energy integration
154 pnew(i) = pnew(i) - psh(i)
155 ENDDO
156
157 ELSEIF(iflag == 2) THEN !--- collocated scheme (law151)
158 ! pressure and sound speed derivatives
159 DO i=1,nel
160 IF (vnew(i) > zero) THEN
161 fac=one
162 fac1=one
163 IF(mu(i) > zero) THEN
164 xx= mu(i)/(one+mu(i))
165 ff=one+(one-half*g0)*mu(i)-half*ga*mu2(i)
166 fg=one-(s1-one+s2*xx+s3*xx*xx)*mu(i)
167 fac=ff/(fg*fg)
168 dff=one-half*g0-ga*mu(i)
169 dfg=one-s1+xx*(-two*s2+xx*(s2-three*s3)+two*s3*xx*xx)
170 fac1=fac*(one+mu(i)*(dff/ff-two*dfg/fg))
171 ENDIF
172 aa=fac*rho0(i)*cc*cc*mu(i)
173 bb=g0+ga*mu(i)
174 pnew(i)=max(aa+bb*espe(i),pmin)*off(i)
175 dpdm(i)=fac1*rho0(i)*cc*cc+pnew(i)*df(i)*df(i)*bb+ga*espe(i)
176 dpde(i)=bb
177 pnew(i)=pnew(i)-psh(i)
178 ENDIF
179 ENDDO
180
181 ENDIF
182
183 RETURN
my_real
#define my_real
Definition cppsort.cpp:32
max
#define max(a, b)
Definition macros.h:21