apply_escape_reflect_periodic_event Subroutine

  subroutine apply_escape_reflect_periodic_event( &
    cfg, event, x, v, alive, escaped, status, redistribution_uniform &
    )
    type(sim_config), intent(in) :: cfg
    type(boundary_event_type), intent(in) :: event
    real(dp), intent(inout) :: x(3), v(3)
    logical, intent(inout) :: alive
    logical, intent(out) :: escaped
    integer(i32), intent(out) :: status
    real(dp), intent(in), optional :: redistribution_uniform(3)

    real(dp) :: x_work(3), v_work(3)
    integer(i32) :: axis, low_face, high_face, low_bc, high_bc
    logical :: alive_work, has_open, has_redistributed_reflect

    escaped = .false.
    status = boundary_event_ok
    if (.not. cfg%use_box .or. .not. valid_event_box(cfg)) then
      status = boundary_event_invalid_geometry
      return
    end if
    if (trim(cfg%open_boundary_model) /= 'escape') then
      status = boundary_event_invalid_geometry
      return
    end if
    if (.not. all(ieee_is_finite(x)) .or. .not. all(ieee_is_finite(v))) then
      status = boundary_event_invalid_geometry
      return
    end if
    if (.not. event%has_event) then
      if (event%face_mask /= 0_i32 .or. event%fraction /= 0.0_dp) status = boundary_event_invalid_geometry
      return
    end if
    if (.not. ieee_is_finite(event%fraction) .or. event%fraction < 0.0_dp .or. event%fraction > 1.0_dp) then
      status = boundary_event_invalid_geometry
      return
    end if
    if (event%face_mask <= 0_i32 .or. iand(event%face_mask, boundary_face_all) /= event%face_mask) then
      status = boundary_event_invalid_geometry
      return
    end if
    if (any(event%face_bc /= boundary_face_conditions(cfg))) then
      status = boundary_event_invalid_geometry
      return
    end if
    do axis = 1_i32, 3_i32
      low_face = 2_i32*axis - 1_i32
      high_face = 2_i32*axis
      if (btest(event%face_mask, low_face - 1_i32) .and. btest(event%face_mask, high_face - 1_i32)) then
        status = boundary_event_invalid_geometry
        return
      end if
    end do

    x_work = x
    v_work = v
    alive_work = alive
    has_open = .false.
    has_redistributed_reflect = .false.
    do axis = 1_i32, 3_i32
      low_face = 2_i32*axis - 1_i32
      high_face = 2_i32*axis
      low_bc = event%face_bc(low_face)
      high_bc = event%face_bc(high_face)
      if (btest(event%face_mask, low_face - 1_i32)) has_open = has_open .or. low_bc == bc_open
      if (btest(event%face_mask, high_face - 1_i32)) has_open = has_open .or. high_bc == bc_open
      if (btest(event%face_mask, low_face - 1_i32)) then
        has_redistributed_reflect = has_redistributed_reflect .or. low_bc == bc_redistributed_reflect
      end if
      if (btest(event%face_mask, high_face - 1_i32)) then
        has_redistributed_reflect = has_redistributed_reflect .or. high_bc == bc_redistributed_reflect
      end if
    end do
    if (has_redistributed_reflect .and. .not. has_open) then
      if (.not. present(redistribution_uniform)) then
        status = boundary_event_invalid_geometry
        return
      end if
      if (.not. valid_redistribution_uniform(redistribution_uniform)) then
        status = boundary_event_invalid_geometry
        return
      end if
    end if

    if (has_open) then
      alive_work = .false.
      escaped = .true.
    else
      do axis = 1_i32, 3_i32
        low_face = 2_i32*axis - 1_i32
        high_face = 2_i32*axis
        low_bc = event%face_bc(low_face)
        high_bc = event%face_bc(high_face)
        if (btest(event%face_mask, low_face - 1_i32)) then
          call apply_surviving_event_axis( &
            axis, low_bc, .false., cfg%box_min(axis), cfg%box_max(axis), x_work, v_work, status &
            )
        else if (btest(event%face_mask, high_face - 1_i32)) then
          call apply_surviving_event_axis( &
            axis, high_bc, .true., cfg%box_min(axis), cfg%box_max(axis), x_work, v_work, status &
            )
        end if
        if (status /= boundary_event_ok) return
      end do
      if (has_redistributed_reflect) then
        call redistribute_event_position(cfg, event%face_mask, redistribution_uniform, x_work)
      end if
    end if
    if (.not. all(ieee_is_finite(x_work)) .or. .not. all(ieee_is_finite(v_work))) then
      escaped = .false.
      status = boundary_event_invalid_geometry
      return
    end if

    x = x_work
    v = v_work
    alive = alive_work
  end subroutine apply_escape_reflect_periodic_event