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

Go to the source code of this file.

Functions/Subroutines

subroutine jwlun51 (time, xl, tburn, uparam, dd, mu, mup1, volume, dvol, v1old, eint1, viscmax, q1, pext, p1, pm1, rho, rho10, mas1, ssp1, qa, qb, bfrac)
subroutine jwl51 (uparam, v1, v1old, mu1, mup1, eint1, p1old, pext, p1, pm1, rho1, rho10, mas1, ssp1, dvdp1, dpdv1, bfrac, v10, flag)

Function/Subroutine Documentation

◆ jwl51()

subroutine jwl51 ( dimension(*), intent(in) uparam,
intent(in) v1,
intent(in) v1old,
intent(inout) mu1,
mup1,
intent(inout) eint1,
intent(in) p1old,
intent(in) pext,
intent(inout) p1,
intent(in) pm1,
intent(inout) rho1,
intent(in) rho10,
intent(in) mas1,
intent(inout) ssp1,
intent(inout) dvdp1,
intent(inout) dpdv1,
intent(in) bfrac,
intent(in) v10,
integer, intent(in) flag )

Definition at line 143 of file jwl51.F.

147C-----------------------------------------------
148C I m p l i c i t T y p e s
149C-----------------------------------------------
150#include "implicit_f.inc"
151C-----------------------------------------------
152C I N P U T O U T P U T A r g u m e n t s
153C-----------------------------------------------
154 my_real,INTENT(IN) :: v1,v1old,
155 . p1old,pext,pm1,
156 . rho10,mas1,
157 . uparam(*),bfrac, v10
158 INTEGER,INTENT(IN) :: FLAG
159 my_real,INTENT(INOUT) :: rho1, mu1, eint1, p1, ssp1, dvdp1,dpdv1
160C-----------------------------------------------
161C L o c a l V a r i a b l e s
162C-----------------------------------------------
163 INTEGER IBFRAC
164 my_real aa,bb,p0,vdet,bhe,b1,b2,w1,r1,r2,r1m,er1m,r2m,er2m,
165 . mup1,dvdp1i,c11,c01,ssp_products,ssp_unreacted,
166 . psol, pgas, psol_min, pgas_min,dpdmu,dpdv1_reacted,dpdv1_unreacted
167C-----------------------------------------------
168 vdet = uparam(42)
169 bhe = uparam(44)
170 b1 = uparam(45)
171 c01 = uparam(49)
172 c11 = uparam(50)
173 b2 = uparam(51)
174 r1 = uparam(52)
175 r2 = uparam(53)
176 w1 = uparam(54)
177 ibfrac = nint(uparam(68))
178C------------------------
179 dvdp1i = dvdp1
180 rho1 = mas1/v1
181 mup1 = rho1/rho10
182 mu1 = mup1 - one
183
184 r1m = r1/mup1
185 r2m = r2/mup1
186 er1m = exp(-r1m)
187 er2m = exp(-r2m)
188
189
190 aa = w1*mup1/v10 !W1/V1 same digits this way
191 aa = aa
192 bb = half*(v1-v1old)
193 IF (flag == 1) eint1 = eint1 - (p1old+pext+pext)*bb
194 p0 = b1*(one-w1/r1m)*er1m + b2*(one-w1/r2m)*er2m
195 IF (flag == 1) THEN
196 p1 = ( p0-pext + aa*eint1 ) / (one+aa*bb)
197 ELSE
198 p1 = p0 - pext + aa * eint1
199 ENDIF
200
201 !--------------------------------!
202 ! Linear and jwl eos !
203 !--------------------------------!
204 psol = c01+c11*mu1 !linear eos relative pressure
205 psol_min = pm1 !p<0 allowed for solid phase. Default : -EP30
206 psol = max(psol,psol_min)
207
208 pgas = p1 !jwl eos relative to pext
209 pgas_min = -pext !p>0 for detonation products
210 pgas = max(pgas,pgas_min)
211
212 p1 = bfrac*pgas + (one-bfrac)*psol
213 IF (flag == 1) eint1 = eint1 - p1*bb
214 IF (flag == 1) eint1 = max(eint1, zero)
215
216 !--------------------------------!
217 ! Sound Speed !
218 !--------------------------------!
219 dpdmu = b1*er1m*( (-w1*mup1/r1) + r1m - w1) + b2*er2m*( (-w1*mup1/r2) + r2m - w1)
220 . + w1*eint1/v1 + (pgas+pext)*w1
221 dpdmu = abs(dpdmu) / mup1
222 ssp_products = sqrt(dpdmu/rho10)
223 ssp_unreacted = sqrt(c11/rho10)
224 ssp1 = (one-bfrac)*ssp_unreacted + bfrac*ssp_products
225
226 !--------------------------------!
227 ! DPDV !
228 !--------------------------------!
229 dpdv1_reacted = -dpdmu*mup1/v1
230 dpdv1_unreacted = -c11*mup1/v1
231 dpdv1 = bfrac*dpdv1_reacted + (one-bfrac)*dpdv1_unreacted
232
233 IF(abs(dpdv1)<em20)THEN
234 dvdp1 = zero
235 ELSE
236 dvdp1 = one/dpdv1
237 ENDIF
238
239 RETURN
my_real
#define my_real
Definition cppsort.cpp:32
jwl
subroutine jwl(iflag, nel, pm, off, eint, mu, mu2, espe, dvol, df, vnew, mat, psh, pnew, dpdm, dpde)
Definition jwl.F:32
max
#define max(a, b)
Definition macros.h:21

◆ jwlun51()

subroutine jwlun51 ( intent(in) time,
intent(in) xl,
intent(in) tburn,
dimension(*), intent(in) uparam,
intent(in) dd,
intent(inout) mu,
mup1,
intent(in) volume,
intent(inout) dvol,
intent(in) v1old,
intent(in) eint1,
intent(inout) viscmax,
intent(inout) q1,
intent(in) pext,
intent(inout) p1,
intent(in) pm1,
intent(inout) rho,
intent(in) rho10,
intent(in) mas1,
intent(inout) ssp1,
intent(in) qa,
intent(in) qb,
intent(inout) bfrac )

Definition at line 28 of file jwl51.F.

32C-----------------------------------------------
33C I m p l i c i t T y p e s
34C-----------------------------------------------
35#include "implicit_f.inc"
36C-----------------------------------------------
37C I N P U T O U T P U T A r g u m e n t s
38C-----------------------------------------------
39 my_real, INTENT(IN) :: time,xl,tburn,dd,
40 . volume,v1old,eint1,
41 . pext,pm1,
42 . rho10,mas1,
43 . uparam(*),qa,qb
44 my_real, INTENT(INOUT) :: ssp1, p1, dvol, rho, bfrac, mu, viscmax, q1
45C-----------------------------------------------
46C L o c a l V a r i a b l e s
47C-----------------------------------------------
48 INTEGER IBFRAC
49 my_real aa,bb,p0,vdet,bhe,b1,b2,w1,r1,r2,r1m,er1m,r2m,er2m,
50 . qal,qbl,dpdmu,mup1,c01,c11,
51 . psol, pgas, psol_min, pgas_min, ssp_unreacted, ssp_reacted
52C-----------------------------------------------
53 vdet = uparam(42)
54 bhe = uparam(44)
55 b1 = uparam(45)
56 c01 = uparam(49)
57 c11 = uparam(50)
58 b2 = uparam(51)
59 r1 = uparam(52)
60 r2 = uparam(53)
61 w1 = uparam(54)
62 ibfrac= nint(uparam(68))
63C
64 IF(r1 == zero) r1=ep30
65 IF(r2 == zero) r2=ep30
66C
67 dvol = volume - v1old
68C
69 !--------------------------------!
70 ! Calculation of BFRAC in [0,1] !
71 !--------------------------------!
72 rho = mas1 / volume
73 IF(bfrac < one) THEN
74 bfrac = zero
75 IF(ibfrac/=1 .AND. time > -tburn) bfrac = vdet*(time+tburn)*two_third/xl
76 IF(ibfrac/=2) bfrac = max( bfrac , bhe * (one - rho10/rho) )
77 IF(bfrac < em03) THEN
78 bfrac = zero
79 ELSEIF(bfrac > one) THEN
80 bfrac = one
81 ENDIF
82 ENDIF
83
84 !--------------------------------!
85 ! SSP & ARTIFICIAL VISCO !
86 !--------------------------------!
87 mup1 = rho/rho10
88 mu = mup1 - one
89 r1m = r1/mup1
90 r2m = r2/mup1
91 er1m = exp(-r1m)
92 er2m = exp(-r2m)
93 aa = w1/volume
94 p0 = b1*(one-w1/r1m)*er1m + b2*(one-w1/r2m)*er2m
95 p1 = p0 + aa*eint1 !total jwl pressure for ssp
96 dpdmu = b1*er1m*( (-w1*mup1/r1) + r1m - w1) + b2*er2m*( (-w1*mup1/r2) + r2m - w1) + w1*eint1/volume +p1*w1
97 dpdmu = abs(dpdmu) / mup1 ! if DPDMU <0 => numerical error during energy integration (increase iteration number or reduce submaterial volume change ratio)
98 ssp_reacted = sqrt(dpdmu/rho10)
99 ssp_unreacted = sqrt(c11/rho10)
100 ssp1 = max(bfrac*ssp_reacted,(one-bfrac)*ssp_unreacted)
101 qal = qa*xl
102 qal = qal*qal
103 qbl = qb*xl
104 viscmax = rho*(qal*max(zero,dd) + qbl*ssp1)
105 q1 = viscmax*max(zero,dd)
106 bb = half*(volume-v1old)
107! EINT1 = EINT1 - (P1OLD+PEXT+PEXT)*BB
108 aa = aa
109 p1 = ( p0-pext + aa*eint1 )! / (ONE+AA*BB)
110
111
112 !--------------------------------!
113 ! Linear and jwl eos !
114 !--------------------------------!
115 psol = c01+c11*mu !linear eos relative pressure
116 psol_min = pm1 !p<0 allowed for solid phase. Default : -EP30
117 psol = max(psol,psol_min)
118
119 pgas = p1 !jwl eos relative to Pext
120 pgas_min = -pext !p>0 for detonation products
121 pgas = max(pgas,pgas_min)
122
123 p1 = bfrac*pgas + (one-bfrac)*psol
124
125 !--------------------------------!
126 ! Update SSP with current state !
127 !--------------------------------!
128 dpdmu = b1*er1m*( (-w1*mup1/r1) + r1m - w1) + b2*er2m*( (-w1*mup1/r2) + r2m - w1) + w1*eint1/volume +(p1+pext)*w1
129 dpdmu = abs(dpdmu) / mup1 ! if DPDMU <0 => numerical error during energy integration (increase iteration number or reduce submaterial volume change ratio)
130 ssp_reacted = sqrt(dpdmu/rho10)
131 ssp_unreacted = sqrt(c11/rho10)
132 ssp1 = max(bfrac*ssp_reacted,(one-bfrac)*ssp_unreacted)
133C
134 RETURN