Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Enable E3SM shared pi in MPAS-Ocean #6485

Merged
merged 11 commits into from
Oct 23, 2024
Merged

Enable E3SM shared pi in MPAS-Ocean #6485

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
11 commits
Select commit Hold shift + click to select a range
6160a60
Add shr_const to analysis mode
cbegeman Jun 21, 2024
2da96ef
Use pi for debug tracer
cbegeman Jun 21, 2024
561ce77
Use pi in cvmix
cbegeman Jun 21, 2024
6223ff4
Use pi in gm
cbegeman Jun 21, 2024
1b92637
Set pi in mpas_ocn_constants
cbegeman Jun 21, 2024
5e7322e
Use pi in conservation earth area
cbegeman Jun 21, 2024
ba7b48e
Set pi in mpas_ocn_constants
cbegeman Jun 21, 2024
58a8d60
Add shr_const to analysis mode
cbegeman Jun 21, 2024
17aed5b
Use MPASO's pi instead of MPAS-framework's pi in AM
cbegeman Jul 1, 2024
5ce53c7
Use MPASO's pi instead of MPAS-framework's pi in shared
cbegeman Jul 1, 2024
b3c8342
Increase precision of MPAS-O standalone pi
cbegeman Jul 1, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
5 changes: 4 additions & 1 deletion components/mpas-ocean/src/analysis_members/Makefile
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,9 @@

OBJS = mpas_ocn_analysis_driver.o

UTILS = shr_kind_mod.o \
shr_const_mod.o

MEMBERS = mpas_ocn_global_stats.o \
mpas_ocn_okubo_weiss.o \
mpas_ocn_layer_volume_weighted_averages.o \
Expand Down Expand Up @@ -35,7 +38,7 @@ MEMBERS = mpas_ocn_global_stats.o \

all: $(OBJS)

mpas_ocn_analysis_driver.o: $(MEMBERS)
mpas_ocn_analysis_driver.o: $(UTILS) $(MEMBERS)

mpas_ocn_okubo_weiss.o: mpas_ocn_okubo_weiss_eigenvalues.o

Expand Down
8 changes: 5 additions & 3 deletions components/mpas-ocean/src/analysis_members/mpas_ocn_conservation_check.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,9 @@

module ocn_conservation_check

use shr_kind_mod, only: SHR_KIND_R8
use shr_const_mod

use mpas_derived_types
use mpas_pool_routines
use mpas_dmpar
Expand Down Expand Up @@ -157,10 +160,9 @@ subroutine ocn_init_conservation_check(domain, err)!{{{
!-----------------------------------------------------------------

! taking PI from SHR_CONST_PI in share/util/shr_const_mod.F90 to match coupler
!real (kind=RKIND), parameter :: piE3SM = 3.14159265358979323846_RKIND ! pi
real (kind=RKIND), parameter :: piE3SM = 3.141592653589793_RKIND
real (kind=RKIND), parameter :: piE3SM = SHR_CONST_PI
! taking earth radius from SHR_CONST_REARTH in share/util/shr_const_mod.F90 to match coupler
real (kind=RKIND), parameter :: earthRadiusE3SM = 6.371229e6_RKIND ! radius of earth, m
real (kind=RKIND), parameter :: earthRadiusE3SM = SHR_CONST_REARTH ! radius of earth, m

err = 0

Expand Down
2 changes: 1 addition & 1 deletion components/mpas-ocean/src/analysis_members/mpas_ocn_harmonic_analysis.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -681,7 +681,7 @@ SUBROUTINE harmonic_analysis_solve(MNP,nfreq,hmat,GLOELV,haff,haface,emagt,phase
REAL(kind=RKIND),ALLOCATABLE :: hap(:),hax(:)
REAL(kind=RKIND),ALLOCATABLE :: ha(:,:)

convrd=180.0_RKIND/pii
convrd=180.0_RKIND/pi

mm = 2*nfreq
ALLOCATE ( PHASEE(nfreq),EMAG(nfreq) )
Expand Down
2 changes: 1 addition & 1 deletion components/mpas-ocean/src/analysis_members/mpas_ocn_lagrangian_particle_tracking.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3833,7 +3833,7 @@ subroutine convert_latlon_from_xyz(lat, lon, x, y, z) !{{{

! ensure range in 0, 2*pi
if (lon < 0.0_RKIND) then
lon = 2*pii + lon
lon = 2*pi + lon
end if

end subroutine convert_latlon_from_xyz !}}}
Expand Down
4 changes: 2 additions & 2 deletions components/mpas-ocean/src/analysis_members/mpas_ocn_okubo_weiss.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -946,8 +946,8 @@ subroutine ocn_compute_eddy_stats(dminfo, block, nVertLevels, nCells, nCellsSolv

! for lat/lon coordinates, convert from radians to degrees for output
if (config_AM_okuboWeiss_use_lat_lon_coords) then
wsPosX = wsPosX *180.0_RKIND/pii
wsPosY = wsPosY *180.0_RKIND/pii
wsPosX = wsPosX *180.0_RKIND/pi
wsPosY = wsPosY *180.0_RKIND/pi
end if

call mpas_timer_stop("stats per proc")
Expand Down
4 changes: 2 additions & 2 deletions components/mpas-ocean/src/analysis_members/mpas_ocn_time_filters.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -62,7 +62,7 @@ module ocn_time_filters
!--------------------------------------------------------------------
#ifdef TIME_FILTERS_DEBUG
integer :: iEdgeOutput = 0, iBlockOutput = 0, iklevel = 1
real (kind=RKIND) :: lonEdgePoint = 10.0_RKIND*pii/180.0_RKIND, latEdgePoint = 30.0_RKIND*pii/180.0_RKIND
real (kind=RKIND) :: lonEdgePoint = 10.0_RKIND*pi/180.0_RKIND, latEdgePoint = 30.0_RKIND*pi/180.0_RKIND
#endif

!***********************************************************************
Expand Down Expand Up @@ -182,7 +182,7 @@ subroutine ocn_init_time_filters(domain, err)!{{{
call mpas_pool_get_subpool(block % structs, 'mesh', statePool)
call mpas_pool_get_array(statePool, 'latEdge', latEdge)
call mpas_pool_get_array(statePool, 'lonEdge', lonEdge)
print *, 'lon = ', 180.0_RKIND/pii*lonEdge(iEdgeOutput), ' lat = ', 180.0_RKIND/pii*latEdge(iEdgeOutput), &
print *, 'lon = ', 180.0_RKIND/pi*lonEdge(iEdgeOutput), ' lat = ', 180.0_RKIND/pi*latEdge(iEdgeOutput), &
' iklevel=',iklevel, ' iEdgeOutput=',iEdgeOutput, ' iBlockOutput = ', iBlockOutput
#endif

Expand Down
1 change: 1 addition & 0 deletions components/mpas-ocean/src/analysis_members/shr_const_mod.F
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
../../../../share/util/shr_const_mod.F90
1 change: 1 addition & 0 deletions components/mpas-ocean/src/analysis_members/shr_kind_mod.F
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
../../../../share/util/shr_kind_mod.F90
4 changes: 2 additions & 2 deletions components/mpas-ocean/src/mode_forward/mpas_ocn_time_integration_si.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2587,7 +2587,7 @@ subroutine ocn_time_integrator_si(domain, dt)!{{{

! Reset tracer2 to 2 in top n layers
! in zonal bands, and 1 outside
lat = latCell(iCell)*180./3.1415
lat = latCell(iCell)*180.0_RKIND/pi
if ( lat>-60.0.and.lat<-55.0 &
.or.lat>-40.0.and.lat<-35.0 &
.or.lat>- 2.5.and.lat< 2.5 &
Expand All @@ -2604,7 +2604,7 @@ subroutine ocn_time_integrator_si(domain, dt)!{{{

! Reset tracer3 to 2 in top n layers
! in zonal bands, and 1 outside
lat = latCell(iCell)*180./3.1415
lat = latCell(iCell)*180.0_RKIND/pi
if ( lat>-55.0.and.lat<-50.0 &
.or.lat>-35.0.and.lat<-30.0 &
.or.lat>-15.0.and.lat<-10.0 &
Expand Down
4 changes: 2 additions & 2 deletions components/mpas-ocean/src/mode_forward/mpas_ocn_time_integration_split.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2225,7 +2225,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{

! Reset tracer2 to 2 in top n layers
! in zonal bands, and 1 outside
lat = latCell(iCell)*180./3.1415
lat = latCell(iCell)*180.0_RKIND/pi
if ( lat>-60.0.and.lat<-55.0 &
.or.lat>-40.0.and.lat<-35.0 &
.or.lat>- 2.5.and.lat< 2.5 &
Expand All @@ -2242,7 +2242,7 @@ subroutine ocn_time_integrator_split(domain, dt)!{{{

! Reset tracer3 to 2 in top n layers
! in zonal bands, and 1 outside
lat = latCell(iCell)*180./3.1415
lat = latCell(iCell)*180.0_RKIND/pi
if ( lat>-55.0.and.lat<-50.0 &
.or.lat>-35.0.and.lat<-30.0 &
.or.lat>-15.0.and.lat<-10.0 &
Expand Down
4 changes: 2 additions & 2 deletions components/mpas-ocean/src/mode_forward/mpas_ocn_time_integration_split_ab2.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2405,7 +2405,7 @@ subroutine ocn_time_integrator_split_ab2(domain, dt)!{{{

! Reset tracer2 to 2 in top n layers
! in zonal bands, and 1 outside
lat = latCell(iCell)*180./3.1415
lat = latCell(iCell)*180.0_RKIND/pi
if ( lat>-60.0.and.lat<-55.0 &
.or.lat>-40.0.and.lat<-35.0 &
.or.lat>- 2.5.and.lat< 2.5 &
Expand All @@ -2422,7 +2422,7 @@ subroutine ocn_time_integrator_split_ab2(domain, dt)!{{{

! Reset tracer3 to 2 in top n layers
! in zonal bands, and 1 outside
lat = latCell(iCell)*180./3.1415
lat = latCell(iCell)*180.0_RKIND/pi
if ( lat>-55.0.and.lat<-50.0 &
.or.lat>-35.0.and.lat<-30.0 &
.or.lat>-15.0.and.lat<-10.0 &
Expand Down
2 changes: 1 addition & 1 deletion components/mpas-ocean/src/shared/Makefile
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -143,7 +143,7 @@ mpas_ocn_tracer_advection_std.o: mpas_ocn_config.o mpas_ocn_mesh.o mpas_ocn_trac

mpas_ocn_tracer_advection_vert.o: mpas_ocn_mesh.o mpas_ocn_config.o

mpas_ocn_tracer_advection_shared.o: mpas_ocn_mesh.o mpas_ocn_config.o
mpas_ocn_tracer_advection_shared.o: mpas_ocn_constants.o mpas_ocn_mesh.o mpas_ocn_config.o

mpas_ocn_tracer_hmix_redi.o: mpas_ocn_constants.o mpas_ocn_config.o

Expand Down
3 changes: 3 additions & 0 deletions components/mpas-ocean/src/shared/mpas_ocn_constants.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -65,6 +65,7 @@ module ocn_constants

real (kind=RKIND), public :: &
T0_Kelvin ,&! zero point for Celsius
pi ,&! pi
mpercm ,&! meters per m
cmperm ,&! m per meter
days_per_second ,&! days per second
Expand Down Expand Up @@ -117,6 +118,7 @@ subroutine ocn_constants_init(configPool, packagePool)!{{{
!-----------------------------------------------------------------------

T0_Kelvin = 273.16_RKIND ! zero point for Celsius
pi = 3.14159265358979323846_RKIND ! pi
rho_air = 1.2_RKIND ! ambient air density (kg/m^3)
rho_sw = config_density0 ! density of salt water (kg/m^3)
rho_fw = 1.0e3_RKIND ! avg. water density (kg/m^3)
Expand Down Expand Up @@ -159,6 +161,7 @@ subroutine ocn_constants_init(configPool, packagePool)!{{{

#ifdef MPAS_ESM_SHR_CONST
T0_Kelvin = SHR_CONST_TKFRZ ! zero point for Celsius
pi = SHR_CONST_PI ! zero point for Celsius
cp_sw = SHR_CONST_CPSW ! erg/g/K
cp_air = SHR_CONST_CPDAIR ! J/kg/K
rho_air = SHR_CONST_RHODAIR ! kg/m^3
Expand Down
4 changes: 2 additions & 2 deletions components/mpas-ocean/src/shared/mpas_ocn_gm.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -484,7 +484,7 @@ subroutine ocn_GM_compute_Bolus_velocity(statePool, &

! Compute the speed of the first baroclinic mode from the Brunt-Vaisala frequency.
cGMphaseSpeed(iEdge) = max(c_min, &
sum_hN/(config_GM_spatially_variable_baroclinic_mode*3.141592_RKIND))
sum_hN/(config_GM_spatially_variable_baroclinic_mode*pi))

end do
!$omp end do
Expand Down Expand Up @@ -615,7 +615,7 @@ subroutine ocn_GM_compute_Bolus_velocity(statePool, &
! for the first-mode (m=1) baroclinic gravity wave speed in m/s
! c_m = 1/pi integral(N dz )

c_mode1 = max(1.0E-5_RKIND, sum_hN/3.141592_RKIND)
c_mode1 = max(1.0E-5_RKIND, sum_hN/pi)

! See Hallberg (2013) https://doi.org/10.1016/j.ocemod.2013.08.007
! just after eqn 4 the defines the deformation
Expand Down
4 changes: 2 additions & 2 deletions components/mpas-ocean/src/shared/mpas_ocn_manufactured_solution.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -217,8 +217,8 @@ subroutine ocn_manufactured_solution_init(domain, err)!{{{

if (.not. config_use_manufactured_solution) return

kx = 2.0_RKIND*pii / config_manufactured_solution_wavelength_x
ky = 2.0_RKIND*pii / config_manufactured_solution_wavelength_y
kx = 2.0_RKIND*pi / config_manufactured_solution_wavelength_x
ky = 2.0_RKIND*pi / config_manufactured_solution_wavelength_y

eta0 = config_manufactured_solution_amplitude

Expand Down
3 changes: 2 additions & 1 deletion components/mpas-ocean/src/shared/mpas_ocn_mesh.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -33,6 +33,7 @@ module ocn_mesh
use mpas_log

use ocn_config
use ocn_constants

implicit none
private
Expand Down Expand Up @@ -1792,7 +1793,7 @@ subroutine ocn_meshScaling() !{{{
! neighboring cells are circles for this calculation.
cellWidth = 2.0_RKIND* &
sqrt((areaCell(cell1) + areaCell(cell2))/ &
2.0_RKIND/pii)
2.0_RKIND/pi)
meshScalingDel2(iEdge) = cellWidth/ &
config_hmix_ref_cell_width
meshScalingDel4(iEdge) = (cellWidth/ &
Expand Down
4 changes: 2 additions & 2 deletions components/mpas-ocean/src/shared/mpas_ocn_tendency.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1453,15 +1453,15 @@ subroutine ocn_tend_freq_filtered_thickness(tendPool, statePool, &
do k = 1, maxLevelCell(iCell)

tend_lowFreqDivergence(k,iCell) = &
-2.0*pii/thickness_filter_timescale_sec &
-2.0*pi/thickness_filter_timescale_sec &
*(lowFreqDivergence(k,iCell) - div_hu(k) &
+ div_hu_btr * layerThickness(k,iCell)/totalThickness)

tend_highFreqThickness(k,iCell) = - div_hu(k) + &
div_hu_btr*layerThickness(k,iCell)/totalThickness + &
lowFreqDivergence(k,iCell) + &
use_highFreqThick_restore* &
(-2.0*pii/highFreqThick_restore_time_sec* &
(-2.0*pi/highFreqThick_restore_time_sec* &
highFreqThickness(k,iCell) )

end do ! vert (k) loop
Expand Down
4 changes: 2 additions & 2 deletions components/mpas-ocean/src/shared/mpas_ocn_tidal_forcing.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -243,8 +243,8 @@ subroutine ocn_tidal_forcing_build_array(domain, meshPool, forcingPool, statePoo
! compute the tidalHeight
if (trim(config_tidal_forcing_model) == 'monochromatic') then
tidalHeight = config_tidal_forcing_monochromatic_amp * &
SIN(2.0_RKIND*pii/config_tidal_forcing_monochromatic_period * daysSinceStartOfSim - &
pii*config_tidal_forcing_monochromatic_phaseLag/180.0_RKIND) - &
SIN(2.0_RKIND*pi/config_tidal_forcing_monochromatic_period * daysSinceStartOfSim - &
pi*config_tidal_forcing_monochromatic_phaseLag/180.0_RKIND) - &
config_tidal_forcing_monochromatic_baseline
elseif (trim(config_tidal_forcing_model) == 'linear') then
tidalHeight = max(config_tidal_forcing_linear_min, &
Expand Down
10 changes: 5 additions & 5 deletions components/mpas-ocean/src/shared/mpas_ocn_tracer_advection_shared.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -24,10 +24,12 @@ module ocn_tracer_advection_shared
use mpas_derived_types
use mpas_hash
use mpas_sort
use mpas_constants
use mpas_geometry_utils
use mpas_log

use ocn_config
use ocn_constants
use ocn_mesh

implicit none
Expand Down Expand Up @@ -370,7 +372,6 @@ subroutine computeDerivTwo(derivTwo, err)!{{{
xec, yec, zec, &! arc bisection coords
thetae_tmp, &! angle
xv1, xv2, yv1, yv2, zv1, zv2, &! vertex cart coords
pii, &! pi math constant
length_scale, &! length scale
cos2t, costsint, sin2t ! trig function temps

Expand Down Expand Up @@ -413,7 +414,6 @@ subroutine computeDerivTwo(derivTwo, err)!{{{

! Initialize derivTwo and pi

pii = 2.*asin(1.0)
derivTwo(:,:,:) = 0.0_RKIND

do iCell = 1, nCellsAll
Expand Down Expand Up @@ -461,9 +461,9 @@ subroutine computeDerivTwo(derivTwo, err)!{{{
end do

if ( zc(1) == 1.0_RKIND) then
theta_abs(iCell) = pii/2.0_RKIND
theta_abs(iCell) = pi/2.0_RKIND
else
theta_abs(iCell) = pii/2.0_RKIND &
theta_abs(iCell) = pi/2.0_RKIND &
- mpas_sphere_angle( xc(1), yc(1), zc(1), &
xc(2), yc(2), zc(2), &
0.0_RKIND, 0.0_RKIND, 1.0_RKIND)
Expand Down Expand Up @@ -509,7 +509,7 @@ subroutine computeDerivTwo(derivTwo, err)!{{{
angle_2d(i) = angleEdge(edgesOnCell(i,iCell))
iEdge = edgesOnCell(i,iCell)
if ( iCell .ne. cellsOnEdge(1,iEdge)) &
angle_2d(i) = angle_2d(i) - pii
angle_2d(i) = angle_2d(i) - pi

xp(i) = dcEdge(edgesOnCell(i,iCell)) * cos(angle_2d(i))
yp(i) = dcEdge(edgesOnCell(i,iCell)) * sin(angle_2d(i))
Expand Down
6 changes: 3 additions & 3 deletions components/mpas-ocean/src/shared/mpas_ocn_vel_forcing_topographic_wave_drag.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -199,7 +199,7 @@ subroutine ocn_vel_forcing_topographic_wave_drag_init(err)!{{{
call mpas_log_write(" Topographic wave drag scheme is: Jayne and St. Laurent")
call mpas_log_write("")
tensorScheme = .false.
kappa = pii/10000.0_RKIND
kappa = pi/10000.0_RKIND

do iEdge = 1, nEdgesAll
kmax = maxLevelEdgeTop(iEdge)
Expand Down Expand Up @@ -241,7 +241,7 @@ subroutine ocn_vel_forcing_topographic_wave_drag_init(err)!{{{

topographicWaveDrag(iEdge) = topographicWaveDragCoeff * gam &
* bed_slope_edges(iEdge)**2 * Nbar * Nb &
/ (8.0_RKIND * omegaM2 * pii**2)
/ (8.0_RKIND * omegaM2 * pi**2)
endif
enddo
else if (config_topographic_wave_drag_scheme.EQ."LGF") then
Expand All @@ -256,7 +256,7 @@ subroutine ocn_vel_forcing_topographic_wave_drag_init(err)!{{{
topographicWaveDrag(iEdge) = topographicWaveDragCoeff &
* (sqrt((topo_buoyancy_N1B(iEdge)**2 - omegaM2**2) &
* (topo_buoyancy_N1V(iEdge)**2 - omegaM2**2))) &
/ (4.0_RKIND*pii*omegaM2)
/ (4.0_RKIND*pi*omegaM2)
normalCoeffTWD(iEdge) = (lonGradEdge(iEdge)**2) &
* (cos(angleEdge(iEdge))**2) + (latGradEdge(iEdge)**2)*(sin(angleEdge(iEdge))**2) &
- 2.0_RKIND*latGradEdge(iEdge)*lonGradEdge(iEdge)*sin(angleEdge(iEdge))*(cos(angleEdge(iEdge)))
Expand Down
2 changes: 1 addition & 1 deletion components/mpas-ocean/src/shared/mpas_ocn_vel_hmix_leith.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -159,7 +159,7 @@ subroutine ocn_vel_hmix_leith_tend(div, relVort, tend, err)!{{{
dcEdgeInv = 1.0_RKIND / dcEdge(iEdge)
dvEdgeInv = 1.0_RKIND / dvEdge(iEdge)

visc2tmp = (leithParam*dxLeith*meshScalingDel2(iEdge)/pii)**3
visc2tmp = (leithParam*dxLeith*meshScalingDel2(iEdge)/pi)**3

do k = minLevelEdgeBot(iEdge), maxLevelEdgeTop(iEdge)

Expand Down
6 changes: 3 additions & 3 deletions components/mpas-ocean/src/shared/mpas_ocn_vel_self_attraction_loading.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -863,8 +863,8 @@ subroutine ocn_vel_self_attraction_loading_init(domain,err)!{{{
! Pre-compute sin and cos of latCell (co-latitude) values
allocate(sinLatCell(nCellsAll), cosLatCell(nCellsAll))
do iCell = 1,nCellsAll
sinLatCell(iCell) = sin(0.5_RKIND*pii-latCell(iCell))
cosLatCell(iCell) = cos(0.5_RKIND*pii-latCell(iCell))
sinLatCell(iCell) = sin(0.5_RKIND*pi-latCell(iCell))
cosLatCell(iCell) = cos(0.5_RKIND*pi-latCell(iCell))
enddo

! Calculate blocking indices
Expand Down Expand Up @@ -2526,7 +2526,7 @@ subroutine associatedLegendrePolynomials(n, m, startIdx, endIdx, l, pmnm2, pmnm1
if (n == m) then

do iCell = startIdx,endIdx
pmnm2(iCell) = sqrt(1.0_RKIND/(4.0_RKIND*pii))*sinLatCell(iCell)**m
pmnm2(iCell) = sqrt(1.0_RKIND/(4.0_RKIND*pi))*sinLatCell(iCell)**m
do i = 1,m
pmnm2(iCell) = pmnm2(iCell)*sqrt(real(2*i+1,RKIND)/real(2*i,RKIND))
enddo
Expand Down
Loading