FrontISTR  5.7.0
Large-scale structural analysis program with finit element method
hecmw_adapt_grid_smooth.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_GRID_SMOOTH
10 !C***
11 !C
12 !C grid smoothing for embedding
13 !C
14 subroutine hecmw_adapt_grid_smooth (hecMESH)
15  use hecmw_util
16  use hecmw_adapt_int_sr
17  use hecmw_adapt_int_sr_rev
18 
19  implicit real*8 (a-h,o-z)
20  integer(kind=kint), dimension(:), allocatable :: WR, WS
21 
22  dimension ndiv(6)
23  integer(kind=kint), dimension(:), allocatable :: NFLAG_INFO
24 
25  type (hecmwST_local_mesh) :: hecMESH
26 
27  !C
28  !C-- INIT.
29  niteradap= 0
30  newtet = 0
31  newprism = 0
32  allocate (nflag_info(hecmesh%PETOT))
33 
34  niteradap_max= 10000
35 
36  90 continue
37  call mpi_barrier (hecmesh%MPI_COMM, ierr)
38 
39  niteradap= niteradap + 1
40  nf0 = 0
41  if ( niteradap .gt. niteradap_max) then
42  call hecmw_adapt_error_exit(hecmesh, 7)
43  endif
44 
45  !C
46  !C +------------+
47  !C | TETRAHEDRA |
48  !C +------------+
49  !C
50  !C RULEs :
51  !C a. READ the PAPER for the PATTERNs
52  !C b. CONSECUTIVE DIRECTIONAL REFINEMENT IS NOT ALLOWED
53  !C***********************************************************************
54 
55  if (hecmesh%my_rank.eq.0) &
56  & write (*,'(" TETRA iteration=", 2i8)') niteradap, newtet
57 
58  newtet= 0
59  do 100 icel0= 1, hecmesh%n_adapt_act_elem_341
60  icel= hecmesh%adapt_act_elem_341(icel0)
61  is= hecmesh%elem_node_index(icel-1)
62  n1= hecmesh%elem_node_item (is+1)
63  n2= hecmesh%elem_node_item (is+2)
64  n3= hecmesh%elem_node_item (is+3)
65  n4= hecmesh%elem_node_item (is+4)
66 
67  call hecmw_adapt_edge_info ( hecmesh, n1, n2, ie1, 1 )
68  call hecmw_adapt_edge_info ( hecmesh, n1, n3, ie2, 1 )
69  call hecmw_adapt_edge_info ( hecmesh, n1, n4, ie3, 1 )
70  call hecmw_adapt_edge_info ( hecmesh, n2, n3, ie4, 1 )
71  call hecmw_adapt_edge_info ( hecmesh, n2, n4, ie5, 1 )
72  call hecmw_adapt_edge_info ( hecmesh, n3, n4, ie6, 1 )
73 
74  ndiv(1)= hecmesh%adapt_iemb(ie1)
75  ndiv(2)= hecmesh%adapt_iemb(ie2)
76  ndiv(3)= hecmesh%adapt_iemb(ie3)
77  ndiv(4)= hecmesh%adapt_iemb(ie4)
78  ndiv(5)= hecmesh%adapt_iemb(ie5)
79  ndiv(6)= hecmesh%adapt_iemb(ie6)
80 
81  ndivsum= ndiv(1)+ndiv(2)+ndiv(3)+ndiv(4)+ndiv(5)+ndiv(6)
82 
83  !C
84  !C +--------------------------+
85  !C | ADJUST the CELL DIVISION |
86  !C +--------------------------+
87  !C===
88  if (ndivsum .eq. 1) goto 95
89  !C
90  !C== 2 edges
91  if (ndivsum .eq. 2) then
92  if ( ( (ndiv(1).eq.1) .and. (ndiv(6).eq.1) ) .or. &
93  & ( (ndiv(2).eq.1) .and. (ndiv(5).eq.1) ) .or. &
94  & ( (ndiv(3).eq.1) .and. (ndiv(4).eq.1) ) ) then
95  hecmesh%adapt_iemb(ie1)= 1
96  hecmesh%adapt_iemb(ie2)= 1
97  hecmesh%adapt_iemb(ie3)= 1
98  hecmesh%adapt_iemb(ie4)= 1
99  hecmesh%adapt_iemb(ie5)= 1
100  hecmesh%adapt_iemb(ie6)= 1
101  newtet = newtet + 1
102  nf0 = 1
103  goto 95
104  !C
105  !C-- PATTERN 1-3-5
106  else if &
107  & ( (ndiv(1).eq.1) .and. (ndiv(3).eq.1) ) then
108  hecmesh%adapt_iemb(ie5)= 1
109  newtet = newtet + 1
110  nf0 = 1
111  goto 95
112  else if &
113  & ( (ndiv(1).eq.1) .and. (ndiv(5).eq.1) ) then
114  hecmesh%adapt_iemb(ie3)= 1
115  newtet = newtet + 1
116  nf0 = 1
117  goto 95
118  else if &
119  & ( (ndiv(3).eq.1) .and. (ndiv(5).eq.1) ) then
120  hecmesh%adapt_iemb(ie1)= 1
121  newtet = newtet + 1
122  nf0 = 1
123  goto 95
124  !C
125  !C-- PATTERN 2-3-6
126  else if &
127  & ( (ndiv(2).eq.1) .and. (ndiv(3).eq.1) ) then
128  hecmesh%adapt_iemb(ie6)= 1
129  newtet = newtet + 1
130  nf0 = 1
131  goto 95
132  else if &
133  & ( (ndiv(2).eq.1) .and. (ndiv(6).eq.1) ) then
134  hecmesh%adapt_iemb(ie3)= 1
135  newtet = newtet + 1
136  nf0 = 1
137  goto 95
138  else if &
139  & ( (ndiv(3).eq.1) .and. (ndiv(6).eq.1) ) then
140  hecmesh%adapt_iemb(ie2)= 1
141  newtet = newtet + 1
142  nf0 = 1
143  goto 95
144  !C
145  !C-- PATTERN 1-2-4
146  else if &
147  & ( (ndiv(1).eq.1) .and. (ndiv(2).eq.1) ) then
148  hecmesh%adapt_iemb(ie4)= 1
149  newtet = newtet + 1
150  nf0 = 1
151  goto 95
152  else if &
153  & ( (ndiv(1).eq.1) .and. (ndiv(4).eq.1) ) then
154  hecmesh%adapt_iemb(ie2)= 1
155  newtet = newtet + 1
156  nf0 = 1
157  goto 95
158  else if &
159  & ( (ndiv(2).eq.1) .and. (ndiv(4).eq.1) ) then
160  hecmesh%adapt_iemb(ie1)= 1
161  newtet = newtet + 1
162  nf0 = 1
163  goto 95
164  !C
165  !C-- PATTERN 4-5-6
166  else if &
167  & ( (ndiv(4).eq.1) .and. (ndiv(5).eq.1) ) then
168  hecmesh%adapt_iemb(ie6)= 1
169  newtet = newtet + 1
170  nf0 = 1
171  goto 95
172  else if &
173  & ( (ndiv(4).eq.1) .and. (ndiv(6).eq.1) ) then
174  hecmesh%adapt_iemb(ie5)= 1
175  newtet = newtet + 1
176  nf0 = 1
177  goto 95
178  else if &
179  & ( (ndiv(5).eq.1) .and. (ndiv(6).eq.1) ) then
180  hecmesh%adapt_iemb(ie4)= 1
181  newtet = newtet + 1
182  nf0 = 1
183  goto 95
184  else
185  hecmesh%adapt_iemb(ie1)= 1
186  hecmesh%adapt_iemb(ie2)= 1
187  hecmesh%adapt_iemb(ie3)= 1
188  hecmesh%adapt_iemb(ie4)= 1
189  hecmesh%adapt_iemb(ie5)= 1
190  hecmesh%adapt_iemb(ie6)= 1
191  newtet = newtet + 1
192  nf0 = 1
193  goto 95
194  endif
195  endif
196  !C==
197 
198  !C
199  !C== 3 edges
200  if (ndivsum .eq. 3) then
201  if ( &
202  & ((ndiv(1).eq.1).and.(ndiv(3).eq.1).and.(ndiv(5).eq.1)) .or.&
203  & ((ndiv(2).eq.1).and.(ndiv(3).eq.1).and.(ndiv(6).eq.1)) .or.&
204  & ((ndiv(1).eq.1).and.(ndiv(2).eq.1).and.(ndiv(4).eq.1)) .or.&
205  & ((ndiv(4).eq.1).and.(ndiv(5).eq.1).and.(ndiv(6).eq.1)) ) &
206  & then
207  goto 95
208  else
209  hecmesh%adapt_iemb(ie1)= 1
210  hecmesh%adapt_iemb(ie2)= 1
211  hecmesh%adapt_iemb(ie3)= 1
212  hecmesh%adapt_iemb(ie4)= 1
213  hecmesh%adapt_iemb(ie5)= 1
214  hecmesh%adapt_iemb(ie6)= 1
215  newtet = newtet + 1
216  nf0 = 1
217  goto 95
218  endif
219  endif
220  !C==
221 
222  !C
223  !C== more than 4 edges
224  if (ndivsum.eq.4 .or. ndivsum.eq.5) then
225  hecmesh%adapt_iemb(ie1)= 1
226  hecmesh%adapt_iemb(ie2)= 1
227  hecmesh%adapt_iemb(ie3)= 1
228  hecmesh%adapt_iemb(ie4)= 1
229  hecmesh%adapt_iemb(ie5)= 1
230  hecmesh%adapt_iemb(ie6)= 1
231  newtet = newtet + 1
232  nf0 = 1
233  goto 95
234  endif
235  !C==
236 
237  !C
238  !C== check the type of PARENT cell
239  95 continue
240 
241  ntyp= hecmesh%adapt_parent_type(icel)
242  if (ntyp.ne.0 .and. ntyp.ne.11 .and. &
243  & ndivsum.ne.0.and.ndivsum.ne.6) then
244  hecmesh%adapt_iemb(ie1)= 1
245  hecmesh%adapt_iemb(ie2)= 1
246  hecmesh%adapt_iemb(ie3)= 1
247  hecmesh%adapt_iemb(ie4)= 1
248  hecmesh%adapt_iemb(ie5)= 1
249  hecmesh%adapt_iemb(ie6)= 1
250  nf0 = 1
251  endif
252 
253  !C==
254  100 continue
255  !C***********************************************************************
256 
257  !C
258  !C +--------+
259  !C | PRISMs |
260  !C +--------+
261  !C
262  !C RULEs :
263  !C a. READ the PAPER for the PATTERNs
264  !C b. CONSECUTIVE DIRECTIONAL REFINEMENT IS NOT ALLOWED
265  !C***********************************************************************
266 
267  if (hecmesh%my_rank.eq.0) &
268  & write (*,'(" PRISM iteration=", 2i8)') niteradap, newprism
269 
270  !C
271  !C-- ADJUST normal-to-surface direction
272  n1= hecmesh%adapt_import_edge_index(hecmesh%n_neighbor_pe)
273  n2= hecmesh%adapt_export_edge_index(hecmesh%n_neighbor_pe)
274  m = max(n1, n2)
275  allocate (ws(m), wr(m))
276 
277  icoum= 0
278  do layer= 1, hecmesh%n_adapt_act_elem_351
279  icou = 0
280  do icel0= 1, hecmesh%n_adapt_act_elem_351
281  icel= hecmesh%adapt_act_elem_351(icel0)
282  is= hecmesh%elem_node_index(icel-1)
283  n1= hecmesh%elem_node_item (is+1)
284  n2= hecmesh%elem_node_item (is+2)
285  n3= hecmesh%elem_node_item (is+3)
286  n4= hecmesh%elem_node_item (is+4)
287  n5= hecmesh%elem_node_item (is+5)
288  n6= hecmesh%elem_node_item (is+6)
289 
290  call hecmw_adapt_edge_info ( hecmesh, n1, n2, ie1, 1 )
291  call hecmw_adapt_edge_info ( hecmesh, n2, n3, ie2, 1 )
292  call hecmw_adapt_edge_info ( hecmesh, n3, n1, ie3, 1 )
293  call hecmw_adapt_edge_info ( hecmesh, n4, n5, ie4, 1 )
294  call hecmw_adapt_edge_info ( hecmesh, n5, n6, ie5, 1 )
295  call hecmw_adapt_edge_info ( hecmesh, n6, n4, ie6, 1 )
296 
297  if (hecmesh%adapt_iemb(ie1).eq.1 .and. &
298  & hecmesh%adapt_iemb(ie4).eq.0) then
299  hecmesh%adapt_iemb(ie4)= 1
300  icou = icou + 1
301  endif
302  if (hecmesh%adapt_iemb(ie2).eq.1 .and. &
303  & hecmesh%adapt_iemb(ie5).eq.0) then
304  hecmesh%adapt_iemb(ie5)= 1
305  icou = icou + 1
306  endif
307  if (hecmesh%adapt_iemb(ie3).eq.1 .and. &
308  & hecmesh%adapt_iemb(ie6).eq.0) then
309  hecmesh%adapt_iemb(ie6)= 1
310  icou = icou + 1
311  endif
312  if (hecmesh%adapt_iemb(ie4).eq.1 .and. &
313  & hecmesh%adapt_iemb(ie1).eq.0) then
314  hecmesh%adapt_iemb(ie1)= 1
315  icou = icou + 1
316  endif
317 
318  if (hecmesh%adapt_iemb(ie5).eq.1 .and. &
319  & hecmesh%adapt_iemb(ie2).eq.0) then
320  hecmesh%adapt_iemb(ie2)= 1
321  icou = icou + 1
322  endif
323  if (hecmesh%adapt_iemb(ie6).eq.1 .and. &
324  & hecmesh%adapt_iemb(ie3).eq.0) then
325  hecmesh%adapt_iemb(ie3)= 1
326  icou = icou + 1
327  endif
328  enddo
329 
330  call mpi_allreduce ( icoum, icoummin, 1, mpi_integer, &
331  & mpi_min, hecmesh%MPI_COMM, ierr)
332 
333  if (icou.eq.0 .and. icoummin.eq.1) exit
334  if (icou.eq.0) icoum= 1
335  if (icou.ne.0) icoum= 0
336 
337  ws= 0
338  wr= 0
339  call hecmw_adapt_int_send_recv_rev &
340  & ( hecmesh%n_adapt_edge, &
341  & hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
342  & hecmesh%adapt_import_edge_index, &
343  & hecmesh%adapt_import_edge_item , &
344  & hecmesh%adapt_export_edge_index, &
345  & hecmesh%adapt_export_edge_item , &
346  & ws, wr, hecmesh%adapt_iemb, hecmesh%MPI_COMM, &
347  & hecmesh%my_rank, 1, m)
348 
349  ws= 0
350  wr= 0
351  call hecmw_adapt_int_send_recv &
352  & ( hecmesh%n_adapt_edge, &
353  & hecmesh%n_neighbor_pe, hecmesh%neighbor_pe, &
354  & hecmesh%adapt_import_edge_index, &
355  & hecmesh%adapt_import_edge_item , &
356  & hecmesh%adapt_export_edge_index, &
357  & hecmesh%adapt_export_edge_item , &
358  & ws, wr, hecmesh%adapt_iemb, hecmesh%MPI_COMM, &
359  & hecmesh%my_rank, 1, m)
360 
361  enddo
362  deallocate (ws, wr)
363 
364  call mpi_barrier (hecmesh%MPI_COMM,ierr)
365 
366  !C
367  !C-- real operation for prismatic region
368  newprism= 0
369  do 110 icel0= 1, hecmesh%n_adapt_act_elem_351
370  icel= hecmesh%adapt_act_elem_351(icel0)
371  is= hecmesh%elem_node_index(icel-1)
372  n1= hecmesh%elem_node_item (is+1)
373  n2= hecmesh%elem_node_item (is+2)
374  n3= hecmesh%elem_node_item (is+3)
375  n4= hecmesh%elem_node_item (is+4)
376  n5= hecmesh%elem_node_item (is+5)
377  n6= hecmesh%elem_node_item (is+6)
378 
379  call hecmw_adapt_edge_info ( hecmesh, n1, n2, ie1, 1 )
380  call hecmw_adapt_edge_info ( hecmesh, n2, n3, ie2, 1 )
381  call hecmw_adapt_edge_info ( hecmesh, n3, n1, ie3, 1 )
382  call hecmw_adapt_edge_info ( hecmesh, n4, n5, ie4, 1 )
383  call hecmw_adapt_edge_info ( hecmesh, n5, n6, ie5, 1 )
384  call hecmw_adapt_edge_info ( hecmesh, n6, n4, ie6, 1 )
385 
386  ndiv(1)= hecmesh%adapt_iemb(ie1)
387  ndiv(2)= hecmesh%adapt_iemb(ie2)
388  ndiv(3)= hecmesh%adapt_iemb(ie3)
389  ndiv(4)= hecmesh%adapt_iemb(ie4)
390  ndiv(5)= hecmesh%adapt_iemb(ie5)
391  ndiv(6)= hecmesh%adapt_iemb(ie6)
392 
393  ndivsum= ndiv(1)+ndiv(2)+ndiv(3)+ndiv(4)+ndiv(5)+ndiv(6)
394 
395  !C
396  !C +--------------------------+
397  !C | ADJUST the CELL DIVISION |
398  !C +--------------------------+
399  !C===
400 
401  !C
402  !C-- 1 edge(s)
403  if (ndivsum .eq. 1) then
404  if (ndiv(1).eq.1) hecmesh%adapt_iemb(ie4)= 1
405  if (ndiv(2).eq.1) hecmesh%adapt_iemb(ie5)= 1
406  if (ndiv(3).eq.1) hecmesh%adapt_iemb(ie6)= 1
407  if (ndiv(4).eq.1) hecmesh%adapt_iemb(ie1)= 1
408  if (ndiv(5).eq.1) hecmesh%adapt_iemb(ie2)= 1
409  if (ndiv(6).eq.1) hecmesh%adapt_iemb(ie3)= 1
410  goto 105
411  endif
412  !C
413  !C-- 2 edges
414  if (ndivsum.eq.2) then
415  if (ndiv(1).eq.1 .and. ndiv(4).eq.1) goto 105
416  if (ndiv(2).eq.1 .and. ndiv(5).eq.1) goto 105
417  if (ndiv(3).eq.1 .and. ndiv(6).eq.1) goto 105
418  hecmesh%adapt_iemb(ie1)= 1
419  hecmesh%adapt_iemb(ie2)= 1
420  hecmesh%adapt_iemb(ie3)= 1
421  hecmesh%adapt_iemb(ie4)= 1
422  hecmesh%adapt_iemb(ie5)= 1
423  hecmesh%adapt_iemb(ie6)= 1
424  goto 105
425  endif
426 
427  !C
428  !C-- >3 edges
429  if (ndivsum.ge.3) then
430  hecmesh%adapt_iemb(ie1)= 1
431  hecmesh%adapt_iemb(ie2)= 1
432  hecmesh%adapt_iemb(ie3)= 1
433  hecmesh%adapt_iemb(ie4)= 1
434  hecmesh%adapt_iemb(ie5)= 1
435  hecmesh%adapt_iemb(ie6)= 1
436  goto 105
437  endif
438  !C
439  !C-- check the type of PARENT cell
440  105 continue
441 
442  ntyp= hecmesh%adapt_parent_type(icel)
443  if (ntyp.ne.0 .and. ntyp.ne.4 .and. &
444  & ndivsum.ne.0.and.ndivsum.ne.6) then
445  hecmesh%adapt_iemb(ie1)= 1
446  hecmesh%adapt_iemb(ie2)= 1
447  hecmesh%adapt_iemb(ie3)= 1
448  hecmesh%adapt_iemb(ie4)= 1
449  hecmesh%adapt_iemb(ie5)= 1
450  hecmesh%adapt_iemb(ie6)= 1
451  nf0 = 1
452  endif
453 
454  !C==
455  110 continue
456  !C***********************************************************************
457 
458  !C
459  !C-- adjust EMBEDDING level
460  call mpi_barrier (hecmesh%MPI_COMM, ierr)
461  call hecmw_adapt_adjemb (hecmesh, nf0)
462 
463  call mpi_barrier (hecmesh%MPI_COMM, ierr)
464  call mpi_gather &
465  & (nf0, 1, mpi_integer, nflag_info, 1, mpi_integer, 0, &
466  & hecmesh%MPI_COMM, ierr)
467 
468  if (hecmesh%my_rank.eq.0) then
469  icou= 0
470  do i= 1, hecmesh%PETOT
471  icou= icou + nflag_info(i)
472  enddo
473  if (icou.ne.0) nf0= 1
474  endif
475 
476  call mpi_bcast (nf0, 1, mpi_integer, 0, hecmesh%MPI_COMM, ierr)
477  if (nf0 .eq. 1) goto 90
478 
479  return
480 end
481 
hecmw_adapt_error_exit
subroutine hecmw_adapt_error_exit(hecMESH, IFLAG)
Adaptive Mesh Refinement.
Definition: hecmw_adapt_error_exit.f90:8
hecmw_adapt_int_sr_rev
Adaptive Mesh Refinement.
Definition: hecmw_adapt_int_sr_rev.f90:7
hecmw_adapt_int_sr
Adaptive Mesh Refinement.
Definition: hecmw_adapt_int_sr.f90:7
hecmw_adapt_int_sr_rev::hecmw_adapt_int_send_recv_rev
subroutine hecmw_adapt_int_send_recv_rev(N, NEIBPETOT, NEIBPE, STACK_EXPORT, NOD_EXPORT, STACK_IMPORT, NOD_IMPORT, WS, WR, X, SOLVER_COMM, my_rank, NB, m)
Definition: hecmw_adapt_int_sr_rev.f90:18
hecmw_util
I/O and Utility.
Definition: hecmw_util_f.F90:7
hecmw_adapt_grid_smooth
subroutine hecmw_adapt_grid_smooth(hecMESH)
Adaptive Mesh Refinement.
Definition: hecmw_adapt_grid_smooth.f90:15
hecmw_adapt_adjemb
subroutine hecmw_adapt_adjemb(hecMESH, NFLAG_INFO)
Adaptive Mesh Refinement.
Definition: hecmw_adapt_adjemb.f90:19
hecmw_adapt_int_sr::hecmw_adapt_int_send_recv
subroutine hecmw_adapt_int_send_recv(N, NEIBPETOT, NEIBPE, STACK_IMPORT, NOD_IMPORT, STACK_EXPORT, NOD_EXPORT, WS, WR, X, SOLVER_COMM, my_rank, NB, m)
Definition: hecmw_adapt_int_sr.f90:18
hecmw_adapt_edge_info
subroutine hecmw_adapt_edge_info(hecMESH, nod1, nod2, iedge, NFLAG)
Adaptive Mesh Refinement.
Definition: hecmw_adapt_edge_info.f90:13