33
35
36
37
38#include "implicit_f.inc"
39
40
41
42#include "mvsiz_p.inc"
43#include "comlock.inc"
44
45
46
47#include "com_xfem1.inc"
48
49
50
51 INTEGER JFT,JLT,NFT,ILEV,IXC(NIXC,*),ELCUTC(2,*),(*),
52 .
53TYPE (XFEM_EDGE_) , DIMENSION(*) :: CRKEDGE
54
55
56
57 INTEGER I,K,K1,K2,K3,K4,KK,p1,p2,NOD1,NOD2,IED0,IED1,IED2,
58 . IEDGE,IEDGE1,IEDGE2,EDGE,EDGE1,EDGE2,IXEL,ILAY,ITRI,ILV,
59 . FAC,OK,SN,ELCUT,ELCRK,IADC1,IADC2,IADC3,IADC4,NX1,NX2,NX3,NX4
60 INTEGER IFI0(4,MVSIZ),POS(2),
61 . IED(4),IXI(4),D(4),DX(8),NX(MVSIZ)
63 . xint,yint,zint,xm,ym,zm,x10,y10,z10,x20,y20,z20,beta
65 . xin(4,mvsiz),yin(4,mvsiz),zin(4,mvsiz),
66 . xx(4,mvsiz),yy(4,mvsiz),zz(4,mvsiz)
67
68 DATA d /2,3,4,1/
69 DATA dx/1,2,3,4,1,2,3,4/
70
71
72
73 ixel = mod(ilev-1, nxel) + 1
74 ilay = (ilev-ixel)/nxel + 1
75 p1 = 0
76 p2
77
78 DO i=jft,jlt
79 xin(1,i) = zero
80 yin(1,i) = zero
81 zin(1,i) = zero
82 xin(2,i) = zero
83 yin(2,i) = zero
84 zin(2,i) = zero
85 xin(3,i) = zero
86 yin(3,i) = zero
87 zin(3,i) = zero
88 xin(4,i) = zero
89 yin(4,i) = zero
90 zin(4,i) = zero
91 END DO
92
93 DO i=jft,jlt
94 elcrk = iel_crk(i+nft)
95 iadc1 = iadc_crk(1,elcrk)
96 iadc2 = iadc_crk(2,elcrk)
97 iadc3 = iadc_crk(3,elcrk)
98 iadc4 = iadc_crk(4,elcrk)
99
101
104
105 ifi0(1,i) = isign(1,ifi0(1,i))
106 ifi0(2,i) = isign(1,ifi0(2,i))
107 ifi0(3,i) = isign(1,ifi0(3,i))
108
109
110
111
112
113 xx(1,i) =
crkavx(ilev)%X(1,iadc1)
114 yy(1,i) =
crkavx(ilev)%X(2,iadc1)
115 zz(1,i) =
crkavx(ilev)%X(3,iadc1)
116
117 xx(2,i) =
crkavx(ilev)%X(1,iadc2)
118 yy(2,i) =
crkavx(ilev)%X(2,iadc2)
119 zz(2,i) =
crkavx(ilev)%X(3,iadc2)
120
121 xx(3,i) =
crkavx(ilev)%X(1,iadc3)
123 zz(3,i) =
crkavx(ilev)%X(3,iadc3)
124
125 xx(4,i) =
crkavx(ilev)%X(1,iadc4)
126 yy(4,i) =
crkavx(ilev)%X(2,iadc4)
127 zz(4,i) =
crkavx(ilev)%X(3,iadc4)
128 END DO
129
130
131
132 DO i=jft,jlt
133 elcrk = iel_crk(i+nft)
135 IF (elcut /= 0) THEN
136 DO k=1,4
137 ied0 = crkedge(ilay)%IEDGEC(k,elcrk)
138 IF (ied0 > 0) THEN
139 edge = xedge4n(k,elcrk)
140 beta = crkedge(ilay)%RATIO(edge)
141 nod1 = nodedge(1,edge)
142 nod2 = nodedge(2,edge)
143 IF (nod1 == ixc(k+1,i+nft) .and.
144 . nod2 == ixc(d(k)+1,i+nft)) THEN
145 p1 = k
146 p2 = d(k)
147 ELSEIF (nod2 == ixc(k+1,i+nft) .and.
148 . nod1 == ixc(d(k)+1,i+nft)) THEN
149 p1 = d(k)
150 p2 = k
151 ENDIF
152 x10 = xx(p1,i)
153 y10 = yy(p1,i)
154 z10 = zz(p1,i)
155 x20 = xx(p2,i)
156 y20 = yy(p2,i)
157 z20 = zz(p2,i)
158 xin(ied0,i) = x10 + beta*(x20-x10)
159 yin(ied0,i) = y10 + beta*(y20-y10)
160 zin(ied0,i) = z10 + beta*(z20-z10)
161 END IF
162 END DO
163 END IF
164 END DO
165
166
167
168
169
170
171 DO i=jft,jlt
172 IF (elcutc(1,i+nft) == 0) cycle
173
174 elcrk = iel_crk(i+nft)
176 IF (elcut == 0) cycle
179 nx2 = dx(nx1+1)
180 nx3 = dx(nx1+2)
181 nx4 = dx(nx1+3)
182
183 IF (itri == 0) THEN
184
185 pos(1) = 0
186 pos(2) = 0
187 fac = 0
188
189 IF (ixel == 1) THEN
190
191 DO k=1,4
192 ied0 = crkedge(ilay)%IEDGEC(k,elcrk)
193 IF (ied0 > 0) THEN
194 fac = fac + 1
195 p1 = k
196 p2 = d(p1)
197 IF (ifi0(p1,i) < 0) pos(fac
198 IF (ifi0(p2,i) < 0) pos(fac) = p2
199 ied(fac) = p1
200 IF (fac == 2) EXIT
201 ENDIF
202 END DO
203
204 IF (pos(1) /= 0 .and. pos(2) /= 0) THEN
205 DO k=1,2
206 iedge = crkedge(ilay)%IEDGEC(ied(k),elcrk)
207 IF (iedge > 0) THEN
208
209 kk =
crkshell(ilev)%XNODEL(pos(k),elcrk)
210 kk = kk -
crknod(ilev)%CRKNUMNODS * (ilev-1)
211 crkavx(ilev)%XX(1,kk) = xin(iedge,i)
212
213 crkavx(ilev)%XX(3,kk) = zin(iedge,i)
214 ENDIF
215 END DO
216 ENDIF
217
218 ELSE IF (ixel == 2) THEN
219
220 DO k=1,4
221 ied0 = crkedge(ilay)%IEDGEC(k,elcrk)
222 IF (ied0 > 0) THEN
223 fac = fac + 1
224 p1 = k
225 p2 = d(p1)
226 IF (ifi0(p1,i) > 0) pos(fac) = p1
227 IF (ifi0(p2,i) > 0) pos(fac) = p2
228 ied(fac) = p1
229 IF (fac == 2) EXIT
230 ENDIF
231 END DO
232
233 IF (pos(1) /= 0 .and. pos(2) /= 0) THEN
234 DO k=1,2
235 iedge = crkedge(ilay)%IEDGEC(ied(k),elcrk)
236 IF (iedge > 0) THEN
237
238 kk =
crkshell(ilev)%XNODEL(pos(k),elcrk)
239 kk = kk -
crknod(ilev)%CRKNUMNODS * (ilev-1)
240 crkavx(ilev)%XX(1,kk) = xin(iedge,i)
241 crkavx(ilev)%XX(2,kk) = yin(iedge,i)
242 crkavx(ilev)%XX(3,kk) = zin(iedge,i)
243 ENDIF
244 END DO
245 ENDIF
246
247 ELSE IF (ixel == 3) THEN
248 ENDIF
249
250
251 ELSEIF (itri < 0) THEN
252
253 ied1 = nx1
254 ied2 = nx4
255 iadc1 = iadc_crk(nx1,elcrk)
256 iadc2 = iadc_crk(nx2,elcrk)
257 iadc3 = iadc_crk(nx3,elcrk)
258 iadc4 = iadc_crk(nx4,elcrk)
259
260 iedge1 = crkedge(ilay)%IEDGEC(ied1,elcrk)
261 iedge2 = crkedge(ilay)%IEDGEC(ied2,elcrk)
262 edge1 = xedge4n(ied1,elcrk)
263 edge2 = xedge4n(ied2,elcrk)
264
265 kk =
crkshell(ilev)%XNODEL(nx1,elcrk)
266 k1 = kk -
crknod(ilev)%CRKNUMNODS * (ilev-1)
267 kk =
crkshell(ilev)%XNODEL(nx2,elcrk)
268 k2 = kk -
crknod(ilev)%CRKNUMNODS * (ilev-1)
269 kk =
crkshell(ilev)%XNODEL(nx3,elcrk)
270 k3 = kk -
crknod(ilev)%CRKNUMNODS * (ilev-1)
271 kk =
crkshell(ilev)%XNODEL(nx4,elcrk)
272 k4 = kk -
crknod(ilev)%CRKNUMNODS * (ilev-1)
273
274
275
276
277
278
279 IF (ixel == 1) THEN
280
281
282
283
284 crkavx(ilev)%XX(1,k2) = xin(iedge1,i)
285 crkavx(ilev)%XX(2,k2) = yin(iedge1,i)
286 crkavx(ilev)%XX(3,k2) = zin(iedge1,i)
287 crkavx(ilev)%XX(1,k3) = xin(iedge2,i)
288 crkavx(ilev)%XX(2,k3) = yin(iedge2,i)
289 crkavx(ilev)%XX(3,k3) = zin(iedge2,i)
293
294 ELSE IF (ixel == 2) THEN
295
296
297
298
299 crkavx(ilev)%XX(1,k1) = xin(iedge2,i)
300 crkavx(ilev)%XX(2,k1) = yin(iedge2,i)
301 crkavx(ilev)%XX(3,k1) = zin(iedge2,i)
305
306 ELSE IF (ixel == 3) THEN
307
308
309
310
311
312 crkavx(ilev)%XX(1,k1) = xin(iedge1,i)
313 crkavx(ilev)%XX(2,k1) = yin
314 crkavx(ilev)%XX(3,k1) = zin(iedge1,i)
315
316 crkavx(ilev)%XX(1,k4) = xin(iedge2,i)
317 crkavx(ilev)%XX(2,k4) = yin(iedge2,i)
318 crkavx(ilev)%XX(3,k4) = zin(iedge2,i)
319
320
321
322
323
324
325
326
327
328
329 END IF
330
331
332 ELSEIF (itri > 0) THEN
333
334 ied1 = nx1
335 ied2 = nx4
336 iadc1 = iadc_crk(nx1,elcrk)
337 iadc2 = iadc_crk(nx2,elcrk)
338 iadc3 = iadc_crk(nx3,elcrk)
339 iadc4 = iadc_crk(nx4,elcrk)
340
341 iedge1 = crkedge(ilay)%IEDGEC(ied1,elcrk)
342 iedge2 = crkedge(ilay)%IEDGEC(ied2,elcrk)
343 edge1 = xedge4n(ied1,elcrk)
344 edge2 = xedge4n(ied2,elcrk)
345
346 kk =
crkshell(ilev)%XNODEL(nx1,elcrk)
347 k1 = kk -
crknod(ilev)%CRKNUMNODS * (ilev-1)
348 kk =
crkshell(ilev)%XNODEL(nx2,elcrk)
349 k2 = kk -
crknod(ilev)%CRKNUMNODS * (ilev-1)
350 kk =
crkshell(ilev)%XNODEL(nx3,elcrk)
351 k3 = kk -
crknod(ilev)%CRKNUMNODS * (ilev-1)
352 kk =
crkshell(ilev)%XNODEL(nx4,elcrk)
353 k4 = kk -
crknod(ilev)%CRKNUMNODS * (ilev-1)
354
355
356
357
358
359
360
361 IF (ixel == 1) THEN
362
363
364
365
366 crkavx(ilev)%XX(1,k1) = xin(iedge1,i)
367 crkavx(ilev)%XX(2,k1) = yin(iedge1,i)
368 crkavx(ilev)%XX(3,k1) = zin(iedge1,i)
372
373
374 ELSE IF (ixel == 2) THEN
375
376
377
378
379 crkavx(ilev)%XX(1,k2) = xin(iedge1,i)
380 crkavx(ilev)%XX(2,k2) = yin(iedge1,i)
381 crkavx(ilev)%XX(3,k2) = zin(iedge1,i)
385 crkavx(ilev)%XX(1,k4) = xin(iedge2,i)
386 crkavx(ilev)%XX(2,k4) = yin(iedge2,i)
387 crkavx(ilev)%XX(3,k4) = zin(iedge2,i)
388
389
390 ELSE IF (ixel == 3) THEN
391
392
393
394
395
396 crkavx(ilev)%XX(1,k1) = xin(iedge2,i)
397 crkavx(ilev)%XX(2,k1) = yin(iedge2,i)
398 crkavx(ilev)%XX(3,k1) = zin(iedge2,i)
399
400
401
402
403 crkavx(ilev)%XX(1,k2) = xin(iedge1,i)
404 crkavx(ilev)%XX(2,k2) = yin(iedge1,i)
405 crkavx(ilev)%XX(3,k2) = zin(iedge1,i)
406
407
408
409
410
411
412
413
414
415
416 END IF
417
418 END IF
419
420 ENDDO
421
422 RETURN
type(xfem_phantom_), dimension(:), allocatable xfem_phantom
type(xfem_nodes_), dimension(:), allocatable crknod
type(xfem_avx_), dimension(:), allocatable crkavx
type(xfem_shell_), dimension(:), allocatable crkshell