FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_adapt_bc_pointer.f90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2! Copyright (c) 2019 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
6
7!C
8!C***
9!C*** hecmw_adapt_bc_pointer
10!C***
11!C
12!C creates NEW BOUNDARY POINTERs
13!C
14subroutine hecmw_adapt_bc_pointer (hecMESH)
15
16 use hecmw_util
17 implicit real*8 (a-h,o-z)
18 integer(kind=kint), dimension(: ), allocatable :: IW1, IW2, IW3
19 integer(kind=kint), dimension(:,:), allocatable :: IUPD
20 type(hecmwst_local_mesh) :: hecMESH
21
22 if (hecmesh%node_group%n_grp.ne.0) then
23 if (hecmesh%node_group%grp_index(hecmesh%node_group%n_grp).ne.0) then
24 !C
25 !C +------------+
26 !C | NODE-GROUP |
27 !C +------------+
28 !C===
29 allocate (iupd(hecmesh%n_adapt_node_cur, &
30 & hecmesh%node_group%n_grp))
31 allocate (iw1(hecmesh%n_adapt_node_cur))
32 allocate (iw2(hecmesh%node_group%n_grp))
33
34 iupd = 0
35 icou = 0
36 iw2 = 0
37
38 do ig= 1, hecmesh%node_group%n_grp
39 icou0= 0
40 iw1= 0
41 do k= hecmesh%node_group%grp_index(ig-1)+1, &
42 & hecmesh%node_group%grp_index(ig)
43 nod = hecmesh%node_group%grp_item(k)
44 icou = icou + 1
45 icou0 = icou0 + 1
46 iw1(nod )= 1
47 iupd(icou0,ig)= nod
48 enddo
49
50 do ie= 1, hecmesh%n_adapt_edge
51 nod1= hecmesh%adapt_edge_node(2*ie-1)
52 nod2= hecmesh%adapt_edge_node(2*ie )
53
54 if (iw1(nod1).eq.1 .and. iw1(nod2).eq.1 .and. &
55 & hecmesh%adapt_iemb(ie ).ne.0) then
56 nod3 = hecmesh%adapt_IWK(ie)
57 icou = icou + 1
58 icou0 = icou0 + 1
59 iupd(icou0,ig)= nod3
60 endif
61 enddo
62 iw2(ig)= icou0
63 enddo
64
65 !C
66 !C-- NEW POINTERs
67 nnn1= hecmesh%node_group%grp_index(hecmesh%node_group%n_grp)
68 deallocate (hecmesh%node_group%grp_item)
69 allocate (hecmesh%node_group%grp_item(icou))
70
71 hecmesh%node_group%grp_index= 0
72 do ig= 1, hecmesh%node_group%n_grp
73 hecmesh%node_group%grp_index(ig)= &
74 & hecmesh%node_group%grp_index(ig-1) + iw2(ig)
75 enddo
76 do ig= 1, hecmesh%node_group%n_grp
77 do k= 1, iw2(ig)
78 kk= hecmesh%node_group%grp_index(ig-1) + k
79 hecmesh%node_group%grp_item(kk)= iupd(k,ig)
80 enddo
81 enddo
82
83 nnn2= hecmesh%node_group%grp_index(hecmesh%node_group%n_grp)
84 ! write (*,'(a, i5, 2i8)') ' node group', hecMESH%my_rank, nnn1, nnn2
85
86 deallocate (iupd, iw1, iw2)
87 !C===
88 endif
89 endif
90
91 if (hecmesh%elem_group%n_grp.ne.0) then
92 if (hecmesh%elem_group%grp_index(hecmesh%elem_group%n_grp).ne.0) then
93 !C
94 !C +---------------+
95 !C | ELEMENT-GROUP |
96 !C +---------------+
97 !C===
98 allocate (iupd(hecmesh%n_elem, &
99 & hecmesh%elem_group%n_grp))
100 allocate (iw1(hecmesh%n_elem))
101 allocate (iw2(hecmesh%elem_group%n_grp))
102
103 iupd = 0
104 icou = 0
105 iw2 = 0
106
107 do ig= 1, hecmesh%elem_group%n_grp
108 icou0= 0
109 iw1= 0
110 do k= hecmesh%elem_group%grp_index(ig-1)+1, &
111 & hecmesh%elem_group%grp_index(ig)
112 icel = hecmesh%elem_group%grp_item(k)
113 icou = icou + 1
114 icou0 = icou0 + 1
115 iw1(icel )= 1
116 iupd(icou0,ig)= icel
117 enddo
118
119 do icel= 1, hecmesh%n_elem
120 if (hecmesh%adapt_type(icel).ne.0.and.iw1(icel).eq.1) then
121 is= hecmesh%adapt_children_index(icel-1) + 1
122 ie= hecmesh%adapt_children_index(icel)
123 ics= hecmesh%adapt_children_local(is)
124
125 if (hecmesh%when_i_was_refined_elem(ics).eq. &
126 & hecmesh%n_adapt) then
127 do k= is, ie
128 if (hecmesh%adapt_children_item(2*k-1).ne.0) then
129 iclocal= hecmesh%adapt_children_local(k)
130 icou = icou + 1
131 icou0 = icou0 + 1
132 iw1(iclocal )= 1
133 iupd(icou0,ig)= iclocal
134 endif
135 enddo
136 endif
137 endif
138 enddo
139 iw2(ig)= icou0
140 enddo
141
142 !C
143 !C-- NEW POINTERs
144 nnn1= hecmesh%elem_group%grp_index(hecmesh%elem_group%n_grp)
145 deallocate (hecmesh%elem_group%grp_item)
146 allocate (hecmesh%elem_group%grp_item(icou))
147
148 hecmesh%elem_group%grp_index= 0
149 do ig= 1, hecmesh%elem_group%n_grp
150 hecmesh%elem_group%grp_index(ig)= &
151 & hecmesh%elem_group%grp_index(ig-1) + iw2(ig)
152 enddo
153 do ig= 1, hecmesh%elem_group%n_grp
154 do k= 1, iw2(ig)
155 kk= hecmesh%elem_group%grp_index(ig-1) + k
156 hecmesh%elem_group%grp_item(kk)= iupd(k,ig)
157 enddo
158 enddo
159 nnn2= hecmesh%elem_group%grp_index(hecmesh%elem_group%n_grp)
160 ! write (*,'(a, i5, 2i8)') ' elem group', hecMESH%my_rank, nnn1, nnn2
161
162 deallocate (iupd, iw1, iw2)
163 !C===
164 endif
165 endif
166
167 if (hecmesh%surf_group%n_grp.ne.0) then
168 if (hecmesh%surf_group%grp_index(hecmesh%surf_group%n_grp).ne.0) then
169 !C
170 !C +---------------+
171 !C | SURFACE-GROUP |
172 !C +---------------+
173 !C===
174 allocate (iupd(2*hecmesh%n_node, &
175 & hecmesh%surf_group%n_grp))
176 allocate (iw1(2*hecmesh%n_elem))
177 allocate (iw2(hecmesh%surf_group%n_grp))
178 allocate (iw3(hecmesh%n_node))
179
180 iupd = 0
181 icou = 0
182 iw2 = 0
183
184 do ig= 1, hecmesh%surf_group%n_grp
185 icou0= 0
186 iw1= 0
187 iw3= 0
188 do k= hecmesh%surf_group%grp_index(ig-1)+1, &
189 & hecmesh%surf_group%grp_index(ig)
190 icel = hecmesh%surf_group%grp_item(2*k-1)
191 isuf = hecmesh%surf_group%grp_item(2*k )
192 icou = icou + 1
193 icou0 = icou0 + 1
194 iw1(icel )= isuf
195 iupd(2*icou0-1,ig)= icel
196 iupd(2*icou0 ,ig)= isuf
197 ip_type= hecmesh%elem_type(icel)
198 if (ip_type.eq.351) then
199 isp= hecmesh%elem_node_index(icel-1)
200 np1= hecmesh%elem_node_item (isp+1)
201 np2= hecmesh%elem_node_item (isp+2)
202 np3= hecmesh%elem_node_item (isp+3)
203 np4= hecmesh%elem_node_item (isp+4)
204 np5= hecmesh%elem_node_item (isp+5)
205 np6= hecmesh%elem_node_item (isp+6)
206
207 if (isuf.eq.1) then
208 iw3(np2)= 1
209 iw3(np3)= 1
210 iw3(np5)= 1
211 iw3(np6)= 1
212 else if (isuf.eq.2) then
213 iw3(np3)= 1
214 iw3(np1)= 1
215 iw3(np6)= 1
216 iw3(np4)= 1
217 else if (isuf.eq.3) then
218 iw3(np1)= 1
219 iw3(np2)= 1
220 iw3(np4)= 1
221 iw3(np5)= 1
222 else if (isuf.eq.4) then
223 iw3(np1)= 1
224 iw3(np2)= 1
225 iw3(np3)= 1
226 else if (isuf.eq.5) then
227 iw3(np4)= 1
228 iw3(np5)= 1
229 iw3(np6)= 1
230 endif
231 endif
232
233 if (ip_type.eq.341) then
234 isp= hecmesh%elem_node_index(icel-1)
235 np1= hecmesh%elem_node_item (isp+1)
236 np2= hecmesh%elem_node_item (isp+2)
237 np3= hecmesh%elem_node_item (isp+3)
238 np4= hecmesh%elem_node_item (isp+4)
239 if (isuf.eq.1) then
240 iw3(np2)= 1
241 iw3(np3)= 1
242 iw3(np4)= 1
243 else if (isuf.eq.2) then
244 iw3(np1)= 1
245 iw3(np3)= 1
246 iw3(np4)= 1
247 else if (isuf.eq.3) then
248 iw3(np1)= 1
249 iw3(np2)= 1
250 iw3(np4)= 1
251 else if (isuf.eq.4) then
252 iw3(np1)= 1
253 iw3(np2)= 1
254 iw3(np3)= 1
255 endif
256 endif
257
258 enddo
259
260 do ie= 1, hecmesh%n_adapt_edge
261 nod1= hecmesh%adapt_edge_node(2*ie-1)
262 nod2= hecmesh%adapt_edge_node(2*ie )
263
264 if (iw3(nod1).eq.1 .and. iw3(nod2).eq.1 .and. &
265 & hecmesh%adapt_iemb(ie).ne.0) then
266 nod3 = hecmesh%adapt_IWK(ie)
267 iw3(nod3)= 1
268 endif
269 enddo
270
271 do kkk= hecmesh%surf_group%grp_index(ig-1)+1, &
272 & hecmesh%surf_group%grp_index(ig)
273 icel= hecmesh%surf_group%grp_item(2*kkk-1)
274 if (hecmesh%adapt_type(icel).ne.0.and.iw1(icel).ne.0) then
275 ip_type= hecmesh%elem_type(icel)
276 is= hecmesh%adapt_children_index(icel-1) + 1
277 ie= hecmesh%adapt_children_index(icel)
278 ics= hecmesh%adapt_children_local(is)
279
280 if (hecmesh%when_i_was_refined_elem(ics).eq. &
281 & hecmesh%n_adapt) then
282 do k= is, ie
283 if (hecmesh%adapt_children_item(2*k-1).ne.0) then
284 iclocal= hecmesh%adapt_children_local(k)
285 if (ip_type.eq.351) then
286 isc= hecmesh%elem_node_index(iclocal-1)
287 nc1= hecmesh%elem_node_item (isc+1)
288 nc2= hecmesh%elem_node_item (isc+2)
289 nc3= hecmesh%elem_node_item (isc+3)
290 nc4= hecmesh%elem_node_item (isc+4)
291 nc5= hecmesh%elem_node_item (isc+5)
292 nc6= hecmesh%elem_node_item (isc+6)
293
294 nsuf1= iw3(nc2) + iw3(nc3) + iw3(nc5) + iw3(nc6)
295 nsuf2= iw3(nc1) + iw3(nc3) + iw3(nc4) + iw3(nc6)
296 nsuf3= iw3(nc1) + iw3(nc2) + iw3(nc4) + iw3(nc5)
297 nsuf4= iw3(nc1) + iw3(nc2) + iw3(nc3)
298 nsuf5= iw3(nc4) + iw3(nc5) + iw3(nc6)
299
300 if (nsuf1.eq.4) then
301 icou = icou + 1
302 icou0 = icou0 + 1
303 iw1(iclocal )= 1
304 iupd(2*icou0-1,ig)= iclocal
305 iupd(2*icou0 ,ig)= 1
306 endif
307
308 if (nsuf2.eq.4) then
309 icou = icou + 1
310 icou0 = icou0 + 1
311 iw1(iclocal )= 2
312 iupd(2*icou0-1,ig)= iclocal
313 iupd(2*icou0 ,ig)= 2
314 endif
315
316 if (nsuf3.eq.4) then
317 icou = icou + 1
318 icou0 = icou0 + 1
319 iw1(iclocal )= 3
320 iupd(2*icou0-1,ig)= iclocal
321 iupd(2*icou0 ,ig)= 3
322 endif
323
324 if (nsuf4.eq.3) then
325 icou = icou + 1
326 icou0 = icou0 + 1
327 iw1(iclocal )= 4
328 iupd(2*icou0-1,ig)= iclocal
329 iupd(2*icou0 ,ig)= 4
330 endif
331
332 if (nsuf5.eq.3) then
333 icou = icou + 1
334 icou0 = icou0 + 1
335 iw1(iclocal )= 5
336 iupd(2*icou0-1,ig)= iclocal
337 iupd(2*icou0 ,ig)= 5
338 endif
339 endif
340
341 if (ip_type.eq.341) then
342 isc= hecmesh%elem_node_index(iclocal-1)
343 nc1= hecmesh%elem_node_item (isc+1)
344 nc2= hecmesh%elem_node_item (isc+2)
345 nc3= hecmesh%elem_node_item (isc+3)
346 nc4= hecmesh%elem_node_item (isc+4)
347
348 nsuf1= iw3(nc2) + iw3(nc3) + iw3(nc4)
349 nsuf2= iw3(nc1) + iw3(nc3) + iw3(nc4)
350 nsuf3= iw3(nc1) + iw3(nc2) + iw3(nc4)
351 nsuf4= iw3(nc1) + iw3(nc2) + iw3(nc3)
352
353 if (nsuf1.eq.3) then
354 icou = icou + 1
355 icou0 = icou0 + 1
356 iw1(iclocal )= 1
357 iupd(2*icou0-1,ig)= iclocal
358 iupd(2*icou0 ,ig)= 1
359 endif
360
361 if (nsuf2.eq.3) then
362 icou = icou + 1
363 icou0 = icou0 + 1
364 iw1(iclocal )= 2
365 iupd(2*icou0-1,ig)= iclocal
366 iupd(2*icou0 ,ig)= 2
367 endif
368
369 if (nsuf3.eq.3) then
370 icou = icou + 1
371 icou0 = icou0 + 1
372 iw1(iclocal )= 3
373 iupd(2*icou0-1,ig)= iclocal
374 iupd(2*icou0 ,ig)= 3
375 endif
376
377 if (nsuf4.eq.3) then
378 icou = icou + 1
379 icou0 = icou0 + 1
380 iw1(iclocal )= 4
381 iupd(2*icou0-1,ig)= iclocal
382 iupd(2*icou0 ,ig)= 4
383 endif
384
385 endif
386 endif
387 enddo
388 endif
389 endif
390 enddo
391 iw2(ig)= icou0
392 enddo
393
394 !C
395 !C-- NEW POINTERs
396 nnn1= hecmesh%surf_group%grp_index(hecmesh%surf_group%n_grp)
397 deallocate (hecmesh%surf_group%grp_item)
398 allocate (hecmesh%surf_group%grp_item(2*icou))
399
400 hecmesh%surf_group%grp_index= 0
401 do ig= 1, hecmesh%surf_group%n_grp
402 hecmesh%surf_group%grp_index(ig)= &
403 & hecmesh%surf_group%grp_index(ig-1) + iw2(ig)
404 enddo
405 do ig= 1, hecmesh%surf_group%n_grp
406 do k= 1, iw2(ig)
407 kk= hecmesh%surf_group%grp_index(ig-1) + k
408 hecmesh%surf_group%grp_item(2*kk-1)= iupd(2*k-1,ig)
409 hecmesh%surf_group%grp_item(2*kk )= iupd(2*k ,ig)
410 enddo
411 enddo
412 nnn2= hecmesh%surf_group%grp_index(hecmesh%surf_group%n_grp)
413 ! write (*,'(a, i5, 2i8)') ' surf group', hecMESH%my_rank, nnn1, nnn2
414
415 deallocate (iupd, iw1, iw2, iw3)
416 !C===
417 endif
418 endif
419
420 return
421end
422
423
424
subroutine hecmw_adapt_bc_pointer(hecmesh)
Adaptive Mesh Refinement.
I/O and Utility.
Definition: hecmw_util_f.F90:7