Neko  0.9.99
A portable framework for high-order spectral element flow simulations
chkp_file.f90
Go to the documentation of this file.
1 ! Copyright (c) 2021-2022, The Neko Authors
2 ! All rights reserved.
3 !
4 ! Redistribution and use in source and binary forms, with or without
5 ! modification, are permitted provided that the following conditions
6 ! are met:
7 !
8 ! * Redistributions of source code must retain the above copyright
9 ! notice, this list of conditions and the following disclaimer.
10 !
11 ! * Redistributions in binary form must reproduce the above
12 ! copyright notice, this list of conditions and the following
13 ! disclaimer in the documentation and/or other materials provided
14 ! with the distribution.
15 !
16 ! * Neither the name of the authors nor the names of its
17 ! contributors may be used to endorse or promote products derived
18 ! from this software without specific prior written permission.
19 !
20 ! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21 ! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22 ! LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
23 ! FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
24 ! COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
25 ! INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
26 ! BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
27 ! LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
28 ! CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
29 ! LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
30 ! ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
31 ! POSSIBILITY OF SUCH DAMAGE.
32 !
35 module chkp_file
36  use generic_file
37  use field_series
38  use checkpoint
39  use num_types
40  use field
41  use dofmap, only: dofmap_t
42  use device_math
43  use utils
44  use space
45  use mesh
46  use math
47  use interpolation
48  use neko_mpi_types
49  use comm
51  implicit none
52  private
53 
55  type, public, extends(generic_file_t) :: chkp_file_t
56  type(space_t), pointer :: chkp_xh
57  type(space_t), pointer :: sim_xh
58  type(interpolator_t) :: space_interp
59  type(global_interpolation_t) :: global_interp
60  logical :: mesh2mesh
61  contains
62  procedure :: read => chkp_file_read
63  procedure :: read_field => chkp_read_field
64  procedure :: write => chkp_file_write
65  end type chkp_file_t
66 
67 contains
68 
70  subroutine chkp_file_write(this, data, t)
71  class(chkp_file_t), intent(inout) :: this
72  class(*), target, intent(in) :: data
73  real(kind=rp), intent(in), optional :: t
74  real(kind=dp) :: time
75  character(len=5) :: id_str
76  character(len=1024) :: fname
77  integer :: ierr, suffix_pos, optional_fields
78  type(field_t), pointer :: u, v, w, p, s
79  type(field_t), pointer :: abx1,abx2
80  type(field_t), pointer :: aby1,aby2
81  type(field_t), pointer :: abz1,abz2
82  type(field_t), pointer :: abs1,abs2
83  type(field_series_t), pointer :: ulag => null()
84  type(field_series_t), pointer :: vlag => null()
85  type(field_series_t), pointer :: wlag => null()
86  type(field_series_t), pointer :: slag => null()
87  real(kind=rp), pointer :: dtlag(:), tlag(:)
88  type(mesh_t), pointer :: msh
89  type(mpi_status) :: status
90  type(mpi_file) :: fh
91  integer (kind=MPI_OFFSET_KIND) :: mpi_offset, byte_offset
92  integer(kind=i8) :: n_glb_dofs, dof_offset
93  logical :: write_lag, write_scalar, write_dtlag, write_scalarlag, write_abvel
94  integer :: i
95 
96  if (present(t)) then
97  time = real(t,dp)
98  else
99  time = 0d0
100  end if
101 
102  select type(data)
103  type is (chkp_t)
104 
105  if ( .not. associated(data%u) .or. &
106  .not. associated(data%v) .or. &
107  .not. associated(data%w) .or. &
108  .not. associated(data%p) ) then
109  call neko_error('Checkpoint not initialized')
110  end if
111 
112  u => data%u
113  v => data%v
114  w => data%w
115  p => data%p
116  msh => u%msh
117 
118  optional_fields = 0
119 
120  if (associated(data%ulag)) then
121  ulag => data%ulag
122  vlag => data%vlag
123  wlag => data%wlag
124  write_lag = .true.
125  optional_fields = optional_fields + 1
126  else
127  write_lag = .false.
128  end if
129 
130  if (associated(data%s)) then
131  s => data%s
132  write_scalar = .true.
133  optional_fields = optional_fields + 2
134  else
135  write_scalar = .false.
136  end if
137 
138  if (associated(data%tlag)) then
139  tlag => data%tlag
140  dtlag => data%dtlag
141  write_dtlag = .true.
142  optional_fields = optional_fields + 4
143  else
144  write_dtlag = .false.
145  end if
146  write_abvel = .false.
147  if (associated(data%abx1)) then
148  abx1 => data%abx1
149  abx2 => data%abx2
150  aby1 => data%aby1
151  aby2 => data%aby2
152  abz1 => data%abz1
153  abz2 => data%abz2
154  optional_fields = optional_fields + 8
155  write_abvel = .true.
156  end if
157  write_scalarlag = .false.
158  if (associated(data%abs1)) then
159  slag => data%slag
160  abs1 => data%abs1
161  abs2 => data%abs2
162  optional_fields = optional_fields + 16
163  write_scalarlag = .true.
164  end if
165  class default
166  call neko_error('Invalid data')
167  end select
168 
169 
170  suffix_pos = filename_suffix_pos(this%fname)
171  write(id_str, '(i5.5)') this%counter
172  fname = trim(this%fname(1:suffix_pos-1))//id_str//'.chkp'
173 
174 
175  dof_offset = int(msh%offset_el, i8) * int(u%Xh%lx * u%Xh%ly * u%Xh%lz, i8)
176  n_glb_dofs = int(u%Xh%lx * u%Xh%ly * u%Xh%lz, i8) * int(msh%glb_nelv, i8)
177 
178  call mpi_file_open(neko_comm, trim(fname), &
179  mpi_mode_wronly + mpi_mode_create, mpi_info_null, fh, ierr)
180  call mpi_file_write_all(fh, msh%glb_nelv, 1, mpi_integer, status, ierr)
181  call mpi_file_write_all(fh, msh%gdim, 1, mpi_integer, status, ierr)
182  call mpi_file_write_all(fh, u%Xh%lx, 1, mpi_integer, status, ierr)
183  call mpi_file_write_all(fh, optional_fields, 1, mpi_integer, status, ierr)
184  call mpi_file_write_all(fh, time, 1, mpi_double_precision, status, ierr)
185 
186 
187  !
188  ! Dump mandatory checkpoint data
189  !
190 
191  byte_offset = 4_i8 * int(mpi_integer_size,i8) + int(mpi_double_precision_size,i8)
192  byte_offset = byte_offset + &
193  dof_offset * int(mpi_real_prec_size, i8)
194  call mpi_file_write_at_all(fh, byte_offset,u%x, u%dof%size(), &
195  mpi_real_precision, status, ierr)
196  mpi_offset = 4_i8 * int(mpi_integer_size,i8) + int(mpi_double_precision_size,i8)
197  mpi_offset = mpi_offset +&
198  n_glb_dofs * int(mpi_real_prec_size, i8)
199 
200  byte_offset = mpi_offset + &
201  dof_offset * int(mpi_real_prec_size, i8)
202  call mpi_file_write_at_all(fh, byte_offset, v%x, v%dof%size(), &
203  mpi_real_precision, status, ierr)
204  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
205 
206  byte_offset = mpi_offset + &
207  dof_offset * int(mpi_real_prec_size, i8)
208  call mpi_file_write_at_all(fh, byte_offset, w%x, w%dof%size(), &
209  mpi_real_precision, status, ierr)
210  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
211 
212  byte_offset = mpi_offset + &
213  dof_offset * int(mpi_real_prec_size, i8)
214  call mpi_file_write_at_all(fh, byte_offset, p%x, p%dof%size(), &
215  mpi_real_precision, status, ierr)
216  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
217 
218  !
219  ! Dump optional payload
220  !
221 
222  if (write_lag) then
223 
224  do i = 1, ulag%size()
225  byte_offset = mpi_offset + &
226  dof_offset * int(mpi_real_prec_size, i8)
227  ! We should not need this extra associate block, ant it works
228  ! great without it for GNU, Intel, NEC and Cray, but throws an
229  ! ICE with NAG.
230  associate(x => ulag%lf(i)%x)
231  call mpi_file_write_at_all(fh, byte_offset, x, &
232  ulag%lf(i)%dof%size(), mpi_real_precision, status, ierr)
233  end associate
234  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
235  end do
236 
237  do i = 1, vlag%size()
238  byte_offset = mpi_offset + &
239  dof_offset * int(mpi_real_prec_size, i8)
240  ! We should not need this extra associate block, ant it works
241  ! great without it for GNU, Intel, NEC and Cray, but throws an
242  ! ICE with NAG.
243  associate(x => vlag%lf(i)%x)
244  call mpi_file_write_at_all(fh, byte_offset, x, &
245  vlag%lf(i)%dof%size(), mpi_real_precision, status, ierr)
246  end associate
247  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
248  end do
249 
250  do i = 1, wlag%size()
251  byte_offset = mpi_offset + &
252  dof_offset * int(mpi_real_prec_size, i8)
253  ! We should not need this extra associate block, ant it works
254  ! great without it for GNU, Intel, NEC and Cray, but throws an
255  ! ICE with NAG.
256  associate(x => wlag%lf(i)%x)
257  call mpi_file_write_at_all(fh, byte_offset, x, &
258  wlag%lf(i)%dof%size(), mpi_real_precision, status, ierr)
259  end associate
260  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
261  end do
262 
263  end if
264 
265  if (write_scalar) then
266  byte_offset = mpi_offset + &
267  dof_offset * int(mpi_real_prec_size, i8)
268  call mpi_file_write_at_all(fh, byte_offset, s%x, p%dof%size(), &
269  mpi_real_precision, status, ierr)
270  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
271  end if
272 
273  if (write_dtlag) then
274  call mpi_file_write_at_all(fh, mpi_offset, tlag, 10, mpi_real_precision, status, ierr)
275  mpi_offset = mpi_offset + 10_i8 * int(mpi_real_prec_size, i8)
276  call mpi_file_write_at_all(fh, mpi_offset, dtlag, 10, mpi_real_precision, status, ierr)
277  mpi_offset = mpi_offset + 10_i8 * int(mpi_real_prec_size, i8)
278  end if
279 
280  if (write_abvel) then
281  byte_offset = mpi_offset + &
282  dof_offset * int(mpi_real_prec_size, i8)
283  call mpi_file_write_at_all(fh, byte_offset, abx1%x, abx1%dof%size(), &
284  mpi_real_precision, status, ierr)
285  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
286  byte_offset = mpi_offset + &
287  dof_offset * int(mpi_real_prec_size, i8)
288  call mpi_file_write_at_all(fh, byte_offset, abx2%x, abx1%dof%size(), &
289  mpi_real_precision, status, ierr)
290  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
291  byte_offset = mpi_offset + &
292  dof_offset * int(mpi_real_prec_size, i8)
293  call mpi_file_write_at_all(fh, byte_offset, aby1%x, abx1%dof%size(), &
294  mpi_real_precision, status, ierr)
295  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
296  byte_offset = mpi_offset + &
297  dof_offset * int(mpi_real_prec_size, i8)
298  call mpi_file_write_at_all(fh, byte_offset, aby2%x, abx1%dof%size(), &
299  mpi_real_precision, status, ierr)
300  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
301  byte_offset = mpi_offset + &
302  dof_offset * int(mpi_real_prec_size, i8)
303  call mpi_file_write_at_all(fh, byte_offset, abz1%x, abx1%dof%size(), &
304  mpi_real_precision, status, ierr)
305  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
306  byte_offset = mpi_offset + &
307  dof_offset * int(mpi_real_prec_size, i8)
308  call mpi_file_write_at_all(fh, byte_offset, abz2%x, abx1%dof%size(), &
309  mpi_real_precision, status, ierr)
310  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
311  end if
312 
313  if (write_scalarlag) then
314  do i = 1, slag%size()
315  byte_offset = mpi_offset + &
316  dof_offset * int(mpi_real_prec_size, i8)
317  ! We should not need this extra associate block, ant it works
318  ! great without it for GNU, Intel, NEC and Cray, but throws an
319  ! ICE with NAG.
320  associate(x => slag%lf(i)%x)
321  call mpi_file_write_at_all(fh, byte_offset, x, &
322  slag%lf(i)%dof%size(), mpi_real_precision, status, ierr)
323  end associate
324  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
325  end do
326 
327  byte_offset = mpi_offset + &
328  dof_offset * int(mpi_real_prec_size, i8)
329  call mpi_file_write_at_all(fh, byte_offset, abs1%x, abx1%dof%size(), &
330  mpi_real_precision, status, ierr)
331  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
332  byte_offset = mpi_offset + &
333  dof_offset * int(mpi_real_prec_size, i8)
334  call mpi_file_write_at_all(fh, byte_offset, abs2%x, abx1%dof%size(), &
335  mpi_real_precision, status, ierr)
336  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
337  end if
338 
339  call mpi_file_close(fh, ierr)
340 
341  this%counter = this%counter + 1
342 
343  end subroutine chkp_file_write
344 
346  subroutine chkp_file_read(this, data)
347  class(chkp_file_t) :: this
348  class(*), target, intent(inout) :: data
349  type(chkp_t), pointer :: chkp
350  character(len=5) :: id_str
351  character(len=1024) :: fname
352  integer :: ierr, suffix_pos
353  type(field_t), pointer :: u, v, w, p, s
354  type(field_series_t), pointer :: ulag => null()
355  type(field_series_t), pointer :: vlag => null()
356  type(field_series_t), pointer :: wlag => null()
357  type(field_series_t), pointer :: slag => null()
358  type(mesh_t), pointer :: msh
359  type(mpi_status) :: status
360  type(mpi_file) :: fh
361  type(field_t), pointer :: abx1,abx2
362  type(field_t), pointer :: aby1,aby2
363  type(field_t), pointer :: abz1,abz2
364  type(field_t), pointer :: abs1,abs2
365  real(kind=rp), allocatable :: x_coord(:,:,:,:)
366  real(kind=rp), allocatable :: y_coord(:,:,:,:)
367  real(kind=rp), allocatable :: z_coord(:,:,:,:)
368  real(kind=rp), pointer :: dtlag(:), tlag(:)
369  integer (kind=MPI_OFFSET_KIND) :: mpi_offset, byte_offset
370  integer(kind=i8) :: n_glb_dofs, dof_offset
371  integer :: glb_nelv, gdim, lx, have_lag, have_scalar, nel, optional_fields, have_dtlag
372  integer :: have_abvel, have_scalarlag
373  logical :: read_lag, read_scalar, read_dtlag, read_abvel, read_scalarlag
374  real(kind=rp) :: tol
375  real(kind=rp) :: center_x, center_y, center_z
376  integer :: i, e
377  type(dofmap_t) :: dof
378 
379  call this%check_exists()
380 
381  select type(data)
382  type is (chkp_t)
383 
384  if ( .not. associated(data%u) .or. &
385  .not. associated(data%v) .or. &
386  .not. associated(data%w) .or. &
387  .not. associated(data%p) ) then
388  call neko_error('Checkpoint not initialized')
389  end if
390 
391  u => data%u
392  v => data%v
393  w => data%w
394  p => data%p
395  this%chkp_Xh => data%previous_Xh
397  if (allocated(data%previous_mesh%elements)) then
398  msh => data%previous_mesh
399  this%mesh2mesh = .true.
400  tol = data%mesh2mesh_tol
401  else
402  msh => u%msh
403  this%mesh2mesh = .false.
404  end if
405 
406  if (associated(data%ulag)) then
407  ulag => data%ulag
408  vlag => data%vlag
409  wlag => data%wlag
410  read_lag = .true.
411  else
412  read_lag = .false.
413  end if
414 
415  if (associated(data%s)) then
416  s => data%s
417  read_scalar = .true.
418  else
419  read_scalar = .false.
420  end if
421  if (associated(data%dtlag)) then
422  dtlag => data%dtlag
423  tlag => data%tlag
424  read_dtlag = .true.
425  else
426  read_dtlag = .false.
427  end if
428  read_abvel = .false.
429  if (associated(data%abx1)) then
430  abx1 => data%abx1
431  abx2 => data%abx2
432  aby1 => data%aby1
433  aby2 => data%aby2
434  abz1 => data%abz1
435  abz2 => data%abz2
436  read_abvel = .true.
437  end if
438  read_scalarlag = .false.
439  if (associated(data%abs1)) then
440  slag => data%slag
441  abs1 => data%abs1
442  abs2 => data%abs2
443  read_scalarlag = .true.
444  end if
445 
446  chkp => data
447 
448  class default
449  call neko_error('Invalid data')
450  end select
451 
452 
453 
454  call mpi_file_open(neko_comm, trim(this%fname), &
455  mpi_mode_rdonly, mpi_info_null, fh, ierr)
456  call mpi_file_read_all(fh, glb_nelv, 1, mpi_integer, status, ierr)
457  call mpi_file_read_all(fh, gdim, 1, mpi_integer, status, ierr)
458  call mpi_file_read_all(fh, lx, 1, mpi_integer, status, ierr)
459  call mpi_file_read_all(fh, optional_fields, 1, mpi_integer, status, ierr)
460  call mpi_file_read_all(fh, chkp%t, 1, mpi_double_precision, status, ierr)
461 
462  have_lag = mod(optional_fields,2)/1
463  have_scalar = mod(optional_fields,4)/2
464  have_dtlag = mod(optional_fields,8)/4
465  have_abvel = mod(optional_fields,16)/8
466  have_scalarlag = mod(optional_fields,32)/16
467 
468  if ( ( glb_nelv .ne. msh%glb_nelv ) .or. &
469  ( gdim .ne. msh%gdim) .or. &
470  ( (have_lag .eq. 0) .and. (read_lag) ) .or. &
471  ( (have_scalar .eq. 0) .and. (read_scalar) ) ) then
472  call neko_error('Checkpoint does not match case')
473  end if
474  nel = msh%nelv
475  this%sim_Xh => u%Xh
476  if (gdim .eq. 3) then
477  call this%chkp_Xh%init(gll, lx, lx, lx)
478  else
479  call this%chkp_Xh%init(gll, lx, lx)
480  end if
481  if (this%mesh2mesh) then
482  call dof%init(msh, this%chkp_Xh)
483  allocate(x_coord(u%Xh%lx,u%Xh%ly,u%Xh%lz,u%msh%nelv))
484  allocate(y_coord(u%Xh%lx,u%Xh%ly,u%Xh%lz,u%msh%nelv))
485  allocate(z_coord(u%Xh%lx,u%Xh%ly,u%Xh%lz,u%msh%nelv))
490  do e = 1, u%dof%msh%nelv
491  center_x = 0d0
492  center_y = 0d0
493  center_z = 0d0
494  do i = 1,u%dof%Xh%lxyz
495  center_x = center_x + u%dof%x(i,1,1,e)
496  center_y = center_y + u%dof%y(i,1,1,e)
497  center_z = center_z + u%dof%z(i,1,1,e)
498  end do
499  center_x = center_x/u%Xh%lxyz
500  center_y = center_y/u%Xh%lxyz
501  center_z = center_z/u%Xh%lxyz
502  do i = 1,u%dof%Xh%lxyz
503  x_coord(i,1,1,e) = u%dof%x(i,1,1,e) - tol*(u%dof%x(i,1,1,e)-center_x)
504  y_coord(i,1,1,e) = u%dof%y(i,1,1,e) - tol*(u%dof%y(i,1,1,e)-center_y)
505  z_coord(i,1,1,e) = u%dof%z(i,1,1,e) - tol*(u%dof%z(i,1,1,e)-center_z)
506  end do
507  end do
508  call this%global_interp%init(dof,tol=tol)
509  call this%global_interp%find_points(x_coord,y_coord,z_coord,u%dof%size())
510  deallocate(x_coord)
511  deallocate(y_coord)
512  deallocate(z_coord)
513  else
514  call this%space_interp%init(this%sim_Xh, this%chkp_Xh)
515  end if
516  dof_offset = int(msh%offset_el, i8) * int(this%chkp_Xh%lxyz, i8)
517  n_glb_dofs = int(this%chkp_Xh%lxyz, i8) * int(msh%glb_nelv, i8)
518 
519  !
520  ! Read mandatory checkpoint data
521  !
522 
523  byte_offset = 4_i8 * int(mpi_integer_size,i8) + int(mpi_double_precision_size,i8)
524  byte_offset = byte_offset + &
525  dof_offset * int(mpi_real_prec_size, i8)
526  call this%read_field(fh, byte_offset, u%x, nel)
527  mpi_offset = 4_i8 * int(mpi_integer_size,i8) + int(mpi_double_precision_size,i8)
528  mpi_offset = mpi_offset +&
529  n_glb_dofs * int(mpi_real_prec_size, i8)
530 
531  byte_offset = mpi_offset + &
532  dof_offset * int(mpi_real_prec_size, i8)
533  call this%read_field(fh, byte_offset, v%x, nel)
534  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
535 
536  byte_offset = mpi_offset + &
537  dof_offset * int(mpi_real_prec_size, i8)
538  call this%read_field(fh, byte_offset, w%x, nel)
539  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
540 
541  byte_offset = mpi_offset + &
542  dof_offset * int(mpi_real_prec_size, i8)
543  call this%read_field(fh, byte_offset, p%x, nel)
544  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
545 
546  !
547  ! Read optional payload
548  !
549 
550  if (read_lag) then
551  do i = 1, ulag%size()
552  byte_offset = mpi_offset + &
553  dof_offset * int(mpi_real_prec_size, i8)
554  call this%read_field(fh, byte_offset, ulag%lf(i)%x, nel)
555  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
556  end do
557 
558  do i = 1, vlag%size()
559  byte_offset = mpi_offset + &
560  dof_offset * int(mpi_real_prec_size, i8)
561  call this%read_field(fh, byte_offset, vlag%lf(i)%x, nel)
562  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
563  end do
564 
565  do i = 1, wlag%size()
566  byte_offset = mpi_offset + &
567  dof_offset * int(mpi_real_prec_size, i8)
568  call this%read_field(fh, byte_offset, wlag%lf(i)%x, nel)
569  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
570  end do
571  end if
572 
573  if (read_scalar) then
574  byte_offset = mpi_offset + &
575  dof_offset * int(mpi_real_prec_size, i8)
576  call this%read_field(fh, byte_offset, s%x, nel)
577  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
578  end if
579 
580  if (read_dtlag .and. have_dtlag .eq. 1) then
581  call mpi_file_read_at_all(fh, mpi_offset, tlag, 10, mpi_real_precision, status, ierr)
582  mpi_offset = mpi_offset + 10_i8 * int(mpi_real_prec_size, i8)
583  call mpi_file_read_at_all(fh, mpi_offset, dtlag, 10, mpi_real_precision, status, ierr)
584  mpi_offset = mpi_offset + 10_i8 * int(mpi_real_prec_size, i8)
585  end if
586 
587  if (read_abvel .and. have_abvel .eq. 1) then
588  byte_offset = mpi_offset + &
589  dof_offset * int(mpi_real_prec_size, i8)
590  call this%read_field(fh, byte_offset, abx1%x, nel)
591  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
592  byte_offset = mpi_offset + &
593  dof_offset * int(mpi_real_prec_size, i8)
594  call this%read_field(fh, byte_offset, abx2%x, nel)
595  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
596  byte_offset = mpi_offset + &
597  dof_offset * int(mpi_real_prec_size, i8)
598  call this%read_field(fh, byte_offset, aby1%x, nel)
599  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
600  byte_offset = mpi_offset + &
601  dof_offset * int(mpi_real_prec_size, i8)
602  call this%read_field(fh, byte_offset, aby2%x, nel)
603  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
604  byte_offset = mpi_offset + &
605  dof_offset * int(mpi_real_prec_size, i8)
606  call this%read_field(fh, byte_offset, abz1%x, nel)
607  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
608  byte_offset = mpi_offset + &
609  dof_offset * int(mpi_real_prec_size, i8)
610  call this%read_field(fh, byte_offset, abz2%x, nel)
611  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
612  end if
613  if (read_scalarlag .and. have_scalarlag .eq. 1) then
614  do i = 1, slag%size()
615  byte_offset = mpi_offset + &
616  dof_offset * int(mpi_real_prec_size, i8)
617  call this%read_field(fh, byte_offset, slag%lf(i)%x, nel)
618  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
619  end do
620  byte_offset = mpi_offset + &
621  dof_offset * int(mpi_real_prec_size, i8)
622  call this%read_field(fh, byte_offset, abs1%x, nel)
623  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
624  byte_offset = mpi_offset + &
625  dof_offset * int(mpi_real_prec_size, i8)
626  call this%read_field(fh, byte_offset, abs2%x, nel)
627  mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
628  end if
629 
630  call mpi_file_close(fh, ierr)
631 
632  call this%global_interp%free()
633  call this%space_interp%free()
634 
635  end subroutine chkp_file_read
636 
637  subroutine chkp_read_field(this, fh, byte_offset, x, nel)
638  class(chkp_file_t) :: this
639  type(mpi_file) :: fh
640  integer(kind=MPI_OFFSET_KIND) :: byte_offset
641  integer, intent(in) :: nel
642  real(kind=rp) :: x(this%sim_Xh%lxyz, nel)
643  real(kind=rp), allocatable :: read_array(:)
644  integer :: nel_stride, frac_space
645  type(mpi_status) :: status
646  integer :: ierr, i, n
647 
648  n = this%chkp_xh%lxyz*nel
649  allocate(read_array(n))
650 
651  call rzero(read_array,n)
652  call mpi_file_read_at_all(fh, byte_offset, read_array, &
653  n, mpi_real_precision, status, ierr)
654  if (this%mesh2mesh) then
655  x = 0.0_rp
656  call this%global_interp%evaluate(x,read_array)
657 
658  else if (this%sim_Xh%lx .ne. this%chkp_Xh%lx) then
659  call this%space_interp%map_host(x, read_array, nel, this%sim_Xh)
660  else
661  do i = 1,n
662  x(i,1) = read_array(i)
663  end do
664  end if
665  deallocate(read_array)
666  end subroutine chkp_read_field
667 
668 end module chkp_file
double real
Definition: device_config.h:12
Defines a checkpoint.
Definition: checkpoint.f90:34
Neko checkpoint file format.
Definition: chkp_file.f90:35
subroutine chkp_file_write(this, data, t)
Write a Neko checkpoint.
Definition: chkp_file.f90:71
subroutine chkp_file_read(this, data)
Load a checkpoint from file.
Definition: chkp_file.f90:347
subroutine chkp_read_field(this, fh, byte_offset, x, nel)
Definition: chkp_file.f90:638
Definition: comm.F90:1
type(mpi_comm) neko_comm
MPI communicator.
Definition: comm.F90:16
type(mpi_datatype) mpi_real_precision
MPI type for working precision of REAL types.
Definition: comm.F90:23
Defines a mapping of the degrees of freedom.
Definition: dofmap.f90:35
Stores a series fields.
Defines a field.
Definition: field.f90:34
Implements global_interpolation given a dofmap.
Routines to interpolate between different spaces.
Definition: math.f90:60
Defines a mesh.
Definition: mesh.f90:34
MPI derived types.
Definition: mpi_types.f90:34
integer, public mpi_double_precision_size
Size of MPI type double precision.
Definition: mpi_types.f90:61
integer, public mpi_real_prec_size
Size of working precision REAL types.
Definition: mpi_types.f90:65
integer, public mpi_integer_size
Size of MPI type integer.
Definition: mpi_types.f90:63
integer, parameter, public i8
Definition: num_types.f90:7
integer, parameter, public dp
Definition: num_types.f90:9
integer, parameter, public rp
Global precision used in computations.
Definition: num_types.f90:12
Defines a function space.
Definition: space.f90:34
Utilities.
Definition: utils.f90:35
pure integer function, public filename_suffix_pos(fname)
Find position (in the string) of a filename's suffix.
Definition: utils.f90:56
Interface for Neko checkpoint files.
Definition: chkp_file.f90:55
A generic file handler.
Implements global interpolation for arbitrary points in the domain.
Interpolation between two space::space_t.
The function space for the SEM solution fields.
Definition: space.f90:62