72 class(*),
target,
intent(in) :: data
73 real(kind=
rp),
intent(in),
optional :: t
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
87 real(kind=
rp),
pointer :: dtlag(:), tlag(:)
88 type(
mesh_t),
pointer :: msh
89 type(mpi_status) :: status
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
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
120 if (
associated(data%ulag))
then
125 optional_fields = optional_fields + 1
130 if (
associated(data%s))
then
132 write_scalar = .true.
133 optional_fields = optional_fields + 2
135 write_scalar = .false.
138 if (
associated(data%tlag))
then
142 optional_fields = optional_fields + 4
144 write_dtlag = .false.
146 write_abvel = .false.
147 if (
associated(data%abx1))
then
154 optional_fields = optional_fields + 8
157 write_scalarlag = .false.
158 if (
associated(data%abs1))
then
162 optional_fields = optional_fields + 16
163 write_scalarlag = .true.
171 write(id_str,
'(i5.5)') this%counter
172 fname = trim(this%fname(1:suffix_pos-1))//id_str//
'.chkp'
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)
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)
192 byte_offset = byte_offset + &
194 call mpi_file_write_at_all(fh, byte_offset,u%x, u%dof%size(), &
197 mpi_offset = mpi_offset +&
200 byte_offset = mpi_offset + &
202 call mpi_file_write_at_all(fh, byte_offset, v%x, v%dof%size(), &
206 byte_offset = mpi_offset + &
208 call mpi_file_write_at_all(fh, byte_offset, w%x, w%dof%size(), &
212 byte_offset = mpi_offset + &
214 call mpi_file_write_at_all(fh, byte_offset, p%x, p%dof%size(), &
224 do i = 1, ulag%size()
225 byte_offset = mpi_offset + &
230 associate(x => ulag%lf(i)%x)
231 call mpi_file_write_at_all(fh, byte_offset, x, &
237 do i = 1, vlag%size()
238 byte_offset = mpi_offset + &
243 associate(x => vlag%lf(i)%x)
244 call mpi_file_write_at_all(fh, byte_offset, x, &
247 mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
250 do i = 1, wlag%size()
251 byte_offset = mpi_offset + &
252 dof_offset * int(mpi_real_prec_size, i8)
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)
260 mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
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)
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)
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)
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)
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)
324 mpi_offset = mpi_offset + n_glb_dofs * int(mpi_real_prec_size, i8)
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)
339 call mpi_file_close(fh, ierr)
341 this%counter = this%counter + 1
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
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
375 real(kind=rp) :: center_x, center_y, center_z
377 type(dofmap_t) :: dof
379 call this%check_exists()
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')
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
403 this%mesh2mesh = .false.
406 if (
associated(data%ulag))
then
415 if (
associated(data%s))
then
419 read_scalar = .false.
421 if (
associated(data%dtlag))
then
429 if (
associated(data%abx1))
then
438 read_scalarlag = .false.
439 if (
associated(data%abs1))
then
443 read_scalarlag = .true.
449 call neko_error(
'Invalid data')
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)
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
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')
476 if (gdim .eq. 3)
then
477 call this%chkp_Xh%init(gll, lx, lx, lx)
479 call this%chkp_Xh%init(gll, lx, lx)
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
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)
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)
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())
514 call this%space_interp%init(this%sim_Xh, this%chkp_Xh)
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)
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)
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)
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)
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)
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)
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)
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)
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)
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)
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)
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)
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)
630 call mpi_file_close(fh, ierr)
632 call this%global_interp%free()
633 call this%space_interp%free()
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
648 n = this%chkp_xh%lxyz*nel
649 allocate(read_array(n))
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
656 call this%global_interp%evaluate(x,read_array)
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)
662 x(i,1) = read_array(i)
665 deallocate(read_array)
Neko checkpoint file format.
subroutine chkp_file_write(this, data, t)
Write a Neko checkpoint.
subroutine chkp_file_read(this, data)
Load a checkpoint from file.
subroutine chkp_read_field(this, fh, byte_offset, x, nel)
type(mpi_comm) neko_comm
MPI communicator.
type(mpi_datatype) mpi_real_precision
MPI type for working precision of REAL types.
Defines a mapping of the degrees of freedom.
Implements global_interpolation given a dofmap.
Routines to interpolate between different spaces.
integer, public mpi_double_precision_size
Size of MPI type double precision.
integer, public mpi_real_prec_size
Size of working precision REAL types.
integer, public mpi_integer_size
Size of MPI type integer.
integer, parameter, public i8
integer, parameter, public dp
integer, parameter, public rp
Global precision used in computations.
Defines a function space.
pure integer function, public filename_suffix_pos(fname)
Find position (in the string) of a filename's suffix.
Interface for Neko checkpoint files.
Implements global interpolation for arbitrary points in the domain.
Interpolation between two space::space_t.
The function space for the SEM solution fields.