MFC:Pre_process  v1.0
m_grid.f90
Go to the documentation of this file.
1 !! __ _______________
2 !! / |/ / ____/ ____/
3 !! / /|_/ / /_ / /
4 !! / / / / __/ / /___
5 !! /_/ /_/_/ \____/
6 !!
7 !! This file is part of MFC.
8 !!
9 !! MFC is the legal property of its developers, whose names
10 !! are listed in the copyright file included with this source
11 !! distribution.
12 !!
13 !! MFC is free software: you can redistribute it and/or modify
14 !! it under the terms of the GNU General Public License as published
15 !! by the Free Software Foundation, either version 3 of the license
16 !! or any later version.
17 !!
18 !! MFC is distributed in the hope that it will be useful,
19 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
20 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21 !! GNU General Public License for more details.
22 !!
23 !! You should have received a copy of the GNU General Public License
24 !! along with MFC (LICENSE).
25 !! If not, see <http://www.gnu.org/licenses/>.
26 
33 
44 MODULE m_grid
45 
46 
47  ! Dependencies =============================================================
48  USE m_derived_types ! Definitions of the derived types
49 
50  USE m_global_parameters ! Global parameters for the code
51 
52  USE m_mpi_proxy ! Message passing interface (MPI) module proxy
53 
54  USE mpi ! Message passing interface (MPI) module
55  ! ==========================================================================
56 
57 
58  IMPLICIT NONE
59 
60 
61  PRIVATE; PUBLIC :: s_initialize_grid_module, &
66 
67  abstract INTERFACE ! ===================================================
68 
69  SUBROUTINE s_generate_abstract_grid(dflt_int) ! ------------------------
70 
71  INTEGER, INTENT(IN) :: dflt_int
72 
73  END SUBROUTINE s_generate_abstract_grid ! ----------------------
74 
75  END INTERFACE ! ========================================================
76 
77  PROCEDURE(s_generate_abstract_grid), POINTER :: s_generate_grid => null()
78 
79  CONTAINS
80 
81 
88  SUBROUTINE s_generate_serial_grid(dflt_int) ! -----------------------------------------
89 
90  INTEGER, INTENT(IN) :: dflt_int
91 
92  ! Generic loop iterator
93  INTEGER :: i,j
94  REAL(KIND(0d0)) :: length
95 
96  ! Grid Generation in the x-direction ===============================
97  dx = (x_domain%end - x_domain%beg) / REAL(m+1, kind(0d0))
98 
99  DO i = 0, m
100  x_cc( i ) = x_domain%beg + 5d-1*dx*REAL(2*i+1, kind(0d0))
101  x_cb(i-1) = x_domain%beg + dx*REAL( i , kind(0d0))
102  END DO
103 
104  x_cb(m) = x_domain%end
105 
106  IF(stretch_x) THEN
107 
108  length = abs(x_cb(m)-x_cb(0))
109  x_cb = x_cb / length
110  x_a = x_a / length
111  x_b = x_b / length
112 
113  DO j = 1, loops_x
114  DO i = -1, m
115  x_cb(i) = x_cb(i) / a_x * &
116  ( a_x + log(cosh( a_x*(x_cb(i) - x_a) )) &
117  + log(cosh( a_x*(x_cb(i) - x_b) )) &
118  - 2d0*log(cosh( 0.5d0*a_x*(x_b - x_a) )) )
119  END DO
120  END DO
121  x_cb = x_cb * length
122 
123  x_cc = (x_cb(0:m) + x_cb(-1:m-1))/2d0
124 
125  dx = minval(x_cb(0:m) - x_cb(-1:m-1))
126  print*, 'Stretched grid: min/max x grid: ', minval(x_cc(:)), maxval(x_cc(:))
127  IF(num_procs > 1) CALL s_mpi_reduce_min(dx)
128 
129  END IF
130  ! ==================================================================
131 
132 
133  ! Grid Generation in the y-direction ===============================
134  IF(n == 0) RETURN
135 
136  IF (grid_geometry == 2.AND.y_domain%beg.eq.0.0d0) THEN
137  !IF (grid_geometry == 2) THEN
138 
139  dy = (y_domain%end - y_domain%beg) / REAL(2*n+1, kind(0d0))
140 
141  y_cc(0) = y_domain%beg + 5d-1*dy
142  y_cb(-1) = y_domain%beg
143 
144  DO i = 1, n
145  y_cc( i ) = y_domain%beg + 2d0*dy*REAL( i , kind(0d0))
146  y_cb(i-1) = y_domain%beg + dy*REAL(2*i-1, kind(0d0))
147  END DO
148 
149  ELSE
150 
151  dy = (y_domain%end - y_domain%beg) / REAL(n+1, kind(0d0))
152 
153  DO i = 0, n
154  y_cc( i ) = y_domain%beg + 5d-1*dy*REAL(2*i+1, kind(0d0))
155  y_cb(i-1) = y_domain%beg + dy*REAL( i , kind(0d0))
156  END DO
157 
158  END IF
159 
160  y_cb(n) = y_domain%end
161 
162  IF(stretch_y) THEN
163 
164  DO j = 1, loops_y
165  DO i = -1, n
166  y_cb(i) = y_cb(i) / a_y * &
167  ( a_y + log(cosh( a_y*(y_cb(i) - y_a) )) &
168  + log(cosh( a_y*(y_cb(i) - y_b) )) &
169  - 2d0*log(cosh( 0.5d0*a_y*(y_b - y_a) )) )
170  END DO
171  END DO
172 
173  y_cc = (y_cb(0:n) + y_cb(-1:n-1))/2d0
174 
175  dy = minval(y_cb(0:n) - y_cb(-1:n-1))
176 
177  IF(num_procs > 1) CALL s_mpi_reduce_min(dy)
178 
179  END IF
180  ! ==================================================================
181 
182 
183  ! Grid Generation in the z-direction ===============================
184  IF(p == 0) RETURN
185 
186  dz = (z_domain%end - z_domain%beg) / REAL(p+1, kind(0d0))
187 
188  DO i = 0, p
189  z_cc( i ) = z_domain%beg + 5d-1*dz*REAL(2*i+1, kind(0d0))
190  z_cb(i-1) = z_domain%beg + dz*REAL( i , kind(0d0))
191  END DO
192 
193  z_cb(p) = z_domain%end
194 
195  IF(stretch_z) THEN
196 
197  DO j = 1, loops_z
198  DO i = -1, p
199  z_cb(i) = z_cb(i) / a_z * &
200  ( a_z + log(cosh( a_z*(z_cb(i) - z_a) )) &
201  + log(cosh( a_z*(z_cb(i) - z_b) )) &
202  - 2d0*log(cosh( 0.5d0*a_z*(z_b - z_a) )) )
203  END DO
204  END DO
205 
206  z_cc = (z_cb(0:p) + z_cb(-1:p-1))/2d0
207 
208  dz = minval(z_cb(0:p) - z_cb(-1:p-1))
209 
210  IF(num_procs > 1) CALL s_mpi_reduce_min(dz)
211 
212  END IF
213  ! ==================================================================
214 
215 
216  END SUBROUTINE s_generate_serial_grid ! ---------------------------------------
217 
218 
219 
226  SUBROUTINE s_generate_parallel_grid(dflt_int) !-------------------------
228  INTEGER, INTENT(IN) :: dflt_int
229 
230  ! Locations of cell boundaries
231  REAL(KIND(0d0)), ALLOCATABLE, DIMENSION(:) :: x_cb_glb, y_cb_glb, z_cb_glb
233 
234 
235  CHARACTER(LEN = path_len + name_len) :: file_loc
237 
238  INTEGER :: ifile, ierr, data_size
239  INTEGER, DIMENSION(MPI_STATUS_SIZE) :: status
240 
241  INTEGER :: i,j
242 
243  ALLOCATE(x_cb_glb(-1:m_glb))
244  ALLOCATE(y_cb_glb(-1:n_glb))
245  ALLOCATE(z_cb_glb(-1:p_glb))
246 
247  ! Grid generation in the x-direction
248  dx = (x_domain%end - x_domain%beg) / REAL(m_glb+1, kind(0d0))
249  DO i = 0, m_glb
250  x_cb_glb(i-1) = x_domain%beg + dx*REAL(i,kind(0d0))
251  END DO
252  x_cb_glb(m_glb) = x_domain%end
253  IF (stretch_x) THEN
254  DO j = 1, loops_x
255  DO i = -1, m_glb
256  x_cb_glb(i) = x_cb_glb(i) / a_x * &
257  (a_x + log(cosh( a_x*(x_cb_glb(i) - x_a))) &
258  + log(cosh( a_x*(x_cb_glb(i) - x_b))) &
259  - 2d0*log(cosh(0.5d0*a_x*(x_b - x_a))) )
260  END DO
261  END DO
262  END IF
263 
264  ! Grid generation in the y-direction
265  IF (n_glb > 0) THEN
266 
267  IF (grid_geometry == 2.AND.y_domain%beg.eq.0.0d0) THEN
268  dy = (y_domain%end - y_domain%beg) / REAL(2*n_glb+1, kind(0d0))
269  y_cb_glb(-1) = y_domain%beg
270  DO i = 1, n_glb
271  y_cb_glb(i-1) = y_domain%beg + dy*REAL(2*i-1, kind(0d0))
272  END DO
273  ELSE
274  dy = (y_domain%end - y_domain%beg) / REAL(n_glb+1, kind(0d0))
275  DO i = 0, n_glb
276  y_cb_glb(i-1) = y_domain%beg + dy*REAL(i, kind(0d0))
277  END DO
278  END IF
279  y_cb_glb(n_glb) = y_domain%end
280  IF (stretch_y) THEN
281  DO j = 1, loops_y
282  DO i = -1, n_glb
283  y_cb_glb(i) = y_cb_glb(i) / a_y * &
284  (a_y + log(cosh( a_y*(y_cb_glb(i) - y_a))) &
285  + log(cosh( a_y*(y_cb_glb(i) - y_b))) &
286  - 2d0*log(cosh(0.5d0*a_y*(y_b - y_a))) )
287  END DO
288  END DO
289  END IF
290 
291  ! Grid generation in the z-direction
292  IF (p_glb > 0) THEN
293  dz = (z_domain%end - z_domain%beg) / REAL(p_glb+1, kind(0d0))
294  DO i = 0, p_glb
295  z_cb_glb(i-1) = z_domain%beg + dz*REAL(i,kind(0d0))
296  END DO
297  z_cb_glb(p_glb) = z_domain%end
298  IF (stretch_z) THEN
299  DO j = 1, loops_z
300  DO i = -1, p_glb
301  z_cb_glb(i) = z_cb_glb(i) / a_z * &
302  (a_z + log(cosh( a_z*(z_cb_glb(i) - z_a))) &
303  + log(cosh( a_z*(z_cb_glb(i) - z_b))) &
304  - 2d0*log(cosh(0.5d0*a_z*(z_b - z_a))) )
305  END DO
306  END DO
307  END IF
308  END IF
309  END IF
310 
311  ! Write cell boundary locations to grid data files
312  file_loc = trim(case_dir) // '/restart_data' // trim(mpiiofs) // 'x_cb.dat'
313  data_size = m_glb+2
314  CALL mpi_file_open(mpi_comm_self,file_loc,ior(mpi_mode_wronly,mpi_mode_create),&
315  mpi_info_int,ifile,ierr)
316  CALL mpi_file_write(ifile,x_cb_glb,data_size,mpi_double_precision,status,ierr)
317  CALL mpi_file_close(ifile,ierr)
318 
319  IF (n > 0) THEN
320  file_loc = trim(case_dir) // '/restart_data' // trim(mpiiofs) // 'y_cb.dat'
321  data_size = n_glb+2
322  CALL mpi_file_open(mpi_comm_self,file_loc,ior(mpi_mode_wronly,mpi_mode_create),&
323  mpi_info_int,ifile,ierr)
324  CALL mpi_file_write(ifile,y_cb_glb,data_size,mpi_double_precision,status,ierr)
325  CALL mpi_file_close(ifile,ierr)
326 
327  IF (p > 0) THEN
328  file_loc = trim(case_dir) // '/restart_data' // trim(mpiiofs) // 'z_cb.dat'
329  data_size = p_glb+2
330  CALL mpi_file_open(mpi_comm_self,file_loc,ior(mpi_mode_wronly,mpi_mode_create),&
331  mpi_info_int,ifile,ierr)
332  CALL mpi_file_write(ifile,z_cb_glb,data_size,mpi_double_precision,status,ierr)
333  CALL mpi_file_close(ifile,ierr)
334  END IF
335  END IF
336 
337  DEALLOCATE(x_cb_glb, y_cb_glb, z_cb_glb)
338 
339  END SUBROUTINE s_generate_parallel_grid ! ------------------------------
340 
343  SUBROUTINE s_initialize_grid_module() ! -----------------------------------
344 
345  IF (parallel_io .NEQV. .true.) THEN
347  ELSE
349  END IF
350 
351  END SUBROUTINE s_initialize_grid_module ! ---------------------------------
352 
354  SUBROUTINE s_finalize_grid_module() ! --------------------------------
355 
356  s_generate_grid => null()
357 
358  END SUBROUTINE s_finalize_grid_module ! ------------------------------
359 
360 END MODULE m_grid
subroutine, public s_finalize_grid_module()
Deallocation procedures for the module.
Definition: m_grid.f90:355
integer mpi_info_int
MPI info for parallel IO with Lustre file systems.
real(kind(0d0)), dimension(:), allocatable y_cc
This module contains all of the parameters characterizing the computational domain, simulation algorithm, initial condition and the stiffened equation of state.
type(bounds_info) z_domain
Locations of the domain bounds in the x-, y- and z-coordinate directions.
type(bounds_info) x_domain
real(kind(0d0)), dimension(:), allocatable y_cb
real(kind(0d0)), dimension(:), allocatable x_cb
integer p_glb
Global number of cells in each direction.
This file contains the definitions of all of the custom-defined types used in the pre-process code...
real(kind(0d0)), dimension(:), allocatable z_cb
Locations of cell-boundaries (cb) in x-, y- and z-directions, respectively.
integer p
Number of cells in the x-, y- and z-coordinate directions.
real(kind(0d0)), dimension(:), allocatable x_cc
subroutine s_mpi_reduce_min(var_loc)
The following subroutine takes the inputted variable and determines its minimum value on the entire c...
integer grid_geometry
Cylindrical coordinates (either axisymmetric or full 3D)
character(len=path_len) case_dir
Case folder location.
character(len=name_len) mpiiofs
integer num_procs
Number of processors.
type(bounds_info) y_domain
subroutine, public s_initialize_grid_module()
Computation of parameters, allocation procedures, and/or any other tasks needed to properly setup the...
Definition: m_grid.f90:344
procedure(s_generate_abstract_grid), pointer, public s_generate_grid
Definition: m_grid.f90:77
subroutine, public s_generate_serial_grid(dflt_int)
The following subroutine generates either a uniform or non-uniform rectilinear grid in serial...
Definition: m_grid.f90:89
real(kind(0d0)), dimension(:), allocatable z_cc
Locations of cell-centers (cc) in x-, y- and z-directions, respectively.
real(kind(0d0)) dz
Minimum cell-widths in the x-, y- and z-coordinate directions.
logical parallel_io
Format of the data files.
This module serves as a proxy to the parameters and subroutines available in the MPI implementation&#39;s...
Definition: m_mpi_proxy.f90:39
logical stretch_z
Grid stretching flags for the x-, y- and z-coordinate directions.
subroutine, public s_generate_parallel_grid(dflt_int)
The following subroutine generates either a uniform or non-uniform rectilinear grid in parallel...
Definition: m_grid.f90:227
This module takes care of creating the rectilinear grid on which the data for the initial condition w...
Definition: m_grid.f90:44