OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
hm_read_eos_gruneisen.F File Reference
#include "implicit_f.inc"
#include "param_c.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine hm_read_eos_gruneisen (iout, pm, unitab, iunit, lsubmodel, imideos, mat_param)

Function/Subroutine Documentation

◆ hm_read_eos_gruneisen()

subroutine hm_read_eos_gruneisen ( integer iout,
pm,
type (unit_type_), intent(in) unitab,
integer iunit,
type(submodel_data), dimension(nsubmod), intent(in) lsubmodel,
integer, intent(in) imideos,
type(matparam_struct_), intent(inout) mat_param )
Parameters
[in]imideospressure shift
[in]imideosinitial pressure

Definition at line 35 of file hm_read_eos_gruneisen.F.

36C-----------------------------------------------
37C M o d u l e s
38C-----------------------------------------------
39 USE unitab_mod
40 USE submodel_mod
41 USE message_mod
42 USE matparam_def_mod, ONLY : matparam_struct_
43C-----------------------------------------------
44C D e s c r i p t i o n
45C-----------------------------------------------
46C reading parameters for
47C MIE-GRUNEISEN EQUATION OF STATE
48C-----------------------------------------------
49C C o m m e n t s
50C-----------------------------------------------
51C RHOI = PM(89) -> provided by /MAT
52C RHOR = PM(01) -> provided by /MAT (can be erased by EOS if present : obsolete)
53C => MU0 = RHO/RHOR-1.
54C PM(31) = P(MU0,E0) -> will be used to initialize diagonal of stress tensor SIG(1:3,*)
55C-----------------------------------------------
56C I m p l i c i t T y p e s
57C-----------------------------------------------
58#include "implicit_f.inc"
59C-----------------------------------------------
60C D u m m y A r g u m e n t s
61C-----------------------------------------------
62 TYPE (UNIT_TYPE_),INTENT(IN) ::UNITAB
63 INTEGER IOUT,IUNIT
64 my_real pm(npropm)
65 TYPE(SUBMODEL_DATA), DIMENSION(NSUBMOD), INTENT(IN) :: LSUBMODEL
66 TYPE(MATPARAM_STRUCT_),INTENT(INOUT) :: MAT_PARAM
67 INTEGER,INTENT(IN) :: IMIDEOS
68C-----------------------------------------------
69C C o m m o n B l o c k s
70C-----------------------------------------------
71#include "param_c.inc"
72C-----------------------------------------------
73C L o c a l V a r i a b l e s
74C-----------------------------------------------
75 my_real :: c, s1, s2, s3, gama0, a, e0, rho0,rhoi,rhor
76 my_real :: mu2,fac1,dpdmu,pp,bb,aa
77 my_real :: mu0, df, ssp0, g0, fac, ff, fg, xx, dff, dfg
78 my_real :: psh !< pressure shift
79 my_real :: p0 !< initial pressure
80 LOGICAL :: IS_ENCRYPTED, IS_AVAILABLE, IS_AVAILABLE_RHO0
81C-----------------------------------------------
82C S o u r c e L i n e s
83C-----------------------------------------------
84 is_encrypted = .false.
85 is_available = .false.
86 is_available_rho0 = .false.
87
88 CALL hm_option_is_encrypted(is_encrypted)
89
90 CALL hm_get_floatv('MAT_C', c, is_available,lsubmodel,unitab)
91 CALL hm_get_floatv('S1', s1, is_available,lsubmodel,unitab)
92 CALL hm_get_floatv('S2', s2, is_available,lsubmodel,unitab)
93 CALL hm_get_floatv('S3', s3, is_available,lsubmodel,unitab)
94
95 CALL hm_get_floatv('GAMMA',gama0, is_available,lsubmodel,unitab)
96 CALL hm_get_floatv('MAT_A', a, is_available,lsubmodel,unitab)
97 CALL hm_get_floatv('R0E0', e0, is_available,lsubmodel,unitab)
98 CALL hm_get_floatv('P0', p0, is_available,lsubmodel,unitab)
99 CALL hm_get_floatv('PSH', psh, is_available,lsubmodel,unitab)
100 CALL hm_get_floatv('Refer_Rho', rho0 ,is_available_rho0,lsubmodel,unitab)
101!-------------------------------------------------------------------------------
102 IF(a == zero) a=gama0
103
104
105 IF(p0 < zero)THEN
106 CALL ancmsg(msgid=67,msgtype=msgerror,anmode=aninfo,i1=imideos,
107 . c1='/EOS/GRUNEISEN',
108 . c2='INITIAL PRESSURE MUST BE STRICTLY POSITIVE (TOTAL PRESSURE). USE PSH PARAMETER TO SHIFT THE PRESSURE')
109 ENDIF
110
111 IF(p0 > zero .AND. e0 /= zero)THEN
112 CALL ancmsg(msgid=67,msgtype=msgerror,anmode=aninfo,i1=imideos,
113 . c1='/EOS/GRUNEISEN',
114 . c2='INITIAL PRESSURE PROVIDED. E0 IS CONSEQUENTLY REDEFINED SUCH AS P(RHO0,E0)=P0')
115 ENDIF
116
117 rhor = pm(1)
118 rhoi = pm(89)
119
120 IF(rho0 > zero) THEN
121 rhor = rho0
122 pm(1)= rho0
123 mat_param%RHO = rho0
124 ELSE
125 rho0=rhor
126 ENDIF
127
128 !COMPUTE INITIAL PRESSURE FOR PM(31)->SIG(1:3,*)
129 !SSP0
130 IF(rhoi == zero)THEN
131 mu0 = zero ! error 683 already displayed
132 ELSE
133 IF(rhor /= zero)THEN
134 mu0 = rhoi/rhor-one
135 ELSE
136 mu0 = zero ! error 683 already displayed
137 ENDIF
138 ENDIF
139 ! Internal Energy (if P0 is provided, calculate E0 consequently)
140 IF(p0 > zero)THEN
141 if(gama0 /= zero)then
142 e0 = (p0-rho0*c*c*mu0)/(gama0+a*mu0)
143 endif
144 ENDIF
145 !Relative volume
146 IF(rhoi /= zero)THEN
147 df = rhor/rhoi
148 ELSE
149 df = zero
150 ENDIF
151 !pressure
152 mu2=mu0*mu0
153 ssp0 = zero
154 g0 = pm(22)
155 rhoi = pm(89)
156 fac=one
157 fac1=one
158 IF(mu0>0)THEN
159 xx= mu0/(one+mu0)
160 ff=one+(one-half*gama0)*mu0-half*a*mu2
161 fg=one-(s1-one+s2*xx+s3*xx*xx)*mu0
162 fac=ff/(fg*fg)
163 dff=one-half*gama0-a*mu0
164 dfg=one-s1+xx*(-two*s2+xx*(s2-three*s3)+two*s3*xx*xx)
165 fac1=fac*(one+mu0*(dff/ff-two*dfg/fg))
166 ENDIF
167 aa=fac*rhor*c*c*mu0
168 bb=gama0+a*mu0
169 pp=max(aa+bb*e0,pm(37))
170 !derivatives and sound speed
171 dpdmu=fac1*rhoi*c*c+pp*df*df*bb+a*e0
172 dpdmu=max(zero,dpdmu)
173 IF(rhor > zero) ssp0 = sqrt((dpdmu + two_third*g0)/rhor)
174
175 ! Storage in Material Buffer PM
176 pm(23) = e0
177 pm(27) = ssp0
178 pm(32) = pm(1)*c*c + a*e0 !BULK = (1+ยต)*dP/dmu (partial derivative at constant E)
179 pm(88) = psh
180 pm(31) = pp - psh
181 pm(104) = pp - psh
182
183 mat_param%EOS%NUPARAM = 6
184 mat_param%EOS%NIPARAM = 0
185 mat_param%EOS%NFUNC = 0
186 mat_param%EOS%NTABLE = 0
187 CALL mat_param%EOS%CONSTRUCT() !allocations
188
189 mat_param%EOS%UPARAM(1) = c
190 mat_param%EOS%UPARAM(2) = s1
191 mat_param%EOS%UPARAM(3) = s2
192 mat_param%EOS%UPARAM(4) = s3
193 mat_param%EOS%UPARAM(5) = gama0
194 mat_param%EOS%UPARAM(6) = a
195 mat_param%EOS%PSH = psh
196 mat_param%EOS%E0 = e0
197 IF (mat_param%THERM%TINI == zero) THEN
198 mat_param%THERM%TINI =three100
199 pm(79) = three100
200 END IF
201
202 WRITE(iout,1000)
203 IF(is_encrypted)THEN
204 WRITE(iout,'(5X,A,//)')'CONFIDENTIAL DATA'
205 ELSE
206 WRITE(iout,1500)c,s1,s2,s3,gama0,a,e0,pp,psh,pp-psh
207 IF(is_available_rho0)WRITE(iout,1501)pm(1)
208 ENDIF
209
210 RETURN
211
212 1000 FORMAT(
213 & 5x,' MIE-GRUNEISEN EOS ',/,
214 & 5x,' ----------------- ',/)
215 1500 FORMAT(
216 & 5x,'C . . . . . . . . . . . . . . . . . . . .=',1pg20.13/,
217 & 5x,'S1. . . . . . . . . . . . . . . . . . . .=',1pg20.13/,
218 & 5x,'S2. . . . . . . . . . . . . . . . . . . .=',1pg20.13/,
219 & 5x,'S3. . . . . . . . . . . . . . . . . . . .=',1pg20.13/,
220 & 5x,'GAMA0 . . . . . . . . . . . . . . . . . .=',1pg20.13/,
221 & 5x,'a . . . . . . . . . . . . . . . . . . . .=',1PG20.13/,
222 & 5X,'initial internal energy per unit volume .=',1PG20.13/,
223 & 5X,'initial pressure . . . . . . . . . . . .=',1PG20.13/,
224 & 5X,'pressure shift . . . . . . . . . . . . .=',1PG20.13/,
225 & 5X,'initial pressure(shifted) . . . . . . .=',1PG20.13)
226 1501 FORMAT(
227 & 5X,'eos reference density . . . . . . . . . .=',1PG20.13)
228
229 RETURN
#define my_real
Definition cppsort.cpp:32
subroutine hm_get_floatv(name, rval, is_available, lsubmodel, unitab)
subroutine hm_option_is_encrypted(is_encrypted)
#define max(a, b)
Definition macros.h:21
subroutine ancmsg(msgid, msgtype, anmode, i1, i2, i3, i4, i5, i6, i7, i8, i9, i10, i11, i12, i13, i14, i15, i16, i17, i18, i19, i20, r1, r2, r3, r4, r5, r6, r7, r8, r9, c1, c2, c3, c4, c5, c6, c7, c8, c9, prmode)
Definition message.F:895