40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
58 USE elbufdef_mod
60
61
62
63#include "implicit_f.inc"
64
65
66
67#include "com01_c.inc"
68#include "com04_c.inc"
69#include "param_c.inc"
70#include "mvsiz_p.inc"
71
72
73
74 INTEGER, INTENT(IN) :: ITYP
75 INTEGER, INTENT(IN) :: IX(NIX,*),IPARG(NPARG,NGROUP),NIX,NG
76 my_real :: wa_l(*),x(3,numnod),evar(mvsiz),vol(mvsiz)
77 TYPE (ELBUF_STRUCT_), DIMENSION(NGROUP), TARGET :: ELBUF_TAB
78TYPE(t_ale_connectivity), INTENT(IN) :: ALE_CONNECTIVITY
79
80
81
82 INTEGER I, J
83 INTEGER :: NFT
84 INTEGER :: MLW
85 INTEGER :: NEL
86 INTEGER IV(6), IE
87 my_real :: grad(3) , xi(8), yi(8), zi(8) , sch
88 INTEGER :: NC(8)
89 my_real :: n(3,6,mvsiz), rho(6), valvois(6),
area(mvsiz), nx, a1, a2
90 INTEGER :: IAD2,
91
92
93
94
95 mlw = iparg(01,ng)
96 nel = iparg(02,ng)
97 nft = iparg(03,ng)
98
99
100
101
102 IF(n2d == 0 .AND. ityp == 1)THEN
103 DO i=1,nel
104 ie=i+nft
105
106 nc(1)=ix(2,ie)
107 nc(2)=ix(3,ie)
108 nc(3)=ix(4,ie)
109 nc(4)=ix(5,ie)
110 nc(5)=ix(6,ie)
111 nc(6)=ix(7,ie)
112 nc(7)=ix(8,ie)
113 nc(8)=ix(9,ie)
114
115 xi(1)=x(1,nc(1))
116 yi(1)=x(2,nc(1))
117 zi(1)=x(3,nc(1))
118
119 xi(2)=x(1,nc(2))
120 yi(2)=x(2,nc(2))
121 zi(2)=x(3,nc(2))
122
123 xi(3)=x(1,nc(3))
124 yi(3)=x(2,nc(3))
125 zi(3)=x(3,nc(3))
126
127 xi(4)=x(1,nc(4))
128 yi(4)=x(2,nc(4))
129 zi(4)=x(3,nc(4))
130
131 xi(5)=x(1,nc(5))
132 yi(5)=x(2,nc(5))
133 zi(5)=x(3,nc(5))
134
135 xi(6)=x(1,nc(6))
136 yi(6)=x(2,nc(6))
137 zi(6)=x(3,nc(6))
138
139 xi(7)=x(1,nc(7))
140 yi(7)=x(2,nc(7))
141 zi(7)=x(3,nc(7))
142
143 xi(8)=x(1,nc(8))
144 yi(8)=x(2,nc(8))
145 zi(8)=x(3,nc(8))
146
147 n(1,1,i)=(yi(3)-yi(1))*(zi(2)-zi(4)) - (zi(3)-zi(1))*(yi(2)-yi(4))
148 n(2,1,i)=(zi(3)-zi(1))*(xi(2)-xi(4)) - (xi(3)-xi(1))*(zi(2)-zi(4))
149 n(3,1,i)=(xi(3)-xi(1))*(yi(2)-yi(4)) - (yi(3)-yi(1))*(xi(2)-xi(4))
150
151 n(1,2,i)=(yi(7)-yi(4))*(zi(3)-zi(8)) - (zi(7)-zi(4))*(yi(3)-yi(8))
152 n(2,2,i)=(zi(7)-zi(4))*(xi(3)-xi(8)) - (xi(7)-xi(4))*(zi(3)-zi(8))
153 n(3,2,i)=(xi(7)-xi(4))*(yi(3)-yi(8)) - (yi(7)-yi(4))*(xi(3)-xi(8))
154
155 n(1,3,i)=(yi(6)-yi(8))*(zi(7)-zi(5)) - (zi(6)-zi(8))*(yi(7)-yi(5))
156 n(2,3,i)=(zi(6)-zi(8))*(xi(7)-xi(5)) - (xi(6)-xi(8))*(zi(7)-zi(5))
157 n(3,3,i)=(xi(6)-xi(8))*(yi(7)-yi(5)) - (yi(6)-yi(8))*(xi(7)-xi(5))
158
159 n(1,4,i)=(yi(2)-yi(5))*(zi(6)-zi(1)) - (zi(2)-zi(5))*(yi(6)-yi(1))
160 n(2,4,i)=(zi(2)-zi(5))*(xi(6)-xi(1)) - (xi(2)-xi(5))*(zi(6)-zi(1))
161 n(3,4,i)=(xi(2)-xi(5))*(yi(6)-yi(1)) - (yi(2)-yi(5))*(xi(6)-xi(1))
162
163 n(1,5,i)=(yi(7)-yi(2))*(zi(6)-zi(3)) - (zi(7)-zi(2))*(yi(6)-yi(3))
164 n(2,5,i)=(zi(7)-zi(2))*(xi(6)-xi(3)) - (xi(7)-xi(2))*(zi(6)-zi(3))
165 n(3,5,i)=(xi(7)-xi(2))*(yi(6)-yi(3)) - (yi(7)-yi(2))*(xi(6)-xi(3))
166
167 n(1,6,i)=(yi(8)-yi(1))*(zi(4)-zi(5)) - (zi(8)-zi(1))*(yi(4)-yi(5))
168 n(2,6,i)=(zi(8)-zi(1))*(xi(4)-xi(5)) - (xi(8)-xi(1))*(zi(4)-zi(5))
169 n(3,6,i)=(xi(8)-xi(1))*(yi(4)-yi(5)) - (yi(8)-yi(1))*(xi(4)-xi(5))
170
171 n(1:3,1,i) = half * n(1:3,1,i)
172 n(1:3,2,i) = half * n(1:3,2,i)
173 n(1:3,3,i) = half * n(1:3,3,i)
174 n(1:3,4,i) = half * n(1:3,4,i)
175 n(1:3,5,i) = half * n(1:3,5,i)
176 n(1:3,6,i) = half * n(1:3,6,i)
177 ENDDO
178 ELSEIF(n2d > 0 .AND. ityp == 2)THEN
179 DO i=1,nel
180 ie=i+nft
181
182 nc(1)=ix(2,ie)
183 nc(2)=ix(3,ie)
184 nc(3)=ix(4,ie)
185 nc(4)=ix(5,ie)
186
187
188 xi(1)=zero
189 yi(1)=x(2,nc(1))
190 zi(1)=x(3,nc(1))
191
192 xi(2)=zero
193 yi(2)=x(2,nc(2))
194 zi(2)=x(3,nc(2))
195
196 xi(3)=zero
197 yi(3)=x(2,nc(3))
198 zi(3)=x(3,nc(3))
199
200 xi(4)=zero
201 yi(4)=x(2,nc(4))
202 zi(4)=x(3,nc(4))
203
204
205 n(1,1,i) = zero
206 n(2,1,i) = (zi(2)-zi(1))
207 n(3,1,i) =-(yi(2)-yi(1))
208
209 n(1,2,i) = zero
210 n(2,2,i) = (zi(3)-zi(2))
211 n(3,2,i) =-(yi(3)-yi(2))
212
213 n(1,3,i) = zero
214 n(2,3,i) = (zi(4)-zi(3))
215 n(3,3,i) =-(yi(4)-yi(3))
216
217 n(1,4,i) = zero
218 n(2,4,i) = (zi(1)-zi(4))
219 n(3,4,i) =-(yi(1)-yi(4))
220
221 n(1:3,5:6,i)=zero
222
223 IF(mlw /= 151)THEN
224 a1 =yi(2)*(zi(3)-zi(4))+yi(3)*(zi(4)-zi(2))+yi(4)*(zi(2)-zi(3))
225 a2 =yi(2)*(zi(4)-zi(1))+yi(4)*(zi(1)-zi(2))+yi(1)*(zi(2)-zi(4))
226 area(i)=(a1+a2)* half
227 ELSE IF (elbuf_tab(ng)%GBUF%G_AREA >= i) THEN
228 area(i)=elbuf_tab(ng)%GBUF%AREA(i)
229 ELSE
231 ENDIF
232 ENDDO
233 ELSEIF(n2d > 0 .AND. ityp == 7)THEN
234 DO i=1,nel
235 ie=i+nft
236
237 nc(1)=ix(2,ie)
238 nc(2)=ix(3,ie)
239 nc(3)=ix(4,ie)
240
241 xi(1)=zero
242 yi(1)=x(2,nc(1))
243 zi(1)=x(3,nc(1))
244
245 xi(2)=zero
246 yi(2)=x(2,nc(2))
247 zi(2)=x(3,nc(2))
248
249 xi(3)=zero
250 yi(3)=x(2,nc(3))
251 zi(3)=x(3,nc(3))
252
253
254 n(1,1,i) = zero
255 n(2,1,i) = (zi(2)-zi(1))
256 n(3,1,i) =-(yi(2)-yi(1))
257
258 n(1,2,i) = zero
259 n(2,2,i) = (zi(3)-zi(2))
260 n(3,2,i) =-(yi(3)-yi(2))
261
262 n(1,3,i) = zero
263 n(2,3,i) = (zi(1)-zi(3))
264 n(3,3,i) =-(yi(1)-yi(3))
265
266 n(1:3,4:6,i)=zero
267
268 IF(mlw /= 151)THEN
269 nx = half * ((yi(2) - yi(1)) * (zi(3) - zi(1)) - (zi(2) - zi(1)) * (yi(3) - yi(1)))
271 ELSE
272 area(i)=elbuf_tab(ng)%GBUF%AREA(i)
273 ENDIF
274 ENDDO
275 ELSE
276 n(1:3,1:6,1:nel) = zero
278 ENDIF
279
280
281
282
283 IF(n2d == 0 .AND. ityp==1)THEN
284 DO i=1,nel
285 ie=i+nft
286 iad2 = ale_connectivity%ee_connect%iad_connect(ie)
287 lgth = ale_connectivity%ee_connect%iad_connect(ie+1)-ale_connectivity%ee_connect%iad_connect(ie)
288 DO j=1,lgth
289 iv(j) = ale_connectivity%ee_connect%connected(iad2 + j - 1)
290 IF(iv(j) > 0)THEN
291 valvois(j
292 ELSE
293 valvois(j) = wa_l(ie)
294 ENDIF
295 ENDDO
296 rho(1) = half*( wa_l(ie) + valvois(1) )
297 rho(2) = half*( wa_l(ie) + valvois(2) )
298 rho(3) = half*( wa_l(ie) + valvois(3) )
299 rho(4) = half*( wa_l(ie) + valvois(4) )
300 rho(5) = half*( wa_l(ie) + valvois(5) )
301 rho(6) = half*( wa_l(ie) + valvois(6) )
302
303
304
305 grad(1:3) = zero
306 grad(1:3) = grad(1:3) + rho(1)*n(1:3,1,i)
307 grad(1:3) = grad(1:3) + rho(2)*n(1:3,2,i)
308 grad(1:3) = grad(1:3) + rho(3)*n(1:3,3,i)
309 grad(1:3) = grad(1:3) + rho(4)*n(1:3,4,i)
310 grad(1:3) = grad(1:3) + rho(5)*n(1:3,5,i)
311 grad(1:3) = grad(1:3) + rho(6)*n(1:3,6,i)
312 grad(1:3) = grad(1:3) / vol(i)
313
314
315
316 sch = sqrt(sum(grad(1:3)*grad(1:3)))
317 sch = exp(-sch)
318 evar(i) = sch
319 enddo
320 ELSEIF(n2d > 0 .AND. ityp == 2)THEN
321 DO i=1,nel
322 ie=i+nft
323 iad2 = ale_connectivity%ee_connect%iad_connect(ie)
324 lgth = ale_connectivity%ee_connect%iad_connect(ie+1)-iad2
325 DO j=1,lgth
326 iv(j)=ale_connectivity%ee_connect%connected(iad2 +
327 IF(iv(j) > 0)THEN
328 valvois(j) = wa_l(iv(j))
329 ELSE
330 valvois(j) = wa_l(ie)
331 ENDIF
332 ENDDO
333 rho(1) = half*( wa_l(ie) + valvois(1) )
334 rho(2) = half*( wa_l(ie) + valvois(2) )
335 rho(3) = half*( wa_l(ie) + valvois(3) )
336 rho(4) = half*( wa_l(ie) + valvois(4) )
337
338
339
340 grad(1:3) = zero
341 grad(2:3) = grad(2:3) + rho(1)*n(2:3,1,i)
342 grad(2:3) = grad(2:3) + rho(2)*n(2:3,2,i)
343 grad(2:3) = grad(2:3) + rho(3)*n(2:3,3,i)
344 grad(2:3) = grad(2:3) + rho(4)*n(2:3,4,i)
345 grad(2:3) = grad(2:3) /
area(i)
346
347
348
349 sch = sqrt(sum(grad(2:3)*grad(2:3)))
350 sch = exp(-sch)
351 evar(i) = sch
352 enddo
353 ELSEIF(n2d > 0 .AND. ityp == 7)THEN
354 DO i=1,nel
355 ie=i+nft
356 iad2 = ale_connectivity%ee_connect%iad_connect(ie)
357 lgth = ale_connectivity%ee_connect%iad_connect(ie+1)-ale_connectivity%ee_connect%iad_connect
358 DO j=1,3
359 iv(j)=ale_connectivity%ee_connect%connected(iad2 + j - 1)
360 IF(iv(j) > 0)THEN
361 valvois(j) = wa_l(iv(j))
362 ELSE
363 valvois(j) = wa_l(ie)
364 ENDIF
365 ENDDO
366 rho(1) = half*( wa_l(ie) + valvois(1) )
367 rho(2) = half*( wa_l(ie) + valvois(2) )
368 rho(3) = half*( wa_l(ie) + valvois(3) )
369
370
371
372 grad(1:3) = zero
373 grad(2:3) = grad(2:3) + rho(1)*n(2:3,1,i)
374 grad(2:3) = grad(2:3) + rho(2)*n(2:3,2,i)
375 grad(2:3) = grad(2:3) + rho(3)*n(2:3,3,i)
376 grad(2:3) = grad(2:3) /
area(i)
377
378
379
380 sch = sqrt(sum(grad(2:3)*grad(2:3)))
381 sch = exp(-sch)
382 evar(i) = sch
383 enddo
384 ELSE
385 evar(1:nel)=one
386 ENDIF
387
subroutine area(d1, x, x2, y, y2, eint, stif0)