FrontISTR  5.7.0
Large-scale structural analysis program with finit element method
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
14 subroutine 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
421 end
422 
423 
424 
hecmw_util
I/O and Utility.
Definition: hecmw_util_f.F90:7
hecmw_util::hecmwst_local_mesh
Definition: hecmw_util_f.F90:234
hecmw_adapt_bc_pointer
subroutine hecmw_adapt_bc_pointer(hecMESH)
Adaptive Mesh Refinement.
Definition: hecmw_adapt_bc_pointer.f90:15