45
46
47
48 USE elbufdef_mod
50
51
52
53#include "implicit_f.inc"
54
55
56
57#include "param_c.inc"
58#include "com01_c.inc"
59#include "com08_c.inc"
60#include "spmd_c.inc"
61
62
63
64 INTEGER NVENT, NELT, NEL, ITTF, ELEM(3,*), IBAGHOL(NIBHOL,*), ITAGEL(*),
65 . ELTG(*), IPARG(NPARG,*), IPM(NPROPMI,*),
66 . MATTG(*), IGROUPC(*), IGROUPTG(*)
67
68
70 . aoutot,
71 . elsout(*), elarea(*), elsini(*), rvolu(*), svent(nvent),
72 . rbaghol(nrbhol,*),
poro(*), p(*), pm(npropm,*), porosity(*)
73 TYPE(ELBUF_STRUCT_), DIMENSION(NGROUP) :: ELBUF_TAB
74
75
76
77
78 INTEGER IEL, K, N1, N2, N3,
79 . IDEF, IVENT, IVENTYP,
80 . IPORT, IPORP, , IPORT1, IPORP1, IPORA1,
81 . NG, IM, NFUNC, MTN,
82 . ILEAKAGE,IBLOCKAGE
83 INTEGER JEL, NFT, NELG
85 . pext, avent, bvent, aout, aout1,
86 . fport, fporp, fpora, fport1, fporp1, fpora1, deri,
87 . scalt, scalp, scals,
area, aini, exten,
88 . ttf, flc, fac, fac1, svtfac, pmean, tt1,
89 . tvent
91 EXTERNAL get_u_func
92
93 pext =rvolu(3)
94 scalt=rvolu(26)
95 scalp=rvolu(27)
96 scals=rvolu(28)
97 ttf =rvolu(60)
98
99 DO ivent=1,nvent
100 svent(ivent)=zero
101 ENDDO
102
103 DO iel=1,nelt
104 IF (itagel(iel)<0) THEN
105 ivent=-itagel(iel)
106 svent(ivent)=svent(ivent)+elarea(iel)
107 ENDIF
108 ENDDO
109
110 DO ivent=1,nvent
111 rbaghol(16,ivent)=zero
112 rbaghol(17,ivent)=zero
113 rbaghol(18,ivent)=zero
114 rbaghol(19,ivent)=zero
115 rbaghol(20,ivent)=zero
116 rbaghol(21,ivent)=zero
117 rbaghol(22,ivent)=zero
118 ENDDO
119
120
121
122 aoutot=zero
123 DO iel=1,nelt
124 elsout(iel)=zero
125 IF (itagel(iel)<0) THEN
127 ivent=-itagel(iel)
128 iventyp=ibaghol(13,ivent)
129 n1=elem(1,iel)
130 n2=elem(2,iel)
131 n3=elem(3,iel)
132 pmean=third*(p(n1)+p(n2)+p(n3))
133
134
135
136 IF(iventyp==0)THEN
137
138 aini=elsini(iel)
140
141 idef=ibaghol(1,ivent)
142 IF (idef==0.OR.idef==2) cycle
143 avent =rbaghol(2,ivent)
144 bvent =rbaghol(6,ivent)
145 tvent =rbaghol(3,ivent)
146
147 iport =ibaghol(3,ivent)
148 iporp =ibaghol(4,ivent)
149 ipora =ibaghol(5,ivent)
150 iport1=ibaghol(6,ivent)
151 iporp1=ibaghol(7,ivent)
152 ipora1=ibaghol(8,ivent)
153
154 fport =rbaghol(7,ivent)
155 fporp =rbaghol(8,ivent)
156 fpora =rbaghol(9,ivent)
157 fport1=rbaghol(10,ivent)
158 fporp1=rbaghol(11,ivent)
159 fpora1=rbaghol(12,ivent)
160
163 tt1=tt-ttf
164 IF (ittf==13) tt1=tt-ttf-tvent
165 IF (ipora/=0) aout=aout*fpora*get_u_func(ipora,exten,deri)
166 IF (iport/=0)aout=aout*fport*get_u_func(iport,tt1*scalt,deri)
167 IF (iporp/=0)aout=aout*fporp*get_u_func(iporp,(pmean-pext)*scalp,deri)
168
169 IF (ipora1/=0) aout1=aout1*fpora1*get_u_func(ipora1,exten,deri)
170 IF (iport1/=0)aout1=aout1*fport1*get_u_func(iport1,tt1*scalt,deri)
171 IF (iporp1/=0)aout1=aout1*fporp1*get_u_func(iporp1,(pmean-pext)*scalp,deri)
172
173
174
175
176 ELSE
177 iblockage=ibaghol(14,ivent)
178 tt1=tt-ttf
179 IF (ittf==13) THEN
180 tvent=rbaghol(3,ivent)
181 tt1=tt-ttf-tvent
182 ENDIF
183 svtfac=zero
184
185 im = mattg(iel)
186 mtn = ipm(2,im)
187 IF (mtn/=19.AND.mtn/=58) cycle
188
189 ileakage = ipm(4,im)
190 nfunc = ipm(10,im)+ipm(6,im)
191 IF(ileakage==0) THEN
192 svtfac=zero
193 ELSEIF(ileakage==1) THEN
194 flc=pm(164,im)
195 fac=pm(165,im)
196 svtfac=flc*fac
197 ELSEIF(ileakage==2.OR.ileakage==3) THEN
198 flc=zero
199 iport=ipm(10+nfunc-1,im)
200 IF(iport > 0) THEN
201 scalt=pm(160,im)
202 fport=pm(164,im)
203 flc=fport*get_u_func(iport,tt1*scalt,deri)
204 ENDIF
205 fac=zero
206 iporp=ipm(10+nfunc-2,im)
207 IF(iporp > 0) THEN
208 scalp=pm(161,im)
209 fporp=pm(165,im)
210 IF(ileakage==2) THEN
211 fac=fporp*get_u_func(iporp,pmean*scalp,deri)
212 ELSE
213 fac=fporp*get_u_func(iporp,(pmean-pext)*scalp,deri)
214 ENDIF
215 ENDIF
216 svtfac=flc*fac
217 ELSEIF(ileakage==4) THEN
218 aini=elsini(iel)
220 ELSEIF(ileakage==5) THEN
221 IF(nspmd > 1) THEN
222 CALL ancmsg(msgid=258,anmode=aninfo,i1=ipm(1,im))
224 ELSE
225 k = eltg(iel)
226 IF(k <= numelcg)THEN
227 ng=igroupc(k)
228 ELSE
229 k=k-numelcg
230 ng=igrouptg(k)
231 ENDIF
232 nelg = iparg(2,ng)
233 nft = iparg(3,ng)
234 jel = k-nft
235 CALL porform5(svtfac,im,ipm,pm,elbuf_tab(ng),p,pext,jel,nelg)
237 ENDIF
238 ELSEIF(ileakage==6) THEN
239 aini=elsini(iel)
241 ENDIF
242
243 IF(intbag==0) THEN
245 aout1= zero
246 ELSE
247 IF(iblockage==1) THEN
248 aout = (one -
poro(iel))*
area*svtfac
249 aout1= zero
250 ELSE
252 fac1=pm(162,im)
253 IF(fac1 == zero) THEN
254 iport=ipm(10+nfunc,im)
255 IF(iport > 0) THEN
256 scalt=pm(160,im)
257 fport=pm(163,im)
258 fac1=fport*get_u_func(iport,tt1*scalt,deri)
259 ENDIF
260 ENDIF
262 ENDIF
263 ENDIF
264 ENDIF
265
266 elsout(iel)=aout+aout1
267 aoutot=aoutot+elsout(iel)
268 rbaghol(16,ivent)=rbaghol(16,ivent)+aout
269 rbaghol(17,ivent)=rbaghol(17,ivent)+aout1
270
271 IF(iel > nel) THEN
272 porosity(iel-nel)=
min(one,elsout(iel)/
area)
273 ENDIF
274 ENDIF
275 ENDDO
276
277 RETURN
subroutine area(d1, x, x2, y, y2, eint, stif0)
subroutine porform4(svtfac, im, ipm, pm, area, area0, p, pext)
subroutine porform5(svtfac, im, ipm, pm, elbuf_str, p, pext, iel, nel)
subroutine porform6(svtfac, im, pm, area, area0, p, pext)
subroutine poro(geo, nodpor, ms, x, v, w, af, am, skew, weight, nporgeo)
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)