!> @file checkpoint_2d.f90
!> @brief Solution-only checkpoint/restart for the 2D Euler solver (MPI-correct),
!!        built on the dimension-blind `checkpoint_io` shared core (R4).
!!
!! Unformatted stream binary, GLOBAL format (rank-count-independent): the
!! solution is gathered to rank 0 and written over the global (nx_g, ny_g) grid.
!! The on-disk format is the unified format (magic 42); scalars iter, t,
!! extents (= [nx_g, ny_g]), neq2d, dt; grid record grid_type + a 2-component
!! fingerprint; the global ub array; then (v6, Q2d) the trailing records,
!! written/read unconditionally at fixed positions and mirroring 1D's v5
!! record ORDER (see checkpoint.f90's module header):
!!   1. the Kahan `t_comp` compensation scalar;
!!   2. the halo width `h` (validated on read against the current scheme's
!!      halo width);
!!   3. the GLOBAL boundary-ghost frame — the four width-h ghost bands
!!      outside the global interior, corners included, in the fixed order
!!      left, right, bottom, top (corners live in the bottom/top bands) —
!!      moved via `gather_boundary_ghosts_to_root` /
!!      `scatter_boundary_ghosts_from_root` (solution_gather_2d.F90). See
!!      the record-block comment in `write_checkpoint_2d` for WHY the frame
!!      must be persisted (multi-stage RK stage-snapshot ghost history).
!!
!! The version number is `checkpoint_io`'s single shared `ckpt_version`
!! constant. F2 bumped it 4→5 to add two 1D-only trailing records (see
!! checkpoint.f90's module header); 2D's on-disk body was byte-for-byte
!! unchanged by that bump, it just stamped version 5 in the header. Q2d bumps
!! it 5→6 to add 2D's own trailing records above (t_comp, then — same
!! initiative, after the restart bit-identity evidence gate — the halo width
!! and ghost-frame records); 1D's on-disk body is unchanged by THIS bump (a
!! freshly WRITTEN 1D checkpoint differs from its v5 form only in the
!! stamped header version). That does NOT make existing v5 files readable,
!! in either dimension — the strict version-equality check in
!! `read_ckpt_header` rejects legacy v4 and now also v5 files (1D or 2D),
!! the same owner-authorized break as before.
!!
!! The mesh is NOT stored — on restart it is rebuilt from grid_file via
!! build_mesh_2d_global/local. The fingerprint (sum and sum-of-squares of cell
!! centres, via par_sum_real over local cells) detects a changed grid even at the
!! same nx/ny. 2D is explicit-only, so there is no BDF2 array (cf. 1D).
!!
!! MPI: write gathers to root, writes on rank 0, agrees via `agree_outcome`
!! (`checkpoint_io`'s collective wrapper over `par_lor`). read opens the global
!! file read-only on every rank, validates, slices this rank's block via
!! decomp_2d % ix_first_global/iy_first_global; one collective decision point
!! via `conclude_read`.
!!
!! Version compatibility: this is the unified v6 format (magic 42, version 6,
!! dim=2, extents=[nx_g, ny_g]). Legacy 2D v1 checkpoints (magic 4242), and
!! now also v4/v5 checkpoints, are NOT readable — `read_ckpt_header`'s
!! magic/strict-version-equality check rejects them cleanly with a clear
!! error (owner-authorized break; no user-generated checkpoints exist).
!! Everything else — the global gather, the pointer-file mechanics, and the
!! MPI contract above — is unchanged.
module checkpoint_2d
  use precision, only: wp
  use solver_state_2d, only: solver_state_2d_t, neq2d
  use logger, only: log_info, log_warn
  use mpi_cart_2d, only: cart_rank_2d
  use solution_gather_2d, only: gather_solution_to_root_2d, &
                                gather_boundary_ghosts_to_root, scatter_boundary_ghosts_from_root
  use mpi_runtime, only: parallel_fatal
  use parallel_reductions, only: par_lor, par_sum_real
  use checkpoint_io, only: ckpt_filename, update_latest_pointer, resolve_latest_pointer, &
                           write_ckpt_header, read_ckpt_header, write_grid_record, &
                           read_grid_record_and_check, agree_outcome, conclude_read
  implicit none
  private
  public :: write_checkpoint_2d, read_checkpoint_2d

contains

  !> Rank-independent geometric fingerprint of the grid: (sum, sum-of-squares) of
  !! cell-centre coordinates. For curvilinear grids it reduces over local cells
  !! (each cell owned once); for uniform grids it is a deterministic scalar combo
  !! identical on every rank.
  function grid_fingerprint_2d(state) result(fp)
    type(solver_state_2d_t), intent(in) :: state
    real(wp) :: fp(2)
    real(wp) :: s1, s2
    integer :: i, j
    ! NOTE: the UNIFORM branch is purely local (no collective); the CURVILINEAR
    ! branch calls par_sum_real (collective). Callers must invoke this function
    ! on ALL ranks at a uniform call site — before any rank can diverge or
    ! early-return.
    if (state % mesh % uniform) then
      fp(1) = state % cfg % x_left + state % cfg % x_right + &
              state % cfg % y_left + state % cfg % y_right
      fp(2) = state % dx + state % dy + &
              real(state % decomp_2d % nx_global + state % decomp_2d % ny_global, wp)
    else
      s1 = 0.0_wp; s2 = 0.0_wp
      do j = 1, state % ny_local
        do i = 1, state % nx_local
          s1 = s1 + state % mesh % xc(i, j) + state % mesh % yc(i, j)
          s2 = s2 + state % mesh % xc(i, j)**2 + state % mesh % yc(i, j)**2
        end do
      end do
      fp(1) = par_sum_real(s1)
      fp(2) = par_sum_real(s2)
    end if
  end function grid_fingerprint_2d

  !> Write a checkpoint: gather to root, write the global file on rank 0, agree
  !! outcome across ranks, update the latest_checkpoint_2d pointer file.
  !!
  !! @param[in] t_comp  Optional Kahan compensation term for `t` (v6, Q2d);
  !!   absent is written as 0, matching a fresh (non-restarted) run's value.
  subroutine write_checkpoint_2d(state, base, t, iter, is_ok, message, written_file, t_comp)
    type(solver_state_2d_t), intent(in) :: state
    character(len=*), intent(in) :: base
    real(wp), intent(in) :: t
    integer, intent(in) :: iter
    logical, intent(out), optional :: is_ok
    character(len=*), intent(out), optional :: message
    !> Path of the file actually written (`<base>_<iter>.bin`), set only on
    !! success; mirrors checkpoint.f90's write_checkpoint contract.
    character(len=*), intent(out), optional :: written_file
    real(wp), intent(in), optional :: t_comp
    character(len=512) :: fname
    integer :: u, info, nx_g, ny_g, h, dstat, wr
    logical :: local_ok, lok
    character(len=256) :: lmsg, gmsg
    real(wp) :: fp(2)
    real(wp) :: t_comp_val
    real(wp), allocatable :: glob(:, :, :)
    real(wp), allocatable :: gh_l(:, :, :), gh_r(:, :, :), gh_b(:, :, :), gh_t(:, :, :)
    logical :: gok

    if (present(is_ok)) is_ok = .true.
    if (present(message)) message = ''
    if (present(written_file)) written_file = ''
    fname = ckpt_filename(base, iter)
    nx_g = state % decomp_2d % nx_global
    ny_g = state % decomp_2d % ny_global
    h = state % halo_width

    ! Collective gather (every rank participates); fingerprint via reductions.
    call gather_solution_to_root_2d(state, glob, gok, gmsg)
    fp = grid_fingerprint_2d(state)
    if (par_lor(.not. gok)) then
      if (present(is_ok)) is_ok = .false.
      if (present(message)) message = 'checkpoint_2d: gather failed'
      return
    end if

    ! v6 (Q2d): gather the GLOBAL boundary-ghost frame to root. FIXED
    ! collective position reached by every rank unconditionally (mirroring
    ! 1D's v5 edge-ghost gather and the tile gather above), not gated on
    ! rank-0-only local_ok, so it cannot desync ranks.
    !
    ! Only root populates and writes the frame, so size the buffers globally
    ! on root and to zero off-root, mirroring gather_solution_to_root_2d's
    ! treatment of `glob`. The COLLECTIVE CALL below stays unconditional --
    ! it is the buffer that is rank-dependent, not the call position -- so
    ! this cannot desync ranks. (Sizing them globally on every rank was a
    ! bounded but real per-rank waste and a remote OOM surface on large grids.)
    wr = cart_rank_2d()
    if (wr == 0) then
      allocate (gh_l(neq2d, h, ny_g), gh_r(neq2d, h, ny_g), &
                gh_b(neq2d, nx_g + 2 * h, h), gh_t(neq2d, nx_g + 2 * h, h), stat=info)
    else
      allocate (gh_l(neq2d, 0, 0), gh_r(neq2d, 0, 0), &
                gh_b(neq2d, 0, 0), gh_t(neq2d, 0, 0), stat=info)
    end if
    if (info /= 0) call parallel_fatal('checkpoint_2d: cannot allocate ghost-frame gather buffers')
    call gather_boundary_ghosts_to_root(state, gh_l, gh_r, gh_b, gh_t)

    local_ok = .true.
    if (wr < 0) local_ok = .false.
    if (local_ok .and. wr == 0) then
      open (newunit=u, file=trim(fname), status='replace', action='write', &
            form='unformatted', access='stream', iostat=info)
      if (info /= 0) then
        call log_warn('checkpoint_2d: cannot open "'//trim(fname)//'" for writing')
        local_ok = .false.
      else
        call write_ckpt_header(u, 2, iter, t, [nx_g, ny_g], neq2d, state % dt, info)
        if (info == 0) call write_grid_record(u, state % cfg % grid_type, fp, info)
        if (info == 0) write (u, iostat=info) glob
        if (info == 0) then
          ! v6 (Q2d): trailing Kahan t_comp scalar, mirroring 1D's v5 F2
          ! record (checkpoint.f90). Absent argument writes 0, matching a
          ! fresh (non-restarted) run's compensation value.
          t_comp_val = 0.0_wp
          if (present(t_comp)) t_comp_val = t_comp
          write (u, iostat=info) t_comp_val
        end if
        ! v6 (Q2d): halo width + GLOBAL boundary-ghost frame, mirroring 1D's
        ! v5 record order. WHY the frame must be persisted (root cause from
        ! the restart bit-identity evidence gate, test_restart_bit_identity_2d):
        ! the multi-stage whole-array RK kernels (stepper_kernels) snapshot
        ! the solution INCLUDING its ghost words (s1 = u) BEFORE the step's
        ! first rhs call refills the boundary ghosts, so each step's
        ! end-of-step ghost words are a weighted blend of that step's fresh
        ! boundary values and the PREVIOUS step's ghost snapshot — a
        ! geometrically decaying history (e.g. TVD-RK3 keeps weight 1/3 per
        ! step; (1/3)^10 ~ 1.7e-5 matched the observed post-restart ghost
        ! divergence) that is NOT reconstructable from the checkpointed
        ! interior. The interior is immune — the BC refill precedes every
        ! residual read, so the stale blend is never consumed by the
        ! dynamics — but the frame records exist to make the WHOLE-ARRAY
        ! restart guarantee (interior AND ghosts bit-identical) hold,
        ! matching 1D's v5 design. Single-stage euler needs no history and
        ! was bit-identical without these records; they are written for
        ! every scheme for one uniform format.
        if (info == 0) write (u, iostat=info) h
        if (info == 0) write (u, iostat=info) gh_l
        if (info == 0) write (u, iostat=info) gh_r
        if (info == 0) write (u, iostat=info) gh_b
        if (info == 0) write (u, iostat=info) gh_t
        if (info /= 0) local_ok = .false.
        close (u, iostat=info)
        if (info /= 0) local_ok = .false.
      end if
    end if
    if (allocated(glob)) then
      deallocate (glob, stat=dstat)
      if (dstat /= 0) local_ok = .false.
    end if
    deallocate (gh_l, gh_r, gh_b, gh_t, stat=dstat)
    if (dstat /= 0) local_ok = .false.

    if (agree_outcome(local_ok, 'checkpoint_2d: global checkpoint write failed (on at least one rank)', &
                      is_ok, message)) return

    local_ok = .true.
    if (wr == 0) then
      call log_info('checkpoint_2d: wrote "'//trim(fname)//'"')
      call update_latest_pointer('latest_checkpoint_2d', fname, 'checkpoint_2d', lok, lmsg)
      if (.not. lok) local_ok = .false.
    end if
    if (agree_outcome(local_ok, 'checkpoint_2d: failed to update latest_checkpoint_2d', &
                      is_ok, message)) return

    if (present(written_file)) written_file = trim(fname)
  end subroutine write_checkpoint_2d

  !> Read a checkpoint and restore state. Every rank reads the global file
  !! read-only, validates magic/version/dim/nx/ny/neq/grid/halo-width, then
  !! slices its interior block and its owned piece of the GLOBAL
  !! boundary-ghost frame (v6; see the write-side record-block comment).
  !!
  !! @param[out] t_comp  Optional: restored Kahan compensation term for `t`
  !!   (v6, Q2d); always read from the file when present, 0 on any failure.
  subroutine read_checkpoint_2d(state, fname, t, iter, is_ok, message, t_comp)
    type(solver_state_2d_t), intent(inout) :: state
    character(len=*), intent(in) :: fname
    real(wp), intent(out) :: t
    integer, intent(out) :: iter
    logical, intent(out), optional :: is_ok
    character(len=*), intent(out), optional :: message
    real(wp), intent(out), optional :: t_comp
    character(len=512) :: actual_file
    integer :: u, info, nx_g, ny_g, ix1, iy1, nxl, nyl, h, h_ck, dstat
    real(wp) :: dt_ck, cur_fp(2)
    real(wp) :: t_comp_val
    logical :: ok, local_ok, opened
    character(len=256) :: err
    real(wp), allocatable :: glob(:, :, :)
    real(wp), allocatable :: gh_l(:, :, :), gh_r(:, :, :), gh_b(:, :, :), gh_t(:, :, :)

    if (present(is_ok)) is_ok = .true.
    if (present(message)) message = ''
    if (present(t_comp)) t_comp = 0.0_wp
    iter = 0; t = 0.0_wp; t_comp_val = 0.0_wp; local_ok = .true.; opened = .false.; err = ''
    nx_g = state % decomp_2d % nx_global
    ny_g = state % decomp_2d % ny_global
    ix1 = state % decomp_2d % ix_first_global
    iy1 = state % decomp_2d % iy_first_global
    nxl = state % nx_local; nyl = state % ny_local
    h = state % halo_width
    cur_fp = grid_fingerprint_2d(state)   ! all ranks (uses reductions; call before any early skip)

    ! v6 (Q2d): allocate the ghost-frame restart buffers UP FRONT, sized by
    ! this rank's own expected halo width and globals (uniform across ranks
    ! by construction), and zero them, so the frame delivery below is always
    ! reachable at its fixed position regardless of any per-rank read failure
    ! (mirrors 1D's v5 edge-ghost discipline).
    allocate (gh_l(neq2d, h, ny_g), gh_r(neq2d, h, ny_g), &
              gh_b(neq2d, nx_g + 2 * h, h), gh_t(neq2d, nx_g + 2 * h, h), stat=info)
    if (info /= 0) call parallel_fatal('checkpoint_2d: cannot allocate ghost-frame restart buffers')
    gh_l = 0.0_wp; gh_r = 0.0_wp; gh_b = 0.0_wp; gh_t = 0.0_wp

    actual_file = resolve_latest_pointer('latest_checkpoint_2d', fname, 'checkpoint_2d', ok, err)
    if (.not. ok) then
      local_ok = .false.
    else
      open (newunit=u, file=trim(actual_file), status='old', action='read', &
            form='unformatted', access='stream', iostat=info)
      if (info /= 0) then
        err = 'checkpoint_2d: cannot open restart file "'//trim(actual_file)//'"'
        local_ok = .false.
      else
        opened = .true.
      end if
    end if

    if (local_ok) then
      call read_ckpt_header(u, 'checkpoint_2d', 2, [nx_g, ny_g], neq2d, iter, t, dt_ck, local_ok, err)
      if (local_ok) state % dt = dt_ck
    end if

    if (local_ok) then
      call read_grid_record_and_check(u, 'checkpoint_2d', state % cfg % grid_type, cur_fp, local_ok, err)
    end if

    if (local_ok) then
      allocate (glob(neq2d, nx_g, ny_g), stat=info)
      if (info /= 0) then; err = 'checkpoint_2d: cannot allocate restart buffer'; local_ok = .false.; end if
    end if
    if (local_ok) then
      read (u, iostat=info) glob
      if (info /= 0) then; err = 'checkpoint_2d: truncated solution'; local_ok = .false.
      else
        state % ub(:, 1:nxl, 1:nyl) = glob(:, ix1:ix1 + nxl - 1, iy1:iy1 + nyl - 1)
      end if
    end if
    if (allocated(glob)) then
      deallocate (glob, stat=dstat)
      if (dstat /= 0) local_ok = .false.
    end if

    ! v6 (Q2d): trailing Kahan t_comp scalar, mirroring 1D's v5 F2 record
    ! (checkpoint.f90). v6 always has it, so read it unconditionally.
    if (local_ok) then
      read (u, iostat=info) t_comp_val
      if (info /= 0) then; err = 'checkpoint_2d: truncated t_comp'; local_ok = .false.; end if
    end if

    ! v6 (Q2d): halo width (validated against this rank's current
    ! reconstruction scheme), then the GLOBAL boundary-ghost frame bands, in
    ! 1D's v5 record order. Every rank reads them from its own handle on the
    ! global file (the module's every-rank-reads contract); WHY they exist is
    ! documented at the write-side record block.
    if (local_ok) then
      read (u, iostat=info) h_ck
      if (info /= 0) then
        err = 'checkpoint_2d: truncated halo width'
        local_ok = .false.
      else if (h_ck /= h) then
        err = 'checkpoint_2d: halo width mismatch — checkpoint written with a different reconstruction scheme'
        local_ok = .false.
      end if
    end if
    if (local_ok) then
      read (u, iostat=info) gh_l
      if (info == 0) read (u, iostat=info) gh_r
      if (info == 0) read (u, iostat=info) gh_b
      if (info == 0) read (u, iostat=info) gh_t
      if (info /= 0) then; err = 'checkpoint_2d: truncated boundary ghost frame'; local_ok = .false.; end if
    end if

    if (opened) then
      close (u, iostat=info)
      if (info /= 0) call log_warn('checkpoint_2d: close failed')
    end if

    ! Fixed frame-delivery position reached by every rank UNCONDITIONALLY
    ! (see the up-front allocation above) — delivers the frame bands read
    ! above into this rank's ub ghosts (zeros, harmlessly, if a preceding
    ! read/validation failed on this rank; the collective conclude_read
    ! below aborts the whole restart in that case).
    call scatter_boundary_ghosts_from_root(state, gh_l, gh_r, gh_b, gh_t)
    deallocate (gh_l, gh_r, gh_b, gh_t, stat=dstat)
    if (dstat /= 0) then
      local_ok = .false.
      err = 'checkpoint_2d: ghost-frame buffer deallocation failed'
    end if

    if (present(t_comp) .and. local_ok) t_comp = t_comp_val

    ! Single collective decision point: every rank proceeds or fails together,
    ! so a failure on one rank cannot leave survivors marching into a collective
    ! the aborted ranks never enter (no restart deadlock).
    if (conclude_read(local_ok, err, 'checkpoint_2d: restart failed on at least one rank', &
                      is_ok, message)) return

    if (cart_rank_2d() == 0) call log_info('checkpoint_2d: resumed from "'//trim(actual_file)//'"')
  end subroutine read_checkpoint_2d

end module checkpoint_2d
