Remove old radius_1d_rj_r and a_r_1d_rj_r set routines

hiromatsui · hiromatsui · commit 2dfe3ac66037 · 2024-12-20T14:42:44.000-08:00
diff --git a/src/Fortran_libraries/PARALLEL_src/CONST_SPH_GRID/copy_sph_1d_global_index.f90 b/src/Fortran_libraries/PARALLEL_src/CONST_SPH_GRID/copy_sph_1d_global_index.f90
@@ -200,16 +200,7 @@ subroutine copy_sph_1d_gl_idx_rj                                  &
         sph_rj%radius_1d_rj_r(i) = s3d_radius%radius_1d_gl(j)
       end do
 !
-!$omp parallel do private(i)
-      do i = 1, sph_rj%nidx_rj(1)
-        sph_rj%a_r_1d_rj_r(i) = one / sph_rj%radius_1d_rj_r(i)
-        sph_rj%ar_1d_rj(i,1) = sph_rj%a_r_1d_rj_r(i)
-        sph_rj%ar_1d_rj(i,2) = sph_rj%ar_1d_rj(i,1)                     &
-     &                        * sph_rj%a_r_1d_rj_r(i)
-        sph_rj%ar_1d_rj(i,3) = sph_rj%ar_1d_rj(i,2)                     &
-     &                        * sph_rj%a_r_1d_rj_r(i)
-      end do
-!$omp end parallel do
+      call set_sph_one_over_radius_rj(sph_rj)
 !
       end subroutine copy_sph_1d_gl_idx_rj
 !
diff --git a/src/Fortran_libraries/PARALLEL_src/SPH_SHELL_src/init_sph_trans.f90 b/src/Fortran_libraries/PARALLEL_src/SPH_SHELL_src/init_sph_trans.f90
@@ -115,7 +115,7 @@ subroutine initialize_legendre_trans                              &
      &   (sph%sph_rtm%nidx_rtm, sph%sph_rlm%nidx_rlm, idx_trns)
 !
       call radial_4_sph_trans                                           &
-     &   (sph%sph_rtp, sph%sph_rtm, sph%sph_rlm, sph%sph_rj)
+     &   (sph%sph_rtp, sph%sph_rtm, sph%sph_rlm)
       call set_mdx_rlm_rtm(sph%sph_params%l_truncation,                 &
      &    sph%sph_rtm%nidx_rtm, sph%sph_rlm%nidx_rlm,                   &
      &    sph%sph_rtm%idx_gl_1d_rtm_m, sph%sph_rlm%idx_gl_1d_rlm_j,     &
diff --git a/src/Fortran_libraries/PARALLEL_src/SPH_SHELL_src/parallel_load_data_4_sph.f90 b/src/Fortran_libraries/PARALLEL_src/SPH_SHELL_src/parallel_load_data_4_sph.f90
@@ -300,6 +300,7 @@ end subroutine sph_rj_index_flags_and_params
       subroutine set_radius_dat_sph_MHD                                 &
      &         (radial_rj_grp, sph_params, sph_rj)
 !
+      use m_spheric_constants
       use set_radius_4_sph_dynamo
 !
       type(group_data), intent(in) :: radial_rj_grp
@@ -315,9 +316,13 @@ subroutine set_radius_dat_sph_MHD                                 &
      &   (sph_rj%nidx_rj(1), sph_rj%radius_1d_rj_r, radial_rj_grp,      &
      &    sph_params%iflag_radial_grid, sph_params%nlayer_ICB,          &
      &    sph_params%nlayer_CMB, sph_params%nlayer_2_center,            &
-     &    sph_rj%ar_1d_rj, sph_rj%r_ele_rj, sph_rj%ar_ele_rj,           &
+     &    sph_rj%r_ele_rj, sph_rj%ar_ele_rj,                            &
      &    sph_params%radius_ICB, sph_params%radius_CMB,                 &
      &    sph_params%R_earth)
+      if(sph_params%iflag_radial_grid .eq. igrid_error) then
+        call calypso_mpi_abort(ierr_sph,                                &
+     &     'Numbedr of radial layers are 0 or negative!')
+      end if
 !
       end subroutine set_radius_dat_sph_MHD
 !
diff --git a/src/Fortran_libraries/PARALLEL_src/SPH_SHELL_src/set_params_sph_trans.f90 b/src/Fortran_libraries/PARALLEL_src/SPH_SHELL_src/set_params_sph_trans.f90
@@ -23,11 +23,10 @@
 !!     &         (nidx_rtm, idx_gl_1d_mphi, mn_rlm)
 !!      subroutine set_sin_theta_rtm                                    &
 !!     &         (nth_rtm, g_colat_rtm, asin_theta_1d_rtm)
-!!      subroutine radial_4_sph_trans(sph_rtp, sph_rtm, sph_rlm, sph_rj)
+!!      subroutine radial_4_sph_trans(sph_rtp, sph_rtm, sph_rlm)
 !!        type(sph_rtp_grid), intent(inout) :: sph_rtp
 !!        type(sph_rtm_grid), intent(inout) :: sph_rtm
 !!        type(sph_rlm_grid), intent(inout) :: sph_rlm
-!!        type(sph_rj_grid), intent(inout) :: sph_rj
 !!@endverbatim
 !
       module set_params_sph_trans
@@ -245,14 +244,13 @@ end subroutine set_sin_theta_rtp
 ! -----------------------------------------------------------------------
 ! -----------------------------------------------------------------------
 !
-      subroutine radial_4_sph_trans(sph_rtp, sph_rtm, sph_rlm, sph_rj)
+      subroutine radial_4_sph_trans(sph_rtp, sph_rtm, sph_rlm)
 !
       use t_spheric_parameter
 !
       type(sph_rtp_grid), intent(inout) :: sph_rtp
       type(sph_rtm_grid), intent(inout) :: sph_rtm
       type(sph_rlm_grid), intent(inout) :: sph_rlm
-      type(sph_rj_grid), intent(inout) :: sph_rj
 !
 !
 !$omp parallel workshare
@@ -262,8 +260,6 @@ subroutine radial_4_sph_trans(sph_rtp, sph_rtm, sph_rlm, sph_rj)
      &      = one / sph_rtm%radius_1d_rtm_r(1:sph_rtm%nidx_rtm(1))
       sph_rlm%a_r_1d_rlm_r(1:sph_rlm%nidx_rlm(1))                       &
      &      = one / sph_rlm%radius_1d_rlm_r(1:sph_rlm%nidx_rlm(1))
-      sph_rj%a_r_1d_rj_r(1:sph_rj%nidx_rj(1))                           &
-     &      = one / sph_rj%radius_1d_rj_r(1:sph_rj%nidx_rj(1))
 !$omp end parallel workshare
 !
       end subroutine radial_4_sph_trans
diff --git a/src/Fortran_libraries/SERIAL_src/SPH_SPECTR_src/copy_sph_rj_mode_4_IO.f90 b/src/Fortran_libraries/SERIAL_src/SPH_SPECTR_src/copy_sph_rj_mode_4_IO.f90
@@ -66,25 +66,6 @@ subroutine copy_sph_node_4_rj_from_IO(sph_IO, sph_rj,             &
      &       = sph_IO%idx_gl_sph(1:sph_rj%nnod_rj,i)
       end do
 !
-!$omp parallel workshare
-      sph_rj%idx_gl_1d_rj_r(1:sph_rj%nidx_rj(1))                        &
-     &       = sph_IO%idx_gl_1(1:sph_rj%nidx_rj(1))
-!
-      sph_rj%radius_1d_rj_r(1:sph_rj%nidx_rj(1))                        &
-     &       = sph_IO%r_gl_1(1:sph_rj%nidx_rj(1))
-      sph_rj%a_r_1d_rj_r(1:sph_rj%nidx_rj(1))                           &
-     &       = one / sph_rj%radius_1d_rj_r(1:sph_rj%nidx_rj(1))
-!
-      sph_rj%ar_1d_rj(1:sph_rj%nidx_rj(1),1)                            &
-     &       = sph_rj%a_r_1d_rj_r(1:sph_rj%nidx_rj(1))
-      sph_rj%ar_1d_rj(1:sph_rj%nidx_rj(1),2)                            &
-     &       = sph_rj%ar_1d_rj(1:sph_rj%nidx_rj(1),1)                   &
-     &        * sph_rj%a_r_1d_rj_r(1:sph_rj%nidx_rj(1))
-      sph_rj%ar_1d_rj(1:sph_rj%nidx_rj(1),3)                            &
-     &       = sph_rj%ar_1d_rj(1:sph_rj%nidx_rj(1),2)                   &
-     &        * sph_rj%a_r_1d_rj_r(1:sph_rj%nidx_rj(1))
-!$omp end parallel workshare
-!
 !$omp parallel workshare
       sph_rj%idx_gl_1d_rj_j(1:sph_rj%nidx_rj(2),1)                      &
      &       = sph_IO%idx_gl_2(1:sph_rj%nidx_rj(2),1)
@@ -93,6 +74,16 @@ subroutine copy_sph_node_4_rj_from_IO(sph_IO, sph_rj,             &
       sph_rj%idx_gl_1d_rj_j(1:sph_rj%nidx_rj(2),3)                      &
      &      = sph_IO%idx_gl_2(1:sph_rj%nidx_rj(2),3)
 !$omp end parallel workshare
+!
+!$omp parallel workshare
+      sph_rj%idx_gl_1d_rj_r(1:sph_rj%nidx_rj(1))                        &
+     &       = sph_IO%idx_gl_1(1:sph_rj%nidx_rj(1))
+!
+      sph_rj%radius_1d_rj_r(1:sph_rj%nidx_rj(1))                        &
+     &       = sph_IO%r_gl_1(1:sph_rj%nidx_rj(1))
+!$omp end parallel workshare
+!
+      call set_sph_one_over_radius_rj(sph_rj)
 !
       end subroutine copy_sph_node_4_rj_from_IO
 !
diff --git a/src/Fortran_libraries/SERIAL_src/SPH_SPECTR_src/m_spheric_constants.f90 b/src/Fortran_libraries/SERIAL_src/SPH_SPECTR_src/m_spheric_constants.f90
@@ -21,6 +21,7 @@ module m_spheric_constants
 !
       use m_precision
 !
+      integer(kind = kint), parameter :: igrid_error = -999
       integer(kind = kint), parameter :: igrid_half_Chebyshev = 3
       integer(kind = kint), parameter :: igrid_Chebyshev =      2
       integer(kind = kint), parameter :: igrid_non_equidist =   1
diff --git a/src/Fortran_libraries/SERIAL_src/SPH_SPECTR_src/set_radial_grid_sph_shell.f90 b/src/Fortran_libraries/SERIAL_src/SPH_SPECTR_src/set_radial_grid_sph_shell.f90
@@ -114,24 +114,32 @@ subroutine set_radial_distance_flag(num_layer, nlayer_ICB,        &
 !
       real(kind = kreal), allocatable :: r_eq(:), r_ch(:), r_hch(:)
 !
+      if(num_layer .le. 0) then
+        iflag_rgrid = igrid_error
+        return
+      end if
 !
       allocate( r_eq(num_layer) )
       allocate( r_ch(num_layer) )
       allocate( r_hch(num_layer) )
 !
-      call set_equi_distance_shell(num_layer, nlayer_ICB, nlayer_CMB,   &
-     &    r_ICB, r_CMB, r_eq)
+      r_eq(1:num_layer) =  0.0d0
+      r_ch(1:num_layer) =  0.0d0
+      r_hch(1:num_layer) = 0.0d0
+!
+      call set_equi_distance_shell(num_layer, nlayer_ICB,               &
+     &                             nlayer_CMB, r_ICB, r_CMB, r_eq)
       call set_chebyshev_distance_shell(num_layer, nlayer_ICB,          &
-     &    nlayer_CMB, r_ICB, r_CMB, r_ch)
+     &                                  nlayer_CMB, r_ICB, r_CMB, r_ch)
       call half_chebyshev_distance_shell(num_layer, nlayer_CMB,         &
-     &    r_CMB, r_hch)
+     &                                   r_CMB, r_hch)
 !
 !
-      diff_eq_max =  abs( r_grid(1) - r_eq(1)) /  r_eq(1)
-      diff_ch_max =  abs( r_grid(1) - r_ch(1)) /  r_ch(1)
-      diff_hch_max = abs( r_grid(1) - r_hch(1)) / r_hch(1)
+      diff_eq_max =  0.0d0
+      diff_ch_max =  0.0d0
+      diff_hch_max = 0.0d0
 !
-      do k = 2, num_layer
+      do k = 1, num_layer
         diff = abs( r_grid(k) - r_eq(k)) / r_eq(k)
         diff_eq_max = max(diff_eq_max,diff)
 !
diff --git a/src/Fortran_libraries/SERIAL_src/SPH_SPECTR_src/set_radius_4_sph_dynamo.f90 b/src/Fortran_libraries/SERIAL_src/SPH_SPECTR_src/set_radius_4_sph_dynamo.f90
@@ -13,15 +13,19 @@
 !!      subroutine set_radius_dat_4_sph_dynamo                          &
 !!     &         (nri, radius_1d_rj_r, radial_rj_grp, iflag_radial_grid,&
 !!     &          nlayer_ICB, nlayer_CMB, nlayer_2_center,              &
-!!     &          ar_1d_rj, r_ele_rj, ar_ele_rj, r_ICB, r_CMB, R_earth)
+!!     &          r_ele_rj, ar_ele_rj, r_ICB, r_CMB, R_earth)
 !!        type(group_data), intent(in) :: radial_rj_grp
-!!***********************************************************************
-!!*
-!!*       ar_1d_rj(k,1)   : 1 / r
-!!*       ar_1d_rj(k,2)   : 1 / r**2
-!!*       ar_1d_rj(k,3)   : 1 / r**3
-!!*
-!!***********************************************************************
+!!        integer(kind = kint), intent(in) :: nri
+!!        real(kind = kreal), intent(in) :: radius_1d_rj_r(nri)
+!!
+!!        integer(kind = kint), intent(inout) :: iflag_radial_grid
+!!        integer(kind = kint), intent(inout) :: nlayer_ICB
+!!        integer(kind = kint), intent(inout) :: nlayer_CMB
+!!        integer(kind = kint), intent(inout) :: nlayer_2_center
+!!        real(kind = kreal), intent(inout) :: r_ele_rj(nri)
+!!        real(kind = kreal), intent(inout) :: ar_ele_rj(nri,3)
+!!        real(kind = kreal), intent(inout) :: R_earth(0:2)
+!!        real(kind = kreal), intent(inout) :: r_ICB, r_CMB
 !!@endverbatim
 !
 !
@@ -47,7 +51,7 @@ module set_radius_4_sph_dynamo
       subroutine set_radius_dat_4_sph_dynamo                            &
      &         (nri, radius_1d_rj_r, radial_rj_grp, iflag_radial_grid,  &
      &          nlayer_ICB, nlayer_CMB, nlayer_2_center,                &
-     &          ar_1d_rj, r_ele_rj, ar_ele_rj, r_ICB, r_CMB, R_earth)
+     &          r_ele_rj, ar_ele_rj, r_ICB, r_CMB, R_earth)
 !
       use set_radial_grid_sph_shell
       use skip_comment_f
@@ -61,7 +65,6 @@ subroutine set_radius_dat_4_sph_dynamo                            &
       integer(kind = kint), intent(inout) :: nlayer_ICB
       integer(kind = kint), intent(inout) :: nlayer_CMB
       integer(kind = kint), intent(inout) :: nlayer_2_center
-      real(kind = kreal), intent(inout) :: ar_1d_rj(nri,3)
       real(kind = kreal), intent(inout) :: r_ele_rj(nri)
       real(kind = kreal), intent(inout) :: ar_ele_rj(nri,3)
       real(kind = kreal), intent(inout) :: R_earth(0:2)
@@ -97,13 +100,8 @@ subroutine set_radius_dat_4_sph_dynamo                            &
 !
       call set_radial_distance_flag(nri, nlayer_ICB, nlayer_CMB,        &
      &    r_ICB, r_CMB, radius_1d_rj_r, iflag_radial_grid)
+      if(iflag_radial_grid .eq. igrid_error) return
 !
-!
-      do k = 1, nri
-        ar_1d_rj(k,1) = one / radius_1d_rj_r(k)
-        ar_1d_rj(k,2) = ar_1d_rj(k,1)**2
-        ar_1d_rj(k,3) = ar_1d_rj(k,1)**3
-      end do
 !
       r_ele_rj(1) = half * radius_1d_rj_r(1)
       do k = 2, nri
diff --git a/src/Fortran_libraries/SERIAL_src/SPH_SPECTR_src/t_spheric_rj_data.f90 b/src/Fortran_libraries/SERIAL_src/SPH_SPECTR_src/t_spheric_rj_data.f90
@@ -17,6 +17,8 @@
 !!      subroutine dealloc_rj_param_smp(sph_rj)
 !!        type(sph_rj_grid), intent(inout) :: sph_rj
 !!
+!!      subroutine set_sph_one_over_radius_rj(sph_rj)
+!!        type(sph_rj_grid), intent(inout) :: sph_rj
 !!      subroutine copy_spheric_rj_data                                 &
 !!     &         (ltr_org, rj_org, ltr_new, rj_new)
 !!        type(sph_rj_grid), intent(in) :: rj_org
@@ -228,6 +230,27 @@ end subroutine dealloc_rj_param_smp
 !
 ! -----------------------------------------------------------------------
 ! -----------------------------------------------------------------------
+!
+      subroutine set_sph_one_over_radius_rj(sph_rj)
+!
+      type(sph_rj_grid), intent(inout) :: sph_rj
+!
+      integer(kind = kint) :: i
+!
+!$omp parallel do private(i)
+      do i = 1, sph_rj%nidx_rj(1)
+        sph_rj%a_r_1d_rj_r(i) = one / sph_rj%radius_1d_rj_r(i)
+        sph_rj%ar_1d_rj(i,1) = sph_rj%a_r_1d_rj_r(i)
+        sph_rj%ar_1d_rj(i,2) = sph_rj%ar_1d_rj(i,1)                     &
+     &                        * sph_rj%a_r_1d_rj_r(i)
+        sph_rj%ar_1d_rj(i,3) = sph_rj%ar_1d_rj(i,2)                     &
+     &                        * sph_rj%a_r_1d_rj_r(i)
+      end do
+!$omp end parallel do
+!
+      end subroutine set_sph_one_over_radius_rj
+!
+! -----------------------------------------------------------------------
 !
       subroutine copy_spheric_rj_data                                   &
      &         (ltr_org, rj_org, ltr_new, rj_new)
diff --git a/src/programs/data_utilities/Rayleigh_link/t_convert_from_rayleigh.f90 b/src/programs/data_utilities/Rayleigh_link/t_convert_from_rayleigh.f90
@@ -192,17 +192,8 @@ subroutine copy_rayleigh_radial_data(ra_rst, org_sph)
         org_sph%sph_rj%radius_1d_rj_r(k) = ra_rst%r_org(kr)
       end do
 !$omp end parallel do
-!$omp parallel do private(k)
-      do k = 1, org_sph%sph_rj%nidx_rj(1)
-        org_sph%sph_rj%a_r_1d_rj_r(k)                                   &
-     &         = one / org_sph%sph_rj%radius_1d_rj_r(k)
-        org_sph%sph_rj%ar_1d_rj(k,1) = org_sph%sph_rj%a_r_1d_rj_r(k)
-        org_sph%sph_rj%ar_1d_rj(k,2) = org_sph%sph_rj%ar_1d_rj(k,1)     &
-     &                                * org_sph%sph_rj%a_r_1d_rj_r(k)
-        org_sph%sph_rj%ar_1d_rj(k,3) = org_sph%sph_rj%ar_1d_rj(k,2)     &
-     &                                * org_sph%sph_rj%a_r_1d_rj_r(k)
-      end do
-!$omp end parallel do
+!
+      call set_sph_one_over_radius_rj(sph_rj)
 !
       end subroutine copy_rayleigh_radial_data
 !

Original file line number	Diff line number	Diff line change
`@@ -200,16 +200,7 @@ subroutine copy_sph_1d_gl_idx_rj &`
`200`	`200`	`sph_rj%radius_1d_rj_r(i) = s3d_radius%radius_1d_gl(j)`
`201`	`201`	`end do`
`202`	`202`	`!`
`203`		`-!$omp parallel do private(i)`
`204`		`- do i = 1, sph_rj%nidx_rj(1)`
`205`		`- sph_rj%a_r_1d_rj_r(i) = one / sph_rj%radius_1d_rj_r(i)`
`206`		`- sph_rj%ar_1d_rj(i,1) = sph_rj%a_r_1d_rj_r(i)`
`207`		`- sph_rj%ar_1d_rj(i,2) = sph_rj%ar_1d_rj(i,1) &`
`208`		`- & * sph_rj%a_r_1d_rj_r(i)`
`209`		`- sph_rj%ar_1d_rj(i,3) = sph_rj%ar_1d_rj(i,2) &`
`210`		`- & * sph_rj%a_r_1d_rj_r(i)`
`211`		`- end do`
`212`		`-!$omp end parallel do`
	`203`	`+ call set_sph_one_over_radius_rj(sph_rj)`
`213`	`204`	`!`
`214`	`205`	`end subroutine copy_sph_1d_gl_idx_rj`
`215`	`206`	`!`
Original file line number	Diff line number	Diff line change
`@@ -115,7 +115,7 @@ subroutine initialize_legendre_trans &`
`115`	`115`	`& (sph%sph_rtm%nidx_rtm, sph%sph_rlm%nidx_rlm, idx_trns)`
`116`	`116`	`!`
`117`	`117`	`call radial_4_sph_trans &`
`118`		`- & (sph%sph_rtp, sph%sph_rtm, sph%sph_rlm, sph%sph_rj)`
	`118`	`+ & (sph%sph_rtp, sph%sph_rtm, sph%sph_rlm)`
`119`	`119`	`call set_mdx_rlm_rtm(sph%sph_params%l_truncation, &`
`120`	`120`	`& sph%sph_rtm%nidx_rtm, sph%sph_rlm%nidx_rlm, &`
`121`	`121`	`& sph%sph_rtm%idx_gl_1d_rtm_m, sph%sph_rlm%idx_gl_1d_rlm_j, &`
Original file line number	Diff line number	Diff line change
`@@ -300,6 +300,7 @@ end subroutine sph_rj_index_flags_and_params`
`300`	`300`	`subroutine set_radius_dat_sph_MHD &`
`301`	`301`	`& (radial_rj_grp, sph_params, sph_rj)`
`302`	`302`	`!`
	`303`	`+ use m_spheric_constants`
`303`	`304`	`use set_radius_4_sph_dynamo`
`304`	`305`	`!`
`305`	`306`	`type(group_data), intent(in) :: radial_rj_grp`
`@@ -315,9 +316,13 @@ subroutine set_radius_dat_sph_MHD &`
`315`	`316`	`& (sph_rj%nidx_rj(1), sph_rj%radius_1d_rj_r, radial_rj_grp, &`
`316`	`317`	`& sph_params%iflag_radial_grid, sph_params%nlayer_ICB, &`
`317`	`318`	`& sph_params%nlayer_CMB, sph_params%nlayer_2_center, &`
`318`		`- & sph_rj%ar_1d_rj, sph_rj%r_ele_rj, sph_rj%ar_ele_rj, &`
	`319`	`+ & sph_rj%r_ele_rj, sph_rj%ar_ele_rj, &`
`319`	`320`	`& sph_params%radius_ICB, sph_params%radius_CMB, &`
`320`	`321`	`& sph_params%R_earth)`
	`322`	`+ if(sph_params%iflag_radial_grid .eq. igrid_error) then`
	`323`	`+ call calypso_mpi_abort(ierr_sph, &`
	`324`	`+ & 'Numbedr of radial layers are 0 or negative!')`
	`325`	`+ end if`
`321`	`326`	`!`
`322`	`327`	`end subroutine set_radius_dat_sph_MHD`
`323`	`328`	`!`
Original file line number	Diff line number	Diff line change
`@@ -23,11 +23,10 @@`
`23`	`23`	`!! & (nidx_rtm, idx_gl_1d_mphi, mn_rlm)`
`24`	`24`	`!! subroutine set_sin_theta_rtm &`
`25`	`25`	`!! & (nth_rtm, g_colat_rtm, asin_theta_1d_rtm)`
`26`		`-!! subroutine radial_4_sph_trans(sph_rtp, sph_rtm, sph_rlm, sph_rj)`
	`26`	`+!! subroutine radial_4_sph_trans(sph_rtp, sph_rtm, sph_rlm)`
`27`	`27`	`!! type(sph_rtp_grid), intent(inout) :: sph_rtp`
`28`	`28`	`!! type(sph_rtm_grid), intent(inout) :: sph_rtm`
`29`	`29`	`!! type(sph_rlm_grid), intent(inout) :: sph_rlm`
`30`		`-!! type(sph_rj_grid), intent(inout) :: sph_rj`
`31`	`30`	`!!@endverbatim`
`32`	`31`	`!`
`33`	`32`	`module set_params_sph_trans`
`@@ -245,14 +244,13 @@ end subroutine set_sin_theta_rtp`
`245`	`244`	`! -----------------------------------------------------------------------`
`246`	`245`	`! -----------------------------------------------------------------------`
`247`	`246`	`!`
`248`		`- subroutine radial_4_sph_trans(sph_rtp, sph_rtm, sph_rlm, sph_rj)`
	`247`	`+ subroutine radial_4_sph_trans(sph_rtp, sph_rtm, sph_rlm)`
`249`	`248`	`!`
`250`	`249`	`use t_spheric_parameter`
`251`	`250`	`!`
`252`	`251`	`type(sph_rtp_grid), intent(inout) :: sph_rtp`
`253`	`252`	`type(sph_rtm_grid), intent(inout) :: sph_rtm`
`254`	`253`	`type(sph_rlm_grid), intent(inout) :: sph_rlm`
`255`		`- type(sph_rj_grid), intent(inout) :: sph_rj`
`256`	`254`	`!`
`257`	`255`	`!`
`258`	`256`	`!$omp parallel workshare`
`@@ -262,8 +260,6 @@ subroutine radial_4_sph_trans(sph_rtp, sph_rtm, sph_rlm, sph_rj)`
`262`	`260`	`& = one / sph_rtm%radius_1d_rtm_r(1:sph_rtm%nidx_rtm(1))`
`263`	`261`	`sph_rlm%a_r_1d_rlm_r(1:sph_rlm%nidx_rlm(1)) &`
`264`	`262`	`& = one / sph_rlm%radius_1d_rlm_r(1:sph_rlm%nidx_rlm(1))`
`265`		`- sph_rj%a_r_1d_rj_r(1:sph_rj%nidx_rj(1)) &`
`266`		`- & = one / sph_rj%radius_1d_rj_r(1:sph_rj%nidx_rj(1))`
`267`	`263`	`!$omp end parallel workshare`
`268`	`264`	`!`
`269`	`265`	`end subroutine radial_4_sph_trans`
Original file line number	Diff line number	Diff line change
`@@ -21,6 +21,7 @@ module m_spheric_constants`
`21`	`21`	`!`
`22`	`22`	`use m_precision`
`23`	`23`	`!`
	`24`	`+ integer(kind = kint), parameter :: igrid_error = -999`
`24`	`25`	`integer(kind = kint), parameter :: igrid_half_Chebyshev = 3`
`25`	`26`	`integer(kind = kint), parameter :: igrid_Chebyshev = 2`
`26`	`27`	`integer(kind = kint), parameter :: igrid_non_equidist = 1`
Original file line number	Diff line number	Diff line change
`@@ -114,24 +114,32 @@ subroutine set_radial_distance_flag(num_layer, nlayer_ICB, &`
`114`	`114`	`!`
`115`	`115`	`real(kind = kreal), allocatable :: r_eq(:), r_ch(:), r_hch(:)`
`116`	`116`	`!`
	`117`	`+ if(num_layer .le. 0) then`
	`118`	`+ iflag_rgrid = igrid_error`
	`119`	`+ return`
	`120`	`+ end if`
`117`	`121`	`!`
`118`	`122`	`allocate( r_eq(num_layer) )`
`119`	`123`	`allocate( r_ch(num_layer) )`
`120`	`124`	`allocate( r_hch(num_layer) )`
`121`	`125`	`!`
`122`		`- call set_equi_distance_shell(num_layer, nlayer_ICB, nlayer_CMB, &`
`123`		`- & r_ICB, r_CMB, r_eq)`
	`126`	`+ r_eq(1:num_layer) = 0.0d0`
	`127`	`+ r_ch(1:num_layer) = 0.0d0`
	`128`	`+ r_hch(1:num_layer) = 0.0d0`
	`129`	`+!`
	`130`	`+ call set_equi_distance_shell(num_layer, nlayer_ICB, &`
	`131`	`+ & nlayer_CMB, r_ICB, r_CMB, r_eq)`
`124`	`132`	`call set_chebyshev_distance_shell(num_layer, nlayer_ICB, &`
`125`		`- & nlayer_CMB, r_ICB, r_CMB, r_ch)`
	`133`	`+ & nlayer_CMB, r_ICB, r_CMB, r_ch)`
`126`	`134`	`call half_chebyshev_distance_shell(num_layer, nlayer_CMB, &`
`127`		`- & r_CMB, r_hch)`
	`135`	`+ & r_CMB, r_hch)`
`128`	`136`	`!`
`129`	`137`	`!`
`130`		`- diff_eq_max = abs( r_grid(1) - r_eq(1)) / r_eq(1)`
`131`		`- diff_ch_max = abs( r_grid(1) - r_ch(1)) / r_ch(1)`
`132`		`- diff_hch_max = abs( r_grid(1) - r_hch(1)) / r_hch(1)`
	`138`	`+ diff_eq_max = 0.0d0`
	`139`	`+ diff_ch_max = 0.0d0`
	`140`	`+ diff_hch_max = 0.0d0`
`133`	`141`	`!`
`134`		`- do k = 2, num_layer`
	`142`	`+ do k = 1, num_layer`
`135`	`143`	`diff = abs( r_grid(k) - r_eq(k)) / r_eq(k)`
`136`	`144`	`diff_eq_max = max(diff_eq_max,diff)`
`137`	`145`	`!`