41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
65 USE matparam_def_mod
68
69
70
71#include "implicit_f.inc"
72
73
74
75#include "units_c.inc"
76#include "com01_c.inc"
77#include "param_c.inc"
78#include "inter22.inc"
79
80
81
82 TYPE (UNIT_TYPE_),INTENT(IN) ::UNITAB
83 my_real,
DIMENSION(NPROPM) ,
INTENT(INOUT) :: pm
84 my_real,
DIMENSION(100) ,
INTENT(INOUT) :: parmat
85 my_real,
DIMENSION(MAXUPARAM) ,
INTENT(INOUT) :: uparam
86 INTEGER, DIMENSION(MAXFUNC) ,INTENT(INOUT) :: IFUNC
87 INTEGER, INTENT(INOUT) :: ISRATE,IMATVIS,NFUNC,MAXFUNC,MAXUPARAM,NUPARAM,NUVAR
88 TYPE(MLAW_TAG_),INTENT(INOUT) :: MTAG
89 INTEGER,INTENT(IN) :: ID
90 CHARACTER(LEN=NCHARTITLE) ,INTENT(IN) :: TITR
91 my_real,
INTENT(INOUT) :: stifint
92 TYPE(SUBMODEL_DATA),INTENT(IN) :: LSUBMODEL(NSUBMOD)
93 INTEGER,INTENT(IN) :: MAT_ID
94 TYPE(MATPARAM_STRUCT_) ,INTENT(INOUT) :: MATPARAM
95
96
97
98 LOGICAL :: IS_AVAILABLE,IS_ENCRYPTED
99 my_real c1, c2, p0, vis, visv, rho10, rho20, a1, pmin, gam, vis2, pshift, visv2, rho0, rhor
100 INTEGER ISOLVER, IRHO
101 character*63 chain1
102
103
104
105
106 israte = 0
107 is_encrypted = .false.
108 is_available = .false.
109
110
111
112
114 IF(trimat==0)trimat = -2
115
117
118 CALL hm_get_floatv(
'MAT_RHO' ,rho0 ,is_available, lsubmodel, unitab)
119 CALL hm_get_floatv(
'Refer_Rho' ,rhor ,is_available, lsubmodel, unitab)
120 CALL hm_get_floatv(
'MAT_PSH' ,pshift ,is_available, lsubmodel, unitab)
121
122
123 CALL hm_get_floatv(
'Lqud_Rho_l' ,rho10 ,is_available, lsubmodel, unitab)
124 CALL hm_get_floatv(
'C_l' ,c1 ,is_available, lsubmodel, unitab)
125 CALL hm_get_floatv(
'ALPHA1' ,a1 ,is_available, lsubmodel, unitab)
126 CALL hm_get_floatv(
'Nu_l' ,vis ,is_available, lsubmodel, unitab)
127 CALL hm_get_floatv(
'Bulk_Ratio_l' ,visv ,is_available, lsubmodel, unitab)
128
129
130 CALL hm_get_floatv(
'Lqud_Rho_g' ,rho20 ,is_available, lsubmodel, unitab)
131 CALL hm_get_floatv(
'Lqud_Gamma_bulk',gam ,is_available, lsubmodel, unitab
133 CALL hm_get_floatv(
'Nu_g' ,vis2 ,is_available, lsubmodel, unitab)
134 CALL hm_get_floatv(
'Bulk_Ratio_g' ,visv2 ,is_available, lsubmodel, unitab)
135
136 irho=0
137 IF(gam*c1==0)irho=int(p0)
138
139
140
141
142 IF(a1<zero .OR. a1>one)THEN
143 chain1 = 'INITIAL MASSIC FRACTION MUST BE BETWEEN 0 AND 1. '
145 . msgtype=msgerror,
146 . anmode=aninfo,
147 . i1=37,
149 . c1='ERROR',
150 . c2=titr,
151 . c3=chain1)
152 a1 = one
153 ENDIF
154
155 IF(rho10<zero .OR. rho20<zero)THEN
156 chain1 = 'INITIAL DENSITY CANNOT BE NEGATIVE. '
158 . msgtype=msgerror,
159 . anmode=aninfo,
160 . i1=37,
162 . c1='ERROR',
163 . c2=titr,
164 . c3=chain1)
165 IF(rho10<zero)rho10 = em20
166 IF(rho20<zero)rho20 = em20
167 ENDIF
168
169
170
171
172
173
174
175
176 isolver = 1
177 IF(int22 > 0) isolver = 2
178
179
180 IF(pshift==zero)THEN
181 pshift=-p0
182 ELSE
183 pshift = -pshift
184 ENDIF
185
186
187
188
189
190 pmin = - p0
191 rho0 = rho10 * a1 + (one-a1)*rho20
192 rhor = rho0
193 pm(01) = rho0
194 pm(89) = rho0
195 pm(91) = rho10
196 pm(31) = p0+pshift
197 pm(104) = p0+pshift
198
199 uparam(1) = vis
200 uparam(2) = two*vis
201 uparam(3) = (visv-uparam(2))*third
202 uparam(4) = c1
203 c2 = - pmin * gam
204 uparam(5) = gam
205 uparam(6) = c1/rho10
206 uparam(7) = gam * p0 / rho20
207 uparam(8) = pmin
208 uparam(9) = p0
209 uparam(10) = a1
210 uparam(11) = rho10
211 uparam(12) = rho20
212 uparam(13) = vis2
213 uparam(14) = two*vis2
214 uparam(15) = (visv2-uparam(14))*third
215 uparam(16) = pshift
216 uparam(17) = isolver
217
218 nuparam = 17
219 ifunc(1) = irho
220 nfunc = 1
221 nuvar = 5
222 stifint = c1
223
224
226
227
229
230
231
232
233 IF(gam*c1>0.)THEN
234 WRITE(iout,1011) trim(titr),mat_id,37
235 WRITE(iout,1000)
236
237 IF(isolver==2)WRITE(iout,1102)
238 IF(is_encrypted)THEN
239 WRITE(iout,'(5X,A,//)')'CONFIDENTIAL DATA'
240 ELSE
241 WRITE(iout,1012) rho0
242 WRITE(iout,1100)rho10,c1,a1,vis,visv,rho20,gam,p0,vis2,visv2,p0,-pshift
243 ENDIF
244 ELSE
245 WRITE(iout,1001)
246 IF(is_encrypted)THEN
247 WRITE(iout,'(5X,A,//)')'CONFIDENTIAL DATA'
248 ELSE
249 WRITE(iout,1101)rho10,rho20,irho
250 END IF
251 END IF
252
253
254 RETURN
255 1000 FORMAT(
256 & 5x,' BI-PHASE LIQUID-GAS LAW',/,
257 & 5x,' -----------------------',//)
258 1011 FORMAT(/
259 & 5x,a,/,
260 & 5x,'MATERIAL NUMBER . . . . . . . . . . . .=',i10/,
261 & 5x,'MATERIAL LAW. . . . . . . . . . . . . .=',i10/)
262 1012 FORMAT(
263 & 5x,'INITIAL DENSITY . . . . . . . . . . . .=',1pg20
264 1001 FORMAT(
265 & 5x,' BOUNDARY FOR BI-PHASE LIQUID-GAS LAW',/,
266 & 5x,' -------------------------------------',//)
267 1100 FORMAT(
268 & 5x,'LIQUID REFERENCE DENSITY. . . . . . . .=',1pg20.13/
269 & 5x,'LIQUID BULK STIFFNESS . . . . . . . . .=',1pg20.13/
270 & 5x,'LIQUID MASS RATIO . . . . . . . . . . .=',1pg20.13/
271 & 5x,' (0. = FULL OF GAZ; 1. = FULL OF LIQUID)'/
272 & 5x,'LIQUID SHEAR KINEMATIC VISCOSITY . . . =',1pg20.13/
273 & 5x,'LIQUID BULK KINEMATIC VISCOSITY . . . =',1pg20.13//
274 & 5x,'GAS REFERENCE DENSITY. . . . . . . .=',1pg20.13/
275 & 5x,'GAS CONSTANT. . . . . . . . . . . . . .=',1pg20.13/
276 & 5x,'GAS REFERENCE PRESSURE. . . . . . . . .=',1pg20.13/
277 & 5x,'GAS SHEAR KINEMATIC VISCOSITY . . . . .=',1pg20.13/
278 & 5x,'GAS BULK KINEMATIC VISCOSITY . . . . .=',1pg20.13//
279 & 5x,'REFERENCE PRESSURE. . . . . . . . . . .=',1pg20.13/
280 & 5x,'PRESSURE SHIFT. . . . . . . . . . . . .=',1pg20.13//)
281 1101 FORMAT(
282 & 5x,'LIQUID REFERENCE DENSITY. . . . . . . .=',1pg20.13/
283 & 5x,'GAS REFERENCE DENSITY. . . . . . . .=',1pg20.13/
284 & 5x,'SCALING FUNCTION FOR GAS DENSITY. . . .=',i10)
285 1102 FORMAT(
286 & 5x,'ITERATIVE NEWTON SOLVER '
subroutine hm_get_floatv(name, rval, is_available, lsubmodel, unitab)
subroutine hm_option_is_encrypted(is_encrypted)
subroutine init_mat_keyword(matparam, keyword)
type(alemuscl_param_) alemuscl_param
integer, parameter nchartitle
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)