All Classes Files Functions Variables Pages
module_init_mesh.f90
Go to the documentation of this file.
1 !=====================================================================================================================================
2 ! EFISPEC3D !
3 ! (Elements FInis SPECtraux 3D) !
4 ! !
5 ! http://efispec.free.fr !
6 ! !
7 ! !
8 ! This file is part of EFISPEC3D !
9 ! Please refer to http://efispec.free.fr if you use it or part of it !
10 ! !
11 ! !
12 !1 ---> French License: CeCILL V2 !
13 ! !
14 ! Copyright BRGM 2009 contributeurs : Florent DE MARTIN !
15 ! David MICHEA !
16 ! Philippe THIERRY !
17 ! !
18 ! Contact: f.demartin at brgm.fr !
19 ! !
20 ! Ce logiciel est un programme informatique servant a resoudre l'equation du !
21 ! mouvement en trois dimensions via une methode des elements finis spectraux. !
22 ! !
23 ! Ce logiciel est regi par la licence CeCILL soumise au droit francais et !
24 ! respectant les principes de diffusion des logiciels libres. Vous pouvez !
25 ! utiliser, modifier et/ou redistribuer ce programme sous les conditions de la !
26 ! licence CeCILL telle que diffusee par le CEA, le CNRS et l'INRIA sur le site !
27 ! "http://www.cecill.info". !
28 ! !
29 ! En contrepartie de l'accessibilite au code source et des droits de copie, de !
30 ! modification et de redistribution accordes par cette licence, il n'est offert !
31 ! aux utilisateurs qu'une garantie limitee. Pour les memes raisons, seule une !
32 ! responsabilite restreinte pese sur l'auteur du programme, le titulaire des !
33 ! droits patrimoniaux et les concedants successifs. !
34 ! !
35 ! A cet egard l'attention de l'utilisateur est attiree sur les risques associes !
36 ! au chargement, a l'utilisation, a la modification et/ou au developpement et a !
37 ! la reproduction du logiciel par l'utilisateur etant donne sa specificite de !
38 ! logiciel libre, qui peut le rendre complexe a manipuler et qui le reserve donc !
39 ! a des developpeurs et des professionnels avertis possedant des connaissances !
40 ! informatiques approfondies. Les utilisateurs sont donc invites a charger et !
41 ! tester l'adequation du logiciel a leurs besoins dans des conditions permettant !
42 ! d'assurer la securite de leurs systemes et ou de leurs donnees et, plus !
43 ! generalement, a l'utiliser et l'exploiter dans les memes conditions de !
44 ! securite. !
45 ! !
46 ! Le fait que vous puissiez acceder a cet en-tete signifie que vous avez pris !
47 ! connaissance de la licence CeCILL et que vous en avez accepte les termes. !
48 ! !
49 ! !
50 !2 ---> International license: GNU GPL V3 !
51 ! !
52 ! EFISPEC3D is a computer program that solves the three-dimensional equations of !
53 ! motion using a finite spectral-element method. !
54 ! !
55 ! Copyright (C) 2009 Florent DE MARTIN !
56 ! !
57 ! Contact: f.demartin at brgm.fr !
58 ! !
59 ! This program is free software: you can redistribute it and/or modify it under !
60 ! the terms of the GNU General Public License as published by the Free Software !
61 ! Foundation, either version 3 of the License, or (at your option) any later !
62 ! version. !
63 ! !
64 ! This program is distributed in the hope that it will be useful, but WITHOUT ANY !
65 ! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A !
66 ! PARTICULAR PURPOSE. See the GNU General Public License for more details. !
67 ! !
68 ! You should have received a copy of the GNU General Public License along with !
69 ! this program. If not, see http://www.gnu.org/licenses/. !
70 ! !
71 ! !
72 !3 ---> Third party libraries !
73 ! !
74 ! EFISPEC3D uses the following of third party libraries: !
75 ! !
76 ! --> METIS 5.1.0 !
77 ! see http://glaros.dtc.umn.edu/gkhome/metis/metis/overview !
78 ! !
79 ! --> Lib_VTK_IO !
80 ! see S. Zaghi's website: https://github.com/szaghi/Lib_VTK_IO !
81 ! !
82 ! --> INTERP_LINEAR !
83 ! see J. Burkardt website: http://people.sc.fsu.edu/~jburkardt/ !
84 ! !
85 ! --> SAC !
86 ! http://ds.iris.edu/files/sac-manual/ !
87 ! !
88 ! --> EXODUS II !
89 ! http://sourceforge.net/projects/exodusii/ !
90 ! !
91 ! --> NETCDF !
92 ! http://www.unidata.ucar.edu/software/netcdf/ !
93 ! !
94 ! --> HDF5 !
95 ! http://www.hdfgroup.org/HDF5/ !
96 ! !
97 ! Some of these libraries are located in directory lib and pre-compiled !
98 ! with intel compiler for x86_64 platform. !
99 ! !
100 ! !
101 !4 ---> Related Articles !
102 ! !
103 ! De Martin, F., Matsushima, M., Kawase, H. (BSSA, 2013) !
104 ! Impact of geometric effects on near-surface Green's functions !
105 ! doi:10.1785/0120130039 !
106 ! !
107 ! Aochi, H., Ducellier, A., Dupros, F., Delatre, M., Ulrich, T., De Martin, F., !
108 ! Yoshimi, M., (Pure Appl. Geophys. 2013) !
109 ! Finite Difference Simulations of Seismic Wave Propagation for the 2007 Mw 6.6 !
110 ! Niigata-ken Chuetsu-Oki Earthquake: Validity of Models and Reliable !
111 ! Input Ground Motion in the Near-Field !
112 ! doi:10.1007/s00024-011-0429-5 !
113 ! !
114 ! De Martin, F. (BSSA, 2011) !
115 ! Verification of a Spectral-Element Method Code for the Southern California !
116 ! Earthquake Center LOH.3 Viscoelastic Case !
117 ! doi:10.1785/0120100305 !
118 ! !
119 !=====================================================================================================================================
120 
123 
126 module mod_init_mesh
127 
128  implicit none
129 
130  private
131 
132  public :: init_mesh
133  private :: init_gll_number
134  private :: propagate_gll_nodes_quad
135  private :: propagate_gll_nodes_face
136  private :: propagate_gll_nodes_edge
137  private :: propagate_gll_nodes_corner
138  private :: init_element
139 
140  contains
141 
142 !
143 !
176 !***********************************************************************************************************************************************************************************
177  subroutine init_mesh()
178 !***********************************************************************************************************************************************************************************
179 
180  use mpi
181 
182  use mod_global_variables, only : &
183  ig_nhexa&
184  ,ig_nhexa_outer&
185  ,ig_nhexa_inner&
186  ,ig_nquad_parax&
187  ,ig_nquad_fsurf&
188  ,ig_mesh_nnode&
189  ,ig_hexa_nnode&
190  ,ig_quad_nnode&
191  ,ig_hexa_gnode_glonum&
192  ,ig_quadp_gnode_glonum&
193  ,ig_quadf_gnode_glonum&
194  ,rg_gnode_x&
195  ,rg_gnode_y&
196  ,rg_gnode_z&
197  ,rg_mesh_xmax&
198  ,rg_mesh_ymax&
199  ,rg_mesh_zmax&
200  ,rg_mesh_xmin&
201  ,rg_mesh_ymin&
202  ,rg_mesh_zmin&
203  ,rg_gll_displacement&
204  ,rg_gll_velocity&
205  ,rg_gll_acceleration&
206  ,cg_prefix&
207  ,cg_myrank&
208  ,cg_ncpu&
209  ,ig_myrank&
210  ,ig_ncpu_neighbor&
211  ,tg_cpu_neighbor&
212  ,ig_medium_type&
213  ,ig_hexa_gll_glonum&
214  ,ig_quadp_gll_glonum&
215  ,ig_quadf_gll_glonum&
216  ,ig_hexa_material_number&
217  ,ig_quadp_neighbor_hexa&
218  ,ig_quadp_neighbor_hexaface&
219  ,ig_ngll&
220  ,ig_ngll_total&
221  ,get_newunit&
222  ,error_stop&
223  ,ig_cpu_neighbor_info&
224  ,ig_nneighbor_all_kind&
225  ,ig_ncpu&
226  ,ig_lst_unit&
227  ,ig_ndof&
228  ,info_all_cpu&
229  ,lg_output_debug_file
230 
231  use mod_init_memory
232 
233  implicit none
234 
235  integer ,parameter :: nface=6
236  integer ,parameter :: nedge=12
237  integer ,parameter :: nnode=8
238 
239  integer, allocatable, dimension(:) :: cpu_neighbor
240  integer :: ios
241  integer :: size_real_t
242  integer :: icpu
243  integer :: ihexa
244  integer :: inode
245  integer :: myunit_debug
246  integer :: ifsurf
247  character(len=255) :: info
248 
249  integer :: igll
250  integer :: jgll
251  integer :: kgll
252 
253 
254  !
255  !
256  !********************************************************************************************************************************
257  !read CUBIT mesh file and partition it with METIS library. Subroutine read_part_cubit_mesh is in file partCubitMesh.c
258  !********************************************************************************************************************************
259  call read_part_cubit_mesh(trim(cg_prefix)//char(0) &
260  ,ig_ncpu &
261  ,ig_myrank &
262  ,size_real_t &
263  ,ig_mesh_nnode &
264  ,ig_ncpu_neighbor &
265  ,ig_nhexa &
266  ,ig_nhexa_outer &
267  ,ig_nquad_parax &
268  ,ig_nquad_fsurf &
269  ,ig_hexa_nnode &
270  ,ig_quad_nnode &
271  ,rg_mesh_xmax &
272  ,rg_mesh_ymax &
273  ,rg_mesh_zmax &
274  ,rg_mesh_xmin &
275  ,rg_mesh_ymin &
276  ,rg_mesh_zmin )
277 
278  ig_nhexa_inner = ig_nhexa - ig_nhexa_outer
279 
280 
281  !
282  !
283  !********************************************************************************************************************************
284  !write to lst file min/max (x,y,z) coordinates of the entire domain (i.e., for all cpu)
285  !********************************************************************************************************************************
286  if (ig_myrank == 0) then
287  write(ig_lst_unit,'(" " ,/,a )') "boundaries of the entire domain"
288  write(ig_lst_unit,'(" -->", a,E14.7)') " xmin = ",rg_mesh_xmin
289  write(ig_lst_unit,'(" -->", a,E14.7)') " xmax = ",rg_mesh_xmax
290  write(ig_lst_unit,'(" -->", a,E14.7)') " ymin = ",rg_mesh_ymin
291  write(ig_lst_unit,'(" -->", a,E14.7)') " ymax = ",rg_mesh_ymax
292  write(ig_lst_unit,'(" -->", a,E14.7)') " zmin = ",rg_mesh_zmin
293  write(ig_lst_unit,'(" -->", a,E14.7)') " zmax = ",rg_mesh_zmax
294  endif
295 
296  !
297  !
298  !********************************************************************************************************************************
299  !read info about hexa and quad (number, number of geometric nodes per hexa, etc.)
300  !********************************************************************************************************************************
301  write(info,'(a)') "geometric nodes"
302  call info_all_cpu(ig_mesh_nnode,info)
303 
304  write(info,'(a)') " hexahedra"
305  call info_all_cpu(ig_nhexa,info)
306 
307  write(info,'(a)') " outer hexahedra"
308  call info_all_cpu(ig_nhexa_outer,info)
309 
310  write(info,'(a)') " inner hexahedra"
311  call info_all_cpu(ig_nhexa_inner,info)
312 
313  write(info,'(a)') "paraxial quadrangles"
314  call info_all_cpu(ig_nquad_parax,info)
315 
316  write(info,'(a)') "free surface quadrangles"
317  call info_all_cpu(ig_nquad_fsurf,info)
318 
319 
320  !
321  !
322  !********************************************************************************************************************************
323  !initialize convention numbering of finite element
324  !********************************************************************************************************************************
325  call init_element(ig_hexa_nnode,ig_quad_nnode)
326 
327 
328  !
329  !
330  !********************************************************************************************************************************
331  !initialize some arrays
332  !********************************************************************************************************************************
333  allocate(tg_cpu_neighbor(ig_ncpu_neighbor),stat=ios)
334  if (ios /= 0) then
335  write(info,'(a)') "Error in subroutine init_mesh while allocating tg_cpu_neighbor"
336  call error_stop(info)
337  endif
338 
339  ios = init_array_int(cpu_neighbor,ig_ncpu_neighbor,"cpu_neighbor") !DAVID: temporary array used to fill tg_cpu_neighbor
340 
341  ios = init_array_int(ig_cpu_neighbor_info,26*ig_nhexa_outer,3,"ig_cpu_neighbor_info") !DAVID: deallocation in subroutine init_mpi_buffers()
342 
343  ios = init_array_int(ig_hexa_gnode_glonum,ig_nhexa,ig_hexa_nnode,"ig_hexa_gnode_glonum")
344 
345  ios = init_array_int(ig_hexa_gll_glonum,ig_nhexa,ig_ngll,ig_ngll,ig_ngll,"ig_hexa_gll_glonum")
346 
347  if (ig_nquad_parax > 0) then
348 
349  ios = init_array_int(ig_quadp_gnode_glonum,ig_nquad_parax,ig_quad_nnode,"ig_quadp_gnode_glonum")
350 
351  ios = init_array_int(ig_quadp_gll_glonum,ig_nquad_parax,ig_ngll,ig_ngll,"ig_quadp_gll_glonum")
352 
353  ios = init_array_int(ig_quadp_neighbor_hexa,ig_nquad_parax,"ig_quadp_neighbor_hexa")
354 
355  ios = init_array_int(ig_quadp_neighbor_hexaface,ig_nquad_parax,"ig_quadp_neighbor_hexaface")
356 
357  else !workaround for C subroutine 'fill_mesh_arrays' when there is no quad parax to avoid segfault at runtime when compiling with fortran option -check all
358 
359  ios = init_array_int(ig_quadp_gnode_glonum,1,ig_quad_nnode,"ig_quadp_gnode_glonum")
360 
361  ios = init_array_int(ig_quadp_gll_glonum,1,ig_ngll,ig_ngll,"ig_quadp_gll_glonum")
362 
363  ios = init_array_int(ig_quadp_neighbor_hexa,1,"ig_quadp_neighbor_hexa")
364 
365  ios = init_array_int(ig_quadp_neighbor_hexaface,1,"ig_quadp_neighbor_hexaface")
366 
367  endif
368 
369  if (ig_nquad_fsurf > 0) then
370 
371  ios = init_array_int(ig_quadf_gnode_glonum,ig_nquad_fsurf,ig_quad_nnode,"ig_quadf_gnode_glonum")
372 
373  ios = init_array_int(ig_quadf_gll_glonum,ig_nquad_fsurf,ig_ngll,ig_ngll,"ig_quadf_gll_glonum")
374 
375  else !workaround for C subroutine 'fill_mesh_arrays' when there is no quad free surface to avoid segfault at runtime when compiling with fortran option -check all
376 
377  ios = init_array_int(ig_quadf_gnode_glonum,1,ig_quad_nnode,"ig_quadf_gnode_glonum")
378 
379  ios = init_array_int(ig_quadf_gll_glonum,1,ig_ngll,ig_ngll,"ig_quadf_gll_glonum")
380 
381  endif
382 
383  ios = init_array_int(ig_hexa_material_number,ig_nhexa,"ig_hexa_material_number")
384 
385  !
386  !
387  !********************************************************************************************************************************
388  !initialize mesh and create GLL nodes numbering
389  !********************************************************************************************************************************
390  if (ig_myrank == 0) write(ig_lst_unit,'(" ",/,a)') "Generating global gll nodes numbering..."
391 
392  call fill_mesh_arrays(ig_ncpu &
393  ,ig_myrank &
394  ,ig_ngll_total &
395  ,ig_nneighbor_all_kind &
396  ,cpu_neighbor &
397  ,ig_cpu_neighbor_info &
398  ,ig_hexa_gnode_glonum &
399  ,ig_quadp_gnode_glonum &
400  ,ig_quadf_gnode_glonum &
401  ,ig_hexa_gll_glonum &
402  ,ig_quadp_gll_glonum &
403  ,ig_quadf_gll_glonum &
404  ,ig_quadp_neighbor_hexa &
405  ,ig_quadp_neighbor_hexaface &
406  ,ig_hexa_material_number )
407 
408  if (ig_myrank == 0) write(ig_lst_unit,'(a)') "Done"
409 
410  !
411  !
412  !********************************************************************************************************************************
413  !cpu_neighbor -> tg_cpu_neighbor (to avoid problem of structure alignment, which can be compiler dependant)
414  !********************************************************************************************************************************
415  do icpu = 1,ig_ncpu_neighbor
416  tg_cpu_neighbor(icpu)%icpu = cpu_neighbor(icpu)
417  enddo
418  deallocate(cpu_neighbor)
419 
420 
421  !
422  !
423  !***********************************************************
424  !initialize displacement, velocity and acceleration arrays
425  !***********************************************************
426  ios = init_array_real(rg_gll_displacement,ig_ngll_total,ig_ndof,"rg_gll_displacement")
427 
428  ios = init_array_real(rg_gll_velocity ,ig_ngll_total,ig_ndof,"rg_gll_velocity")
429 
430  ios = init_array_real(rg_gll_acceleration,ig_ngll_total,ig_ndof,"rg_gll_acceleration")
431 
432 
433  !
434  !
435  !******************************
436  !write gll info in *.lst file
437  !******************************
438  write(info,'(a)') "gll nodes"
439  call info_all_cpu(ig_ngll_total,info)
440 
441 
442  !
443  !
444  !*************************************************
445  !debug mode : write array ig_hexa_gll_glonum
446  !*************************************************
447  if (lg_output_debug_file) then
448  open(unit=get_newunit(myunit_debug),file="debug."//trim(cg_prefix)//".global.gll."//trim(cg_myrank))
449  do ihexa = 1,ig_nhexa
450  write(unit=myunit_debug,fmt='(a,i10)') "hexa ",ihexa
451  do igll = 1,ig_ngll
452  write(unit=myunit_debug,fmt='(a,i10)') "igll ",igll
453  do jgll = 1,ig_ngll
454  write(unit=myunit_debug,fmt='(10I10)') (ig_hexa_gll_glonum(kgll,jgll,igll,ihexa),kgll=1,ig_ngll)
455  enddo
456  enddo
457  enddo
458  close(myunit_debug)
459  endif
460 
461  !
462  !
463  !*************************************************
464  !debug mode : write free surface geom nodes
465  !*************************************************
466  if (lg_output_debug_file) then
467  open(unit=get_newunit(myunit_debug),file="debug."//trim(cg_prefix)//".freesurface.geom."//trim(cg_myrank))
468  do ifsurf = 1,ig_nquad_fsurf
469  write(unit=myunit_debug,fmt='(a,I10)') "Quad ",ifsurf
470  do inode = 1,ig_quad_nnode
471  write(unit=myunit_debug,fmt='(2(E14.7,1X))') rg_gnode_x(ig_quadf_gnode_glonum(inode,ifsurf))&
472  ,rg_gnode_y(ig_quadf_gnode_glonum(inode,ifsurf))
473  enddo
474  enddo
475  close(myunit_debug)
476  endif
477 
478  if (ig_nquad_parax == 0) then
479 
480  deallocate(ig_quadp_gnode_glonum)
481 
482  deallocate(ig_quadp_gll_glonum)
483 
484  deallocate(ig_quadp_neighbor_hexa)
485 
486  deallocate(ig_quadp_neighbor_hexaface)
487 
488  endif
489 
490  if (ig_nquad_fsurf == 0) then
491 
492  deallocate(ig_quadf_gnode_glonum)
493 
494  deallocate(ig_quadf_gll_glonum)
495 
496  endif
497 
498  return
499 !***********************************************************************************************************************************************************************************
500  end subroutine init_mesh
501 !***********************************************************************************************************************************************************************************
502 
503 !
504 !
509 !***********************************************************************************************************************************************************************************
510  subroutine init_gll_number(ihexa,ngll_total)
511 !***********************************************************************************************************************************************************************************
512 
513  use mpi
514 
515  use mod_global_variables, only : &
516  ig_hexa_gll_glonum&
517  ,ig_ngll
518 
519  implicit none
520 
521  integer, intent(in ) :: ihexa
522  integer, intent(inout) :: ngll_total
523 
524  integer :: k
525  integer :: l
526  integer :: m
527 
528  do k = 1,ig_ngll !along zeta
529  do l = 1,ig_ngll !along eta
530  do m = 1,ig_ngll !along xi
531 
532  if (ig_hexa_gll_glonum(m,l,k,ihexa) == 0) then !GLL has not been numbered yet
533 
534  ngll_total = ngll_total + 1
535  ig_hexa_gll_glonum(m,l,k,ihexa) = ngll_total
536 
537  endif
538 
539  enddo
540  enddo
541  enddo
542 
543  return
544 !***********************************************************************************************************************************************************************************
545  end subroutine init_gll_number
546 !***********************************************************************************************************************************************************************************
547 
548 !
549 !
556 !***********************************************************************************************************************************************************************************
557  subroutine propagate_gll_nodes_quad(ihexa,iface,iquad,global_gll_of_quad,number_of_quad)
558 !***********************************************************************************************************************************************************************************
559 
560  use mpi
561 
562  use mod_global_variables, only : ig_ngll&
563  ,error_stop&
564  ,ig_hexa_gll_glonum
565 
566  implicit none
567 
568  integer, intent(in) :: ihexa
569  integer, intent(in) :: iface
570  integer, intent(in) :: iquad
571  integer, intent(in) :: number_of_quad
572  integer, intent(inout), dimension(IG_NGLL,IG_NGLL,number_of_quad) :: global_gll_of_quad
573 
574  integer :: k,l,m
575 
576  character(len=255) :: info
577 
578  select case(iface)
579  case(1)
580  do l = 1,ig_ngll
581  do m = 1,ig_ngll
582  global_gll_of_quad(m,l,iquad) = ig_hexa_gll_glonum(m,l,1,ihexa)
583  enddo
584  enddo
585  return
586 
587  case(2)
588  do k = 1,ig_ngll
589  do m = 1,ig_ngll
590  global_gll_of_quad(k,m,iquad) = ig_hexa_gll_glonum(m,1,k,ihexa)
591  enddo
592  enddo
593  return
594 
595  case(3)
596  do k = 1,ig_ngll
597  do l = 1,ig_ngll
598  global_gll_of_quad(k,l,iquad) = ig_hexa_gll_glonum(ig_ngll,l,k,ihexa)
599  enddo
600  enddo
601  return
602 
603  case(4)
604  do k = 1,ig_ngll
605  do m = 1,ig_ngll
606  global_gll_of_quad(m,k,iquad) = ig_hexa_gll_glonum(m,ig_ngll,k,ihexa)
607  enddo
608  enddo
609  return
610 
611  case(5)
612  do k = 1,ig_ngll
613  do l = 1,ig_ngll
614  global_gll_of_quad(l,k,iquad) = ig_hexa_gll_glonum(1,l,k,ihexa)
615  enddo
616  enddo
617  return
618 
619  case(6)
620  do l = 1,ig_ngll
621  do m = 1,ig_ngll
622  global_gll_of_quad(l,m,iquad) = ig_hexa_gll_glonum(m,l,ig_ngll,ihexa)
623  enddo
624  enddo
625  return
626 
627  case default
628  write(info,'(a)') "error in subroutine propagate_gll_nodes_quad. Invalid face number"
629  call error_stop(info)
630 
631  end select
632 
633 !***********************************************************************************************************************************************************************************
634  end subroutine propagate_gll_nodes_quad
635 !***********************************************************************************************************************************************************************************
636 
637 !
638 !
645 !***********************************************************************************************************************************************************************************
646  subroutine propagate_gll_nodes_face(ihexa_old,iface_old,ihexa_new,iface_new,icoty_new)
647 !***********************************************************************************************************************************************************************************
648 
649  use mpi
650 
651  use mod_global_variables, only : ig_ngll&
652  ,ig_hexa_gll_glonum
653 
654  implicit none
655 
656  integer, intent(in) :: ihexa_new
657  integer, intent(in) :: ihexa_old
658  integer, intent(in) :: iface_new
659  integer, intent(in) :: iface_old
660  integer, intent(in) :: icoty_new
661 
662  integer :: ext_vec
663  integer :: ext_sgn
664  integer :: int_vec
665  integer :: int_sgn
666  integer :: fixed_val
667  integer :: i_ext
668  integer :: i_int
669  integer :: i_ext_rev
670  integer :: i_int_rev
671  integer :: i_ext_target
672  integer :: i_int_target
673  integer :: igll_source
674  integer :: mid_gll
675 
676 
677 ! icoty_new :
678 ! 0 . 0 .SX.VX.VX.SI.VI.VI
679 ! 7 6 5 4 3 2 1 0 th bit
680 ! 128 64 32 16 8 4 2 1
681 
682  ext_vec = ishft(iand(24, icoty_new), -3)
683  ext_sgn = ishft(iand(32, icoty_new), -5)
684 
685  int_vec = iand(3, icoty_new)
686  int_sgn = ishft(iand(4, icoty_new), -2)
687 
688  mid_gll = ceiling(ig_ngll/2.0)
689  fixed_val = 0
690  select case(iface_new)
691  case(1)
692  if (ig_hexa_gll_glonum(mid_gll,mid_gll,1,ihexa_new) /= 0) return
693  fixed_val = 1
694  case(2)
695  if (ig_hexa_gll_glonum(mid_gll,1,mid_gll,ihexa_new) /= 0) return
696  fixed_val = 1
697  case(3)
698  if (ig_hexa_gll_glonum(ig_ngll,mid_gll,mid_gll,ihexa_new) /= 0) return
699  fixed_val = ig_ngll
700  case(4)
701  if (ig_hexa_gll_glonum(mid_gll,ig_ngll,mid_gll,ihexa_new) /= 0) return
702  fixed_val = ig_ngll
703  case(5)
704  if (ig_hexa_gll_glonum(1,mid_gll,mid_gll,ihexa_new) /= 0) return
705  fixed_val = 1
706  case(6)
707  if (ig_hexa_gll_glonum(mid_gll,mid_gll,ig_ngll,ihexa_new) /= 0) return
708  fixed_val = ig_ngll
709  end select
710 
711  do i_ext = 1,ig_ngll
712 
713  i_ext_rev = (ig_ngll+1)-i_ext
714 
715  do i_int = 1,ig_ngll
716  i_int_rev = (ig_ngll+1)-i_int
717  select case(iface_old)
718  case(1)
719  igll_source = ig_hexa_gll_glonum(i_ext, i_int, 1, ihexa_old)
720  case(2)
721  igll_source = ig_hexa_gll_glonum(i_int, 1, i_ext, ihexa_old)
722  case(3)
723  igll_source = ig_hexa_gll_glonum(ig_ngll, i_int_rev, i_ext_rev, ihexa_old)
724  case(4)
725  igll_source = ig_hexa_gll_glonum(i_ext_rev, ig_ngll, i_int_rev, ihexa_old)
726  case(5)
727  igll_source = ig_hexa_gll_glonum(1, i_ext, i_int, ihexa_old)
728  case(6)
729  igll_source = ig_hexa_gll_glonum(i_int_rev, i_ext_rev, ig_ngll, ihexa_old)
730  end select
731 
732  if (ext_sgn == 1) then
733  i_ext_target = i_ext
734  else
735  i_ext_target = i_ext_rev
736  endif
737  if (int_sgn == 1) then
738  i_int_target = i_int
739  else
740  i_int_target = i_int_rev
741  endif
742 
743  select case(ext_vec)
744  case(1)
745  select case(int_vec)
746  case(2)
747  if (ig_hexa_gll_glonum(fixed_val, i_int_target, i_ext_target, ihexa_new) /= igll_source&
748  .and. ig_hexa_gll_glonum(fixed_val, i_int_target, i_ext_target, ihexa_new) /= 0) write(*,*) 'propagate_gll_nodes_face() ', ihexa_old,iface_old,ihexa_new,iface_new
749  ig_hexa_gll_glonum(fixed_val, i_int_target, i_ext_target, ihexa_new) = igll_source
750  case(3)
751  if (ig_hexa_gll_glonum(i_int_target, fixed_val, i_ext_target, ihexa_new) /= igll_source&
752  .and. ig_hexa_gll_glonum(i_int_target, fixed_val, i_ext_target, ihexa_new) /= 0) write(*,*) 'propagate_gll_nodes_face() ', ihexa_old,iface_old,ihexa_new,iface_new
753  ig_hexa_gll_glonum(i_int_target, fixed_val, i_ext_target, ihexa_new) = igll_source
754  end select
755  case(2)
756  select case(int_vec)
757  case(1)
758  if (ig_hexa_gll_glonum(fixed_val, i_ext_target, i_int_target, ihexa_new) /= igll_source&
759  .and. ig_hexa_gll_glonum(fixed_val, i_ext_target, i_int_target, ihexa_new) /= 0) write(*,*) 'propagate_gll_nodes_face() ', ihexa_old,iface_old,ihexa_new,iface_new
760  ig_hexa_gll_glonum(fixed_val, i_ext_target, i_int_target, ihexa_new) = igll_source
761  case(3)
762  if (ig_hexa_gll_glonum(i_int_target, i_ext_target, fixed_val, ihexa_new) /= igll_source&
763  .and. ig_hexa_gll_glonum(i_int_target, i_ext_target, fixed_val, ihexa_new) /= 0) write(*,*) 'propagate_gll_nodes_face() ', ihexa_old,iface_old,ihexa_new,iface_new
764  ig_hexa_gll_glonum(i_int_target, i_ext_target, fixed_val, ihexa_new) = igll_source
765  end select
766  case(3)
767  select case(int_vec)
768  case(1)
769  if (ig_hexa_gll_glonum(i_ext_target, fixed_val, i_int_target, ihexa_new) /= igll_source&
770  .and. ig_hexa_gll_glonum(i_ext_target, fixed_val, i_int_target, ihexa_new) /= 0) write(*,*) 'propagate_gll_nodes_face() ', ihexa_old,iface_old,ihexa_new,iface_new
771  ig_hexa_gll_glonum(i_ext_target, fixed_val, i_int_target, ihexa_new) = igll_source
772  case(2)
773  if (ig_hexa_gll_glonum(i_ext_target, i_int_target, fixed_val, ihexa_new) /= igll_source&
774  .and. ig_hexa_gll_glonum(i_ext_target, i_int_target, fixed_val, ihexa_new) /= 0) write(*,*) 'propagate_gll_nodes_face() ', ihexa_old,iface_old,ihexa_new,iface_new
775  ig_hexa_gll_glonum(i_ext_target, i_int_target, fixed_val, ihexa_new) = igll_source
776  end select
777  end select
778  enddo
779  enddo
780 
781  return
782 !***********************************************************************************************************************************************************************************
783  end subroutine propagate_gll_nodes_face
784 !***********************************************************************************************************************************************************************************
785 
786 !
787 !
794 !***********************************************************************************************************************************************************************************
795  subroutine propagate_gll_nodes_edge(ihexa_old,iedge_old,ihexa_new,iedge_new,icoty_new)
796 !***********************************************************************************************************************************************************************************
797 
798  use mpi
799 
800  use mod_global_variables, only : ig_ngll&
801  ,ig_myrank&
802  ,ig_hexa_gll_glonum&
803  ,lg_output_debug_file
804 
805  implicit none
806 
807  integer, intent(in) :: ihexa_new
808  integer, intent(in) :: ihexa_old
809  integer, intent(in) :: iedge_new
810  integer, intent(in) :: iedge_old
811  integer, intent(in) :: icoty_new
812 
813  integer, pointer :: pedge_old(:)
814  integer, pointer :: pedge_new(:)
815  integer :: i, j
816 
817  select case(iedge_old)
818  case(1)
819  pedge_old => ig_hexa_gll_glonum(1:ig_ngll,1,1,ihexa_old)
820  case(2)
821  pedge_old => ig_hexa_gll_glonum(ig_ngll,1:ig_ngll,1,ihexa_old)
822  case(3)
823  pedge_old => ig_hexa_gll_glonum(1:ig_ngll,ig_ngll,1,ihexa_old)
824  case(4)
825  pedge_old => ig_hexa_gll_glonum(1,1:ig_ngll,1,ihexa_old)
826  case(5)
827  pedge_old => ig_hexa_gll_glonum(1:ig_ngll,1,ig_ngll,ihexa_old)
828  case(6)
829  pedge_old => ig_hexa_gll_glonum(ig_ngll,1:ig_ngll,ig_ngll,ihexa_old)
830  case(7)
831  pedge_old => ig_hexa_gll_glonum(1:ig_ngll,ig_ngll,ig_ngll,ihexa_old)
832  case(8)
833  pedge_old => ig_hexa_gll_glonum(1,1:ig_ngll,ig_ngll,ihexa_old)
834  case(9)
835  pedge_old => ig_hexa_gll_glonum(1,1,1:ig_ngll,ihexa_old)
836  case(10)
837  pedge_old => ig_hexa_gll_glonum(ig_ngll,1,1:ig_ngll,ihexa_old)
838  case(11)
839  pedge_old => ig_hexa_gll_glonum(ig_ngll,ig_ngll,1:ig_ngll,ihexa_old)
840  case(12)
841  pedge_old => ig_hexa_gll_glonum(1,ig_ngll,1:ig_ngll,ihexa_old)
842  end select
843 
844  select case(iedge_new)
845  case(1)
846  pedge_new => ig_hexa_gll_glonum(1:ig_ngll,1,1,ihexa_new)
847  case(2)
848  pedge_new => ig_hexa_gll_glonum(ig_ngll,1:ig_ngll,1,ihexa_new)
849  case(3)
850  pedge_new => ig_hexa_gll_glonum(1:ig_ngll,ig_ngll,1,ihexa_new)
851  case(4)
852  pedge_new => ig_hexa_gll_glonum(1,1:ig_ngll,1,ihexa_new)
853  case(5)
854  pedge_new => ig_hexa_gll_glonum(1:ig_ngll,1,ig_ngll,ihexa_new)
855  case(6)
856  pedge_new => ig_hexa_gll_glonum(ig_ngll,1:ig_ngll,ig_ngll,ihexa_new)
857  case(7)
858  pedge_new => ig_hexa_gll_glonum(1:ig_ngll,ig_ngll,ig_ngll,ihexa_new)
859  case(8)
860  pedge_new => ig_hexa_gll_glonum(1,1:ig_ngll,ig_ngll,ihexa_new)
861  case(9)
862  pedge_new => ig_hexa_gll_glonum(1,1,1:ig_ngll,ihexa_new)
863  case(10)
864  pedge_new => ig_hexa_gll_glonum(ig_ngll,1,1:ig_ngll,ihexa_new)
865  case(11)
866  pedge_new => ig_hexa_gll_glonum(ig_ngll,ig_ngll,1:ig_ngll,ihexa_new)
867  case(12)
868  pedge_new => ig_hexa_gll_glonum(1,ig_ngll,1:ig_ngll,ihexa_new)
869  end select
870 
871  do i = 1,ig_ngll
872  if (icoty_new == 1) then
873  j = i
874  else
875  j = ig_ngll+1-i
876  endif
877  if ( (pedge_new(j) /= pedge_old(i) .and. pedge_new(j) /= 0) .or. pedge_old(i) == 0 ) write(*,*) 'propagate_gll_nodes_edge',ihexa_old,iedge_old,ihexa_new,iedge_new
878  pedge_new(j) = pedge_old(i)
879  enddo
880 
881  return
882 !***********************************************************************************************************************************************************************************
883  end subroutine propagate_gll_nodes_edge
884 !***********************************************************************************************************************************************************************************
885 
886 !
887 !
893 !***********************************************************************************************************************************************************************************
894  subroutine propagate_gll_nodes_corner(ihexa_old,icorner_old,ihexa_new,icorner_new)
895 !***********************************************************************************************************************************************************************************
896 
897  use mpi
898 
899  use mod_global_variables, only : ig_ngll&
900  ,ig_myrank&
901  ,ig_hexa_gll_glonum&
902  ,lg_output_debug_file
903 
904  implicit none
905 
906  integer, intent(in) :: ihexa_new
907  integer, intent(in) :: ihexa_old
908  integer, intent(in) :: icorner_new
909  integer, intent(in) :: icorner_old
910 
911  integer, pointer :: pcorner_old
912  integer, pointer :: pcorner_new
913 
914  select case(icorner_old)
915  case(1)
916  pcorner_old => ig_hexa_gll_glonum(1,1,1,ihexa_old)
917  case(2)
918  pcorner_old => ig_hexa_gll_glonum(ig_ngll,1,1,ihexa_old)
919  case(3)
920  pcorner_old => ig_hexa_gll_glonum(ig_ngll,ig_ngll,1,ihexa_old)
921  case(4)
922  pcorner_old => ig_hexa_gll_glonum(1,ig_ngll,1,ihexa_old)
923  case(5)
924  pcorner_old => ig_hexa_gll_glonum(1,1,ig_ngll,ihexa_old)
925  case(6)
926  pcorner_old => ig_hexa_gll_glonum(ig_ngll,1,ig_ngll,ihexa_old)
927  case(7)
928  pcorner_old => ig_hexa_gll_glonum(ig_ngll,ig_ngll,ig_ngll,ihexa_old)
929  case(8)
930  pcorner_old => ig_hexa_gll_glonum(1,ig_ngll,ig_ngll,ihexa_old)
931  end select
932 
933  select case(icorner_new)
934  case(1)
935  pcorner_new => ig_hexa_gll_glonum(1,1,1,ihexa_new)
936  case(2)
937  pcorner_new => ig_hexa_gll_glonum(ig_ngll,1,1,ihexa_new)
938  case(3)
939  pcorner_new => ig_hexa_gll_glonum(ig_ngll,ig_ngll,1,ihexa_new)
940  case(4)
941  pcorner_new => ig_hexa_gll_glonum(1,ig_ngll,1,ihexa_new)
942  case(5)
943  pcorner_new => ig_hexa_gll_glonum(1,1,ig_ngll,ihexa_new)
944  case(6)
945  pcorner_new => ig_hexa_gll_glonum(ig_ngll,1,ig_ngll,ihexa_new)
946  case(7)
947  pcorner_new => ig_hexa_gll_glonum(ig_ngll,ig_ngll,ig_ngll,ihexa_new)
948  case(8)
949  pcorner_new => ig_hexa_gll_glonum(1,ig_ngll,ig_ngll,ihexa_new)
950  end select
951  if ( (pcorner_new /= pcorner_old .and. pcorner_new /= 0) .or. pcorner_old == 0 ) write(*,*) 'propagate_gll_nodes_corner',ihexa_old,icorner_old,ihexa_new,icorner_new
952  pcorner_new = pcorner_old
953 
954  return
955 !***********************************************************************************************************************************************************************************
956  end subroutine propagate_gll_nodes_corner
957 !***********************************************************************************************************************************************************************************
958 
959 !
960 !
972 !***********************************************************************************************************************************************************************************
973  subroutine init_element(ilnnhe,ilnnqu)
974 !***********************************************************************************************************************************************************************************
975 
976  use mpi
977 
978  use mod_global_variables, only : ig_hexa_gnode_xiloc&
979  ,ig_hexa_gnode_etloc&
980  ,ig_hexa_gnode_zeloc&
981  ,ig_line_nnode&
982  ,ig_hexa_node2gll&
983  ,rg_gnode_abscissa&
984  ,ig_quad_gnode_xiloc&
985  ,ig_quad_gnode_etloc&
986  ,rg_gnode_abscissa_dist&
987  ,ig_ngll&
988  ,error_stop
989 
990  use mod_init_memory
991 
992  implicit none
993 
994  integer, intent(in) :: ilnnhe
995  integer, intent(in) :: ilnnqu
996 
997  integer :: ios
998  integer :: i
999  integer :: j
1000 
1001  character(len=255) :: info
1002  !
1003  !
1004  !**********************************************************************************
1005  !fill local position of geometric node of hexa8 or hexa27 on xi, eta and zeta axis
1006  !**********************************************************************************
1007  ios = init_array_int(ig_hexa_gnode_xiloc,ilnnhe,"ig_hexa_gnode_xiloc")
1008 
1009  ios = init_array_int(ig_hexa_gnode_etloc,ilnnhe,"ig_hexa_gnode_etloc")
1010 
1011  ios = init_array_int(ig_hexa_gnode_zeloc,ilnnhe,"ig_hexa_gnode_zeloc")
1012 
1013  if (ilnnhe == 8) then
1014 
1015  ig_line_nnode = 2
1016 !
1017 !------->fill rg_gnode_abscissa: local coordinate of geometric nodes
1018  ios = init_array_real(rg_gnode_abscissa,ig_line_nnode,"rg_gnode_abscissa")
1019 
1020  rg_gnode_abscissa(1) = +1.0
1021  rg_gnode_abscissa(2) = -1.0
1022 !
1023 !------->local position of geometric nodes
1024  ig_hexa_gnode_xiloc(1) = 1
1025  ig_hexa_gnode_etloc(1) = 1
1026  ig_hexa_gnode_zeloc(1) = 1
1027  ig_hexa_gnode_xiloc(2) = 2
1028  ig_hexa_gnode_etloc(2) = 1
1029  ig_hexa_gnode_zeloc(2) = 1
1030  ig_hexa_gnode_xiloc(3) = 2
1031  ig_hexa_gnode_etloc(3) = 2
1032  ig_hexa_gnode_zeloc(3) = 1
1033  ig_hexa_gnode_xiloc(4) = 1
1034  ig_hexa_gnode_etloc(4) = 2
1035  ig_hexa_gnode_zeloc(4) = 1
1036  ig_hexa_gnode_xiloc(5) = 1
1037  ig_hexa_gnode_etloc(5) = 1
1038  ig_hexa_gnode_zeloc(5) = 2
1039  ig_hexa_gnode_xiloc(6) = 2
1040  ig_hexa_gnode_etloc(6) = 1
1041  ig_hexa_gnode_zeloc(6) = 2
1042  ig_hexa_gnode_xiloc(7) = 2
1043  ig_hexa_gnode_etloc(7) = 2
1044  ig_hexa_gnode_zeloc(7) = 2
1045  ig_hexa_gnode_xiloc(8) = 1
1046  ig_hexa_gnode_etloc(8) = 2
1047  ig_hexa_gnode_zeloc(8) = 2
1048 
1049  elseif (ilnnhe == 27) then
1050 
1051  ig_line_nnode = 3
1052 !
1053 !------->fill rg_gnode_abscissa: local coordinate of geometric nodes
1054  ios = init_array_real(rg_gnode_abscissa,ig_line_nnode,"rg_gnode_abscissa")
1055 
1056  rg_gnode_abscissa(1) = +1.0
1057  rg_gnode_abscissa(2) = 0.0
1058  rg_gnode_abscissa(3) = -1.0
1059 !
1060 !------->local position of geometric nodes
1061  ig_hexa_gnode_xiloc( 1) = 1
1062  ig_hexa_gnode_etloc( 1) = 1
1063  ig_hexa_gnode_zeloc( 1) = 1
1064  ig_hexa_gnode_xiloc( 2) = 3
1065  ig_hexa_gnode_etloc( 2) = 1
1066  ig_hexa_gnode_zeloc( 2) = 1
1067  ig_hexa_gnode_xiloc( 3) = 3
1068  ig_hexa_gnode_etloc( 3) = 3
1069  ig_hexa_gnode_zeloc( 3) = 1
1070  ig_hexa_gnode_xiloc( 4) = 1
1071  ig_hexa_gnode_etloc( 4) = 3
1072  ig_hexa_gnode_zeloc( 4) = 1
1073  ig_hexa_gnode_xiloc( 5) = 1
1074  ig_hexa_gnode_etloc( 5) = 1
1075  ig_hexa_gnode_zeloc( 5) = 3
1076  ig_hexa_gnode_xiloc( 6) = 3
1077  ig_hexa_gnode_etloc( 6) = 1
1078  ig_hexa_gnode_zeloc( 6) = 3
1079  ig_hexa_gnode_xiloc( 7) = 3
1080  ig_hexa_gnode_etloc( 7) = 3
1081  ig_hexa_gnode_zeloc( 7) = 3
1082  ig_hexa_gnode_xiloc( 8) = 1
1083  ig_hexa_gnode_etloc( 8) = 3
1084  ig_hexa_gnode_zeloc( 8) = 3
1085  ig_hexa_gnode_xiloc( 9) = 2
1086  ig_hexa_gnode_etloc( 9) = 1
1087  ig_hexa_gnode_zeloc( 9) = 1
1088  ig_hexa_gnode_xiloc(10) = 3
1089  ig_hexa_gnode_etloc(10) = 2
1090  ig_hexa_gnode_zeloc(10) = 1
1091  ig_hexa_gnode_xiloc(11) = 2
1092  ig_hexa_gnode_etloc(11) = 3
1093  ig_hexa_gnode_zeloc(11) = 1
1094  ig_hexa_gnode_xiloc(12) = 1
1095  ig_hexa_gnode_etloc(12) = 2
1096  ig_hexa_gnode_zeloc(12) = 1
1097  ig_hexa_gnode_xiloc(13) = 2
1098  ig_hexa_gnode_etloc(13) = 1
1099  ig_hexa_gnode_zeloc(13) = 3
1100  ig_hexa_gnode_xiloc(14) = 3
1101  ig_hexa_gnode_etloc(14) = 2
1102  ig_hexa_gnode_zeloc(14) = 3
1103  ig_hexa_gnode_xiloc(15) = 2
1104  ig_hexa_gnode_etloc(15) = 3
1105  ig_hexa_gnode_zeloc(15) = 3
1106  ig_hexa_gnode_xiloc(16) = 1
1107  ig_hexa_gnode_etloc(16) = 2
1108  ig_hexa_gnode_zeloc(16) = 3
1109  ig_hexa_gnode_xiloc(17) = 1
1110  ig_hexa_gnode_etloc(17) = 1
1111  ig_hexa_gnode_zeloc(17) = 2
1112  ig_hexa_gnode_xiloc(18) = 3
1113  ig_hexa_gnode_etloc(18) = 1
1114  ig_hexa_gnode_zeloc(18) = 2
1115  ig_hexa_gnode_xiloc(19) = 3
1116  ig_hexa_gnode_etloc(19) = 3
1117  ig_hexa_gnode_zeloc(19) = 2
1118  ig_hexa_gnode_xiloc(20) = 1
1119  ig_hexa_gnode_etloc(20) = 3
1120  ig_hexa_gnode_zeloc(20) = 2
1121  ig_hexa_gnode_xiloc(21) = 2
1122  ig_hexa_gnode_etloc(21) = 2
1123  ig_hexa_gnode_zeloc(21) = 1
1124  ig_hexa_gnode_xiloc(22) = 2
1125  ig_hexa_gnode_etloc(22) = 1
1126  ig_hexa_gnode_zeloc(22) = 2
1127  ig_hexa_gnode_xiloc(23) = 3
1128  ig_hexa_gnode_etloc(23) = 2
1129  ig_hexa_gnode_zeloc(23) = 2
1130  ig_hexa_gnode_xiloc(24) = 2
1131  ig_hexa_gnode_etloc(24) = 3
1132  ig_hexa_gnode_zeloc(24) = 2
1133  ig_hexa_gnode_xiloc(25) = 1
1134  ig_hexa_gnode_etloc(25) = 2
1135  ig_hexa_gnode_zeloc(25) = 2
1136  ig_hexa_gnode_xiloc(26) = 2
1137  ig_hexa_gnode_etloc(26) = 2
1138  ig_hexa_gnode_zeloc(26) = 3
1139  ig_hexa_gnode_xiloc(27) = 2
1140  ig_hexa_gnode_etloc(27) = 2
1141  ig_hexa_gnode_zeloc(27) = 2
1142  else
1143  write(info,'(a)') "error in init_element: invalid number of geometrical nodes for hexa"
1144  call error_stop(info)
1145  endif
1146 
1147 
1148  !
1149  !
1150  !******************************************************************
1151  !init rg_gnode_abscissa_dist
1152  !******************************************************************
1153  ios = init_array_real(rg_gnode_abscissa_dist,ig_line_nnode,ig_line_nnode,"rg_gnode_abscissa_dist")
1154 
1155  do i = 1,ig_line_nnode
1156  do j = 1,ig_line_nnode
1157  rg_gnode_abscissa_dist(j,i) = 0.0
1158  if (i /= j) rg_gnode_abscissa_dist(j,i) = 1.0/(rg_gnode_abscissa(i) - rg_gnode_abscissa(j))
1159  enddo
1160  enddo
1161 
1162  !
1163  !
1164  !******************************************************************
1165  !init local position of geometric node of quad4 or quad9 on xi, eta
1166  !******************************************************************
1167  ios = init_array_int(ig_quad_gnode_xiloc,ilnnqu,"ig_quad_gnode_xiloc")
1168 
1169  ios = init_array_int(ig_quad_gnode_etloc,ilnnqu,"ig_quad_gnode_etloc")
1170 
1171  if (ilnnqu == 4) then
1172  ig_quad_gnode_xiloc(1) = 1
1173  ig_quad_gnode_etloc(1) = 1
1174  ig_quad_gnode_xiloc(2) = 2
1175  ig_quad_gnode_etloc(2) = 1
1176  ig_quad_gnode_xiloc(3) = 2
1177  ig_quad_gnode_etloc(3) = 2
1178  ig_quad_gnode_xiloc(4) = 1
1179  ig_quad_gnode_etloc(4) = 2
1180  elseif (ilnnqu == 9) then
1181  ig_quad_gnode_xiloc(1) = 1
1182  ig_quad_gnode_etloc(1) = 1
1183  ig_quad_gnode_xiloc(2) = 3
1184  ig_quad_gnode_etloc(2) = 1
1185  ig_quad_gnode_xiloc(3) = 3
1186  ig_quad_gnode_etloc(3) = 3
1187  ig_quad_gnode_xiloc(4) = 1
1188  ig_quad_gnode_etloc(4) = 3
1189  ig_quad_gnode_xiloc(5) = 2
1190  ig_quad_gnode_etloc(5) = 1
1191  ig_quad_gnode_xiloc(6) = 3
1192  ig_quad_gnode_etloc(6) = 2
1193  ig_quad_gnode_xiloc(7) = 2
1194  ig_quad_gnode_etloc(7) = 3
1195  ig_quad_gnode_xiloc(8) = 1
1196  ig_quad_gnode_etloc(8) = 2
1197  ig_quad_gnode_xiloc(9) = 2
1198  ig_quad_gnode_etloc(9) = 2
1199  else
1200  write(info,'(a)') "error in init_element: invalid number of geometrical nodes for quad"
1201  call error_stop(info)
1202  endif
1203  !
1204  !
1205  !***********************************************************************
1206  !indirection from geom nodes to GLL nodes (see docs/Convention.html)
1207  !***********************************************************************
1208  ! xi eta zeta
1209  ! m - l - k
1210  !
1211  !geom node
1212  ! 1 1 1 1
1213  ! 2 IG_NGLL 1 1
1214  ! 3 IG_NGLL IG_NGLL 1
1215  ! 4 1 IG_NGLL 1
1216  ! 5 1 1 IG_NGLL
1217  ! 6 IG_NGLL 1 IG_NGLL
1218  ! 7 IG_NGLL IG_NGLL IG_NGLL
1219  ! 8 1 IG_NGLL IG_NGLL
1220  !
1221 
1222  ig_hexa_node2gll(1,1) = 1
1223  ig_hexa_node2gll(2,1) = 1
1224  ig_hexa_node2gll(3,1) = 1
1225 
1226  ig_hexa_node2gll(1,2) = 1
1227  ig_hexa_node2gll(2,2) = 1
1228  ig_hexa_node2gll(3,2) = ig_ngll
1229 
1230  ig_hexa_node2gll(1,3) = 1
1231  ig_hexa_node2gll(2,3) = ig_ngll
1232  ig_hexa_node2gll(3,3) = ig_ngll
1233 
1234  ig_hexa_node2gll(1,4) = 1
1235  ig_hexa_node2gll(2,4) = ig_ngll
1236  ig_hexa_node2gll(3,4) = 1
1237 
1238  ig_hexa_node2gll(1,5) = ig_ngll
1239  ig_hexa_node2gll(2,5) = 1
1240  ig_hexa_node2gll(3,5) = 1
1241 
1242  ig_hexa_node2gll(1,6) = ig_ngll
1243  ig_hexa_node2gll(2,6) = 1
1244  ig_hexa_node2gll(3,6) = ig_ngll
1245 
1246  ig_hexa_node2gll(1,7) = ig_ngll
1247  ig_hexa_node2gll(2,7) = ig_ngll
1248  ig_hexa_node2gll(3,7) = ig_ngll
1249 
1250  ig_hexa_node2gll(1,8) = ig_ngll
1251  ig_hexa_node2gll(2,8) = ig_ngll
1252  ig_hexa_node2gll(3,8) = 1
1253 
1254  return
1255 !***********************************************************************************************************************************************************************************
1256  end subroutine init_element
1257 !***********************************************************************************************************************************************************************************
1258 
1259 end module mod_init_mesh
mod_init_mesh
This module contains subroutines to read mesh files and creates GLL nodes global numbering in cpu myr...
Definition: module_init_mesh.f90:126
mod_init_memory
This module contains subroutines to allocate arrays and to compute an approximation of the total RAM ...
Definition: module_init_memory.f90:126
mod_init_mesh::propagate_gll_nodes_corner
subroutine, private propagate_gll_nodes_corner(ihexa_old, icorner_old, ihexa_new, icorner_new)
This subroutine propagates existing GLL numbering of hexahedron elements to corner (i...
Definition: module_init_mesh.f90:894
mod_init_mesh::propagate_gll_nodes_face
subroutine, private propagate_gll_nodes_face(ihexa_old, iface_old, ihexa_new, iface_new, icoty_new)
This subroutine propagates existing GLL numbering of hexahedron elements to the face of its neighbors...
Definition: module_init_mesh.f90:646
mod_global_variables
This module defines all global variables of EFISPEC3D. Scalar variables are initialized directly in t...
Definition: module_global_variables.f90:128
mod_init_memory::init_array_real
Interface init_array_real to redirect allocation to n-rank arrays.
Definition: module_init_memory.f90:148
mod_init_mesh::propagate_gll_nodes_quad
subroutine, private propagate_gll_nodes_quad(ihexa, iface, iquad, global_gll_of_quad, number_of_quad)
This subroutine creates the local to global GLL indirection array for quadrangle elements by propagat...
Definition: module_init_mesh.f90:557
mod_init_mesh::init_element
subroutine, private init_element(ilnnhe, ilnnqu)
subroutine to set up convention of hexahedron and quadrangle elements
Definition: module_init_mesh.f90:973
mod_init_mesh::init_mesh
subroutine, public init_mesh()
This subroutine reads mesh files *.inp for cpu myrank and creates GLL numbering of hexahedron and qua...
Definition: module_init_mesh.f90:177
mod_init_mesh::propagate_gll_nodes_edge
subroutine, private propagate_gll_nodes_edge(ihexa_old, iedge_old, ihexa_new, iedge_new, icoty_new)
This subroutine propagates existing GLL numbering of hexahedron elements to the edge of its neighbors...
Definition: module_init_mesh.f90:795
mod_init_memory::init_array_int
Interface init_array_int to redirect allocation to n-rank arrays.
Definition: module_init_memory.f90:160
mod_init_mesh::init_gll_number
subroutine, private init_gll_number(ihexa, ngll_total)
This subroutine increments GLL numbering of hexahedron elements in cpu myrank (see variable mod_globa...
Definition: module_init_mesh.f90:510