RUN_PART_A Subroutine

subroutine RUN_PART_A(N_FAIL)

Part A: DECOMPOSE FIELD_KNOWN_TOR, confirm exact recovery of the eight known unit a_tor coefficients and near-zero everywhere else, and confirm a_rad/a_pol are near-zero at every mode (the direct purity check). Writes examples/vsh_decomposition_tor/known_modes.dat.

a_tor coefficient error or the max a_rad/a_pol purity leakage exceeds TOL_A.

Arguments

Type IntentOptional Attributes Name
integer(kind=i4), intent(inout) :: N_FAIL

Running failure counter, incremented if either the max

Calls

proc~~run_part_a~2~~CallsGraph proc~run_part_a~2 RUN_PART_A proc~decompose~2 DECOMPOSE proc~run_part_a~2->proc~decompose~2 proc~ylm_index YLM_INDEX proc~run_part_a~2->proc~ylm_index dacos dacos proc~decompose~2->dacos interface~dot DOT proc~decompose~2->interface~dot interface~pvsh_pol_all PVSH_POL_ALL proc~decompose~2->interface~pvsh_pol_all interface~pvsh_rad_all PVSH_RAD_ALL proc~decompose~2->interface~pvsh_rad_all interface~pvsh_tor_all PVSH_TOR_ALL proc~decompose~2->interface~pvsh_tor_all interface~shglq SHGLQ proc~decompose~2->interface~shglq proc~dot_dp DOT_DP interface~dot->proc~dot_dp proc~pvsh_pol_all_dp PVSH_POL_ALL_DP interface~pvsh_pol_all->proc~pvsh_pol_all_dp proc~pvsh_rad_all_dp PVSH_RAD_ALL_DP interface~pvsh_rad_all->proc~pvsh_rad_all_dp proc~pvsh_tor_all_dp PVSH_TOR_ALL_DP interface~pvsh_tor_all->proc~pvsh_tor_all_dp proc~shglq_dp SHGLQ_DP interface~shglq->proc~shglq_dp proc~pvsh_pol_all_dp->proc~ylm_index interface~vsh_core VSH_CORE proc~pvsh_pol_all_dp->interface~vsh_core proc~pvsh_rad_all_dp->interface~vsh_core proc~pvsh_tor_all_dp->proc~ylm_index proc~pvsh_tor_all_dp->interface~vsh_core dcos dcos proc~shglq_dp->dcos proc~vsh_core_dp VSH_CORE_DP interface~vsh_core->proc~vsh_core_dp proc~vsh_core_dp->proc~ylm_index proc~vsh_core_dp->dcos dsin dsin proc~vsh_core_dp->dsin interface~assoc_legendre_norm_all ASSOC_LEGENDRE_NORM_ALL proc~vsh_core_dp->interface~assoc_legendre_norm_all interface~ddx_assoc_legendre_norm_all DDX_ASSOC_LEGENDRE_NORM_ALL proc~vsh_core_dp->interface~ddx_assoc_legendre_norm_all proc~plm_index PLM_INDEX proc~vsh_core_dp->proc~plm_index proc~assoc_legendre_norm_all_dp ASSOC_LEGENDRE_NORM_ALL_DP interface~assoc_legendre_norm_all->proc~assoc_legendre_norm_all_dp proc~ddx_assoc_legendre_norm_all_dp DDX_ASSOC_LEGENDRE_NORM_ALL_DP interface~ddx_assoc_legendre_norm_all->proc~ddx_assoc_legendre_norm_all_dp proc~assoc_legendre_norm_all_dp->proc~plm_index proc~ddx_assoc_legendre_norm_all_dp->proc~plm_index

Called by

proc~~run_part_a~2~~CalledByGraph proc~run_part_a~2 RUN_PART_A program~vsh_decomposition_tor VSH_DECOMPOSITION_TOR program~vsh_decomposition_tor->proc~run_part_a~2

Source Code

  SUBROUTINE RUN_PART_A(N_FAIL)
  IMPLICIT NONE
  INTEGER(KIND=i4), INTENT(INOUT) :: N_FAIL
  INTEGER(KIND=i4), PARAMETER :: LMAX_QUAD_A = 30
  INTEGER(KIND=i4), PARAMETER :: NPHI_A = 32
  COMPLEX(KIND=dp) :: A_RAD((LMAX_A+1)**2), A_POL((LMAX_A+1)**2), &
                       A_TOR((LMAX_A+1)**2)
  REAL(KIND=dp) :: EXPECTED, MAX_ERR, MAX_PURITY
  INTEGER(KIND=i4) :: L, M, IDX, K

  CALL DECOMPOSE(FIELD_KNOWN_TOR, LMAX_A, LMAX_QUAD_A, NPHI_A, &
                 A_RAD, A_POL, A_TOR)

  OPEN(UNIT=73, FILE="./examples/vsh_decomposition_tor/known_modes.dat")
  WRITE(73,'(A)') '# L  M  Re(a_rad)  Im(a_rad)  Re(a_pol)  Im(a_pol)  '// &
                  'Re(a_tor)  Im(a_tor)  expected'
  MAX_ERR = 0.0_dp
  MAX_PURITY = 0.0_dp
  DO L = 0, LMAX_A
    DO M = -L, L
      IDX = YLM_INDEX(L, M)
      EXPECTED = 0.0_dp
      DO K = 1, N_KNOWN
        IF (L_KNOWN(K) == L .AND. M_KNOWN(K) == M) EXPECTED = 1.0_dp
      ENDDO
      MAX_ERR = MAX(MAX_ERR, ABS(A_TOR(IDX)-EXPECTED))
      MAX_PURITY = MAX(MAX_PURITY, ABS(A_RAD(IDX)), ABS(A_POL(IDX)))
      WRITE(73,*) L, M, DBLE(A_RAD(IDX)), AIMAG(A_RAD(IDX)), &
                  DBLE(A_POL(IDX)), AIMAG(A_POL(IDX)), &
                  DBLE(A_TOR(IDX)), AIMAG(A_TOR(IDX)), EXPECTED
    ENDDO
  ENDDO
  CLOSE(73)

  IF (MAX_ERR > TOL_A) THEN
    N_FAIL = N_FAIL + 1
    WRITE(*,'(A,ES10.3)') "FAIL  known_mode_recovery  max_err=", MAX_ERR
  ELSE
    WRITE(*,'(A,ES10.3)') "PASS  known_mode_recovery  max_err=", MAX_ERR
  END IF

  IF (MAX_PURITY > TOL_A) THEN
    N_FAIL = N_FAIL + 1
    WRITE(*,'(A,ES10.3)') "FAIL  toroidal_purity  max_leakage=", MAX_PURITY
  ELSE
    WRITE(*,'(A,ES10.3)') "PASS  toroidal_purity  max_leakage=", MAX_PURITY
  END IF
  END SUBROUTINE RUN_PART_A