All Classes Files Functions Variables Pages
module_gll_value.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 
127 
128  use mpi
129 
130  implicit none
131 
132  private
133 
134  public :: get_hexa_gll_value
135  public :: get_quad_gll_value
136  public :: get_node_gll_value
137 
138 !
141  module procedure get_node_gll_value_x
142  module procedure get_node_gll_value_xyz
143  end interface
144 
145  contains
146 
147 !
148 !
153 !***********************************************************************************************************************************************************************************
154  subroutine get_hexa_gll_value(gll_all_values,gll_global_number,gll_selected_values)
155 !***********************************************************************************************************************************************************************************
156 
157  use mod_global_variables, only :&
158  ig_ndof&
159  ,ig_ngll&
160  ,ig_ngll_total
161 
162  implicit none
163 
164  real , intent( in), dimension(IG_NDOF,ig_ngll_total) :: gll_all_values
165  integer, intent( in), dimension( IG_NGLL,IG_NGLL,IG_NGLL) :: gll_global_number
166  real , intent(out), dimension(IG_NDOF,IG_NGLL,IG_NGLL,IG_NGLL) :: gll_selected_values
167 
168  integer :: igll
169  integer :: k
170  integer :: l
171  integer :: m
172 
173  do k = 1,ig_ngll
174  do l = 1,ig_ngll
175  do m = 1,ig_ngll
176 
177  igll = gll_global_number(m,l,k)
178 
179  gll_selected_values(1,m,l,k) = gll_all_values(1,igll)
180  gll_selected_values(2,m,l,k) = gll_all_values(2,igll)
181  gll_selected_values(3,m,l,k) = gll_all_values(3,igll)
182 
183  enddo
184  enddo
185  enddo
186 
187  return
188 
189 !***********************************************************************************************************************************************************************************
190  end subroutine get_hexa_gll_value
191 !***********************************************************************************************************************************************************************************
192 
193 !
194 !
199 !***********************************************************************************************************************************************************************************
200  subroutine get_quad_gll_value(gll_all_values,gll_global_number,gll_selected_values)
201 !***********************************************************************************************************************************************************************************
202 
203  use mod_global_variables, only :&
204  ig_ndof&
205  ,ig_ngll&
206  ,ig_ngll_total
207 
208  implicit none
209 
210  real , intent( in), dimension(IG_NDOF,ig_ngll_total) :: gll_all_values
211  integer, intent( in), dimension( IG_NGLL,IG_NGLL) :: gll_global_number
212  real , intent(out), dimension(IG_NDOF,IG_NGLL,IG_NGLL) :: gll_selected_values
213 
214  integer :: igll
215  integer :: k
216  integer :: l
217 
218  do k = 1,ig_ngll
219  do l = 1,ig_ngll
220 
221  igll = gll_global_number(l,k)
222 
223  gll_selected_values(1,l,k) = gll_all_values(1,igll)
224  gll_selected_values(2,l,k) = gll_all_values(2,igll)
225  gll_selected_values(3,l,k) = gll_all_values(3,igll)
226 
227  enddo
228  enddo
229 
230  return
231 
232 !***********************************************************************************************************************************************************************************
233  end subroutine get_quad_gll_value
234 !***********************************************************************************************************************************************************************************
235 
236 !
237 !
242 !***********************************************************************************************************************************************************************************
243  subroutine get_node_gll_value_x(gll_all_values,gll_global_number,gll_selected_values)
244 !***********************************************************************************************************************************************************************************
245 
246  use mod_global_variables, only :&
247  ig_ngll_total
248 
249  implicit none
250 
251  real , intent( in), dimension(ig_ngll_total) :: gll_all_values
252  integer, intent( in), dimension(:) :: gll_global_number
253  real , intent(out), dimension(:) :: gll_selected_values
254 
255  integer :: igll
256  integer :: k
257 
258  do k = 1,size(gll_global_number)
259 
260  igll = gll_global_number(k)
261 
262  gll_selected_values(k) = gll_all_values(igll)
263 
264  enddo
265 
266  return
267 
268 !***********************************************************************************************************************************************************************************
269  end subroutine get_node_gll_value_x
270 !***********************************************************************************************************************************************************************************
271 
272 !
273 !
278 !***********************************************************************************************************************************************************************************
279  subroutine get_node_gll_value_xyz(gll_all_values,gll_global_number,gll_selected_values)
280 !***********************************************************************************************************************************************************************************
281 
282  use mod_global_variables, only :&
283  ig_ndof&
284  ,ig_ngll_total
285 
286  implicit none
287 
288  real , intent( in), dimension(IG_NDOF,ig_ngll_total) :: gll_all_values
289  integer, intent( in), dimension(:) :: gll_global_number
290  real , intent(out), dimension(:,:) :: gll_selected_values
291 
292  integer :: igll
293  integer :: k
294 
295  do k = 1,size(gll_global_number)
296 
297  igll = gll_global_number(k)
298 
299  gll_selected_values(1,k) = gll_all_values(1,igll)
300  gll_selected_values(2,k) = gll_all_values(2,igll)
301  gll_selected_values(3,k) = gll_all_values(3,igll)
302 
303  enddo
304 
305  return
306 
307 !***********************************************************************************************************************************************************************************
308  end subroutine get_node_gll_value_xyz
309 !***********************************************************************************************************************************************************************************
310 
311 end module mod_gll_value
subroutine, public get_hexa_gll_value(gll_all_values, gll_global_number, gll_selected_values)
subroutine to get hexahedron element GLL values from an array ordered by global GLL nodes numbering...
This module defines all global variables of EFISPEC3D. Scalar variables are initialized directly in t...
Interface for redirection to subroutines mod_gll_value::get_node_gll_value_x or mod_gll_value::get_no...
subroutine get_node_gll_value_xyz(gll_all_values, gll_global_number, gll_selected_values)
subroutine to get GLL node values from an array ordered by global GLL nodes numbering.
subroutine get_node_gll_value_x(gll_all_values, gll_global_number, gll_selected_values)
subroutine to get GLL node values from an array ordered by global GLL nodes numbering.
This module contains subroutines to get GLL nodes values from global GLL nodes numbering.
subroutine, public get_quad_gll_value(gll_all_values, gll_global_number, gll_selected_values)
subroutine to get quadrangle element GLL values from an array ordered by global GLL nodes numbering...