redi_interpolate_position Function

public pure function redi_interpolate_position(dRhoNeg, Pneg, dRhoPos, Ppos) result(pos)

Non-dimensional position in [0,1] where the interpolated density difference is zero. Guards the vanished/inverted (Ppos==Pneg) and degenerate (dRhoPos==dRhoNeg) cases device-safely (clamped values, no host I/O).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: dRhoNeg

negative density difference
real(kind=wp), intent(in) :: Pneg

position of the negative difference
real(kind=wp), intent(in) :: dRhoPos

positive density difference
real(kind=wp), intent(in) :: Ppos

position of the positive difference

Return Value real(kind=wp)


Called by

proc~~redi_interpolate_position~~CalledByGraph proc~redi_interpolate_position redi_interpolate_position proc~redi_neutral_positions_continuous redi_neutral_positions_continuous proc~redi_neutral_positions_continuous->proc~redi_interpolate_position proc~redi_face_coeffs redi_face_coeffs proc~redi_face_coeffs->proc~redi_neutral_positions_continuous proc~redi_calc_coeffs_x redi_calc_coeffs_x proc~redi_calc_coeffs_x->proc~redi_face_coeffs proc~redi_calc_coeffs_y redi_calc_coeffs_y proc~redi_calc_coeffs_y->proc~redi_face_coeffs proc~redi_calc_coeffs redi_calc_coeffs proc~redi_calc_coeffs->proc~redi_calc_coeffs_x proc~redi_calc_coeffs->proc~redi_calc_coeffs_y proc~ocean_dyn_step_split ocean_dyn_step_split proc~ocean_dyn_step_split->proc~redi_calc_coeffs

Source Code

   pure function redi_interpolate_position(dRhoNeg, Pneg, dRhoPos, Ppos) result(pos)
      !$acc routine seq
      real(wp), intent(in) :: dRhoNeg  !! negative density difference
      real(wp), intent(in) :: Pneg     !! position of the negative difference
      real(wp), intent(in) :: dRhoPos  !! positive density difference
      real(wp), intent(in) :: Ppos     !! position of the positive difference
      real(wp) :: pos

      if ((Ppos > Pneg) .and. (dRhoPos - dRhoNeg >= 0.0_wp)) then
         if (dRhoPos - dRhoNeg > 0.0_wp) then
            pos = min(1.0_wp, max(0.0_wp, -dRhoNeg/(dRhoPos - dRhoNeg)))
         else  ! dRhoPos - dRhoNeg == 0
            if (dRhoNeg > 0.0_wp) then
               pos = 0.0_wp
            else if (dRhoNeg < 0.0_wp) then
               pos = 1.0_wp
            else
               pos = 0.5_wp
            end if
         end if
      else if (Ppos == Pneg) then  ! vanished or inverted layers
         pos = 0.5_wp
      else  ! (Ppos < Pneg) .or. (dRhoNeg > dRhoPos): MOM6 errors here; device-safe clamp
         pos = 0.5_wp
      end if
   end function redi_interpolate_position