ASSOC_LEGENDRE_NORM_ALL_DP Subroutine

private subroutine ASSOC_LEGENDRE_NORM_ALL_DP(PNORM, LMAX, X)

All -normalized associated Legendre functions , (i.e. the real-valued radial part of SSH, for ), using the Holmes & Featherstone (2002) modified forward-column recurrence. The normalization is folded directly into the recurrence coefficients, keeping every intermediate value rather than letting itself grow/shrink combinatorially like the unnormalized ASSOC_LEGENDRE_ALL/Bonnet recurrence does. This is the recurrence SSH_ALL and every batch VSH routine build on.

Warning

Numerical stability: verified safe for at any . Beyond that, safety depends on : this recurrence has a distinct failure mode near the classical turning point of the associated Legendre ODE () -- for within about 30-56 degrees of either pole, values can blow up by 100+ orders of magnitude once exceeds a -dependent onset (as low as in the worst case, degrees). Near the poles ( degrees) or equator ( degrees) it remains accurate at least to (untested beyond that). This supersedes an earlier, overly optimistic "stable to " claim that didn't account for the -dependence. Full methodology, per-degree onset table, and a (\ell,\theta) heatmap: STABILITY_FINDINGS.md at the repo root (regenerable via the opt-in VSH_BUILD_STABILITY CMake target).

PLM_INDEX(l,m).

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(out) :: PNORM((LMAX+1)*(LMAX+2)/2)
integer(kind=i4), intent(in) :: LMAX

Maximum degree, .
real(kind=dp), intent(in) :: X

Argument, (typically ).

Calls

proc~~assoc_legendre_norm_all_dp~~CallsGraph proc~assoc_legendre_norm_all_dp ASSOC_LEGENDRE_NORM_ALL_DP proc~plm_index PLM_INDEX proc~assoc_legendre_norm_all_dp->proc~plm_index

Called by

proc~~assoc_legendre_norm_all_dp~~CalledByGraph proc~assoc_legendre_norm_all_dp ASSOC_LEGENDRE_NORM_ALL_DP interface~assoc_legendre_norm_all ASSOC_LEGENDRE_NORM_ALL interface~assoc_legendre_norm_all->proc~assoc_legendre_norm_all_dp proc~bench_legendre BENCH_LEGENDRE proc~bench_legendre->interface~assoc_legendre_norm_all proc~ssh_all_dp SSH_ALL_DP proc~ssh_all_dp->interface~assoc_legendre_norm_all proc~vsh_core_dp VSH_CORE_DP proc~vsh_core_dp->interface~assoc_legendre_norm_all program~convergence_main CONVERGENCE_MAIN program~convergence_main->interface~assoc_legendre_norm_all interface~ssh_all SSH_ALL program~convergence_main->interface~ssh_all proc~batch_ssh_cons BATCH_SSH_CONS program~convergence_main->proc~batch_ssh_cons interface~vsh_pol_dn_all VSH_POL_DN_ALL program~convergence_main->interface~vsh_pol_dn_all interface~vsh_pol_up_all VSH_POL_UP_ALL program~convergence_main->interface~vsh_pol_up_all proc~batch_vsh_pol_dn_cons BATCH_VSH_POL_DN_CONS program~convergence_main->proc~batch_vsh_pol_dn_cons proc~batch_vsh_pol_up_cons BATCH_VSH_POL_UP_CONS program~convergence_main->proc~batch_vsh_pol_up_cons interface~vsh_tor_all VSH_TOR_ALL program~convergence_main->interface~vsh_tor_all proc~batch_vsh_tor_cons BATCH_VSH_TOR_CONS program~convergence_main->proc~batch_vsh_tor_cons program~stability_main STABILITY_MAIN program~stability_main->interface~assoc_legendre_norm_all interface~ssh_all->proc~ssh_all_dp interface~vsh_core VSH_CORE interface~vsh_core->proc~vsh_core_dp program~bench_main BENCH_MAIN program~bench_main->proc~bench_legendre proc~bench_ssh BENCH_SSH program~bench_main->proc~bench_ssh proc~bench_vsh_pol_dn BENCH_VSH_POL_DN program~bench_main->proc~bench_vsh_pol_dn proc~bench_vsh_pol_up BENCH_VSH_POL_UP program~bench_main->proc~bench_vsh_pol_up proc~bench_vsh_tor BENCH_VSH_TOR program~bench_main->proc~bench_vsh_tor proc~batch_ssh_cons->interface~ssh_all proc~bench_ssh->interface~ssh_all proc~grad_ssh_all_dp GRAD_SSH_ALL_DP proc~grad_ssh_all_dp->interface~vsh_core proc~l_ssh_all_dp L_SSH_ALL_DP proc~l_ssh_all_dp->interface~vsh_core proc~pvsh_pol_all_dp PVSH_POL_ALL_DP proc~pvsh_pol_all_dp->interface~vsh_core proc~pvsh_rad_all_dp PVSH_RAD_ALL_DP proc~pvsh_rad_all_dp->interface~vsh_core proc~pvsh_tor_all_dp PVSH_TOR_ALL_DP proc~pvsh_tor_all_dp->interface~vsh_core proc~vsh_pol_dn_all_dp VSH_POL_DN_ALL_DP proc~vsh_pol_dn_all_dp->interface~vsh_core proc~vsh_pol_up_all_dp VSH_POL_UP_ALL_DP proc~vsh_pol_up_all_dp->interface~vsh_core interface~grad_ssh_all GRAD_SSH_ALL interface~grad_ssh_all->proc~grad_ssh_all_dp interface~l_ssh_all L_SSH_ALL interface~l_ssh_all->proc~l_ssh_all_dp interface~pvsh_pol_all PVSH_POL_ALL interface~pvsh_pol_all->proc~pvsh_pol_all_dp interface~pvsh_rad_all PVSH_RAD_ALL interface~pvsh_rad_all->proc~pvsh_rad_all_dp interface~pvsh_tor_all PVSH_TOR_ALL interface~pvsh_tor_all->proc~pvsh_tor_all_dp interface~vsh_pol_dn_all->proc~vsh_pol_dn_all_dp interface~vsh_pol_up_all->proc~vsh_pol_up_all_dp program~main MAIN program~main->proc~batch_ssh_cons proc~batch_grad_ssh_cons BATCH_GRAD_SSH_CONS program~main->proc~batch_grad_ssh_cons proc~batch_l_ssh_cons BATCH_L_SSH_CONS program~main->proc~batch_l_ssh_cons proc~batch_pvsh_pol_cons BATCH_PVSH_POL_CONS program~main->proc~batch_pvsh_pol_cons proc~batch_pvsh_rad_cons BATCH_PVSH_RAD_CONS program~main->proc~batch_pvsh_rad_cons proc~batch_pvsh_tor_cons BATCH_PVSH_TOR_CONS program~main->proc~batch_pvsh_tor_cons program~main->proc~batch_vsh_pol_dn_cons program~main->proc~batch_vsh_pol_up_cons proc~pvsh_pol_tor_ortho PVSH_POL_TOR_ORTHO program~main->proc~pvsh_pol_tor_ortho proc~vsh_pol_inversion VSH_POL_INVERSION program~main->proc~vsh_pol_inversion program~main->proc~batch_vsh_tor_cons proc~batch_grad_ssh_cons->interface~grad_ssh_all proc~batch_l_ssh_cons->interface~l_ssh_all proc~batch_pvsh_pol_cons->interface~pvsh_pol_all proc~batch_pvsh_rad_cons->interface~pvsh_rad_all proc~batch_pvsh_tor_cons->interface~pvsh_tor_all proc~batch_vsh_pol_dn_cons->interface~vsh_pol_dn_all proc~batch_vsh_pol_up_cons->interface~vsh_pol_up_all proc~bench_vsh_pol_dn->interface~vsh_pol_dn_all proc~bench_vsh_pol_up->interface~vsh_pol_up_all proc~decompose DECOMPOSE proc~decompose->interface~pvsh_pol_all proc~decompose->interface~pvsh_rad_all proc~decompose->interface~pvsh_tor_all proc~decompose~2 DECOMPOSE proc~decompose~2->interface~pvsh_pol_all proc~decompose~2->interface~pvsh_rad_all proc~decompose~2->interface~pvsh_tor_all proc~pvsh_pol_tor_ortho->interface~pvsh_pol_all proc~pvsh_pol_tor_ortho->interface~pvsh_tor_all proc~reconstruct RECONSTRUCT proc~reconstruct->interface~pvsh_pol_all proc~reconstruct->interface~pvsh_rad_all proc~reconstruct->interface~pvsh_tor_all proc~reconstruct~2 RECONSTRUCT proc~reconstruct~2->interface~pvsh_pol_all proc~reconstruct~2->interface~pvsh_rad_all proc~reconstruct~2->interface~pvsh_tor_all proc~vsh_pol_inversion->interface~pvsh_pol_all proc~vsh_pol_inversion->interface~pvsh_rad_all proc~vsh_pol_inversion->interface~vsh_pol_dn_all proc~vsh_pol_inversion->interface~vsh_pol_up_all proc~vsh_tor_all_dp VSH_TOR_ALL_DP proc~vsh_tor_all_dp->interface~pvsh_tor_all interface~vsh_tor_all->proc~vsh_tor_all_dp proc~run_part_a RUN_PART_A proc~run_part_a->proc~decompose proc~run_part_a~2 RUN_PART_A proc~run_part_a~2->proc~decompose~2 proc~run_part_b RUN_PART_B proc~run_part_b->proc~decompose proc~run_part_b->proc~reconstruct proc~run_part_b~2 RUN_PART_B proc~run_part_b~2->proc~decompose~2 proc~run_part_b~2->proc~reconstruct~2 proc~batch_vsh_tor_cons->interface~vsh_tor_all proc~bench_vsh_tor->interface~vsh_tor_all program~vsh_decomposition VSH_DECOMPOSITION program~vsh_decomposition->proc~run_part_a program~vsh_decomposition->proc~run_part_b program~vsh_decomposition_tor VSH_DECOMPOSITION_TOR program~vsh_decomposition_tor->proc~run_part_a~2 program~vsh_decomposition_tor->proc~run_part_b~2

Source Code

  SUBROUTINE ASSOC_LEGENDRE_NORM_ALL_DP(PNORM, LMAX, X)
  IMPLICIT NONE
  INTEGER(KIND=i4), INTENT(IN)  :: LMAX
  REAL(KIND=dp),    INTENT(IN)  :: X
  REAL(KIND=dp),    INTENT(OUT) :: PNORM((LMAX+1)*(LMAX+2)/2)
  INTEGER(KIND=i4) :: L, M
  REAL(KIND=dp)    :: SINTH

  SINTH = SQRT((1.0_dp - X) * (1.0_dp + X))

  PNORM(PLM_INDEX(0, 0)) = 1.0_dp / DSQRT(4.0_dp*pi)

  DO M = 0, LMAX
    IF (M .GT. 0) &
      PNORM(PLM_INDEX(M, M)) = &
        -DSQRT((2*M+1.d0)/(2*M)) * SINTH * &
        PNORM(PLM_INDEX(M-1, M-1))
    IF (M .LT. LMAX) &
      PNORM(PLM_INDEX(M+1, M)) = &
        DSQRT(2*M+3.d0) * X * PNORM(PLM_INDEX(M, M))
    DO L = M+2, LMAX
      PNORM(PLM_INDEX(L, M)) = &
        DSQRT((4.d0*L**2-1.d0)/(L**2-M**2)) * X * &
        PNORM(PLM_INDEX(L-1, M)) - &
        DSQRT((2*L+1.d0)*(L-M-1.d0)*(L+M-1.d0)/ &
              ((2*L-3.d0)*(L**2-M**2))) * &
        PNORM(PLM_INDEX(L-2, M))
    END DO
  END DO
  END SUBROUTINE ASSOC_LEGENDRE_NORM_ALL_DP