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

Jump to bottom

User/z1l/remove edge arrays #41

Merged
merged 3 commits into from
Aug 26, 2014
Merged

User/z1l/remove edge arrays #41

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
f3725a7
Implement Openmp in MOM_PressureForce. Block data structure is used i…
Zhi-Liang Jul 28, 2014
d401e43
Merge remote-tracking branch 'origin' into user/z1l/openmp
Zhi-Liang Jul 28, 2014
37a2eff
remove type edgeValueArrays and edgeSlopeArrays for openmp
Zhi-Liang Jul 31, 2014
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
23 changes: 4 additions & 19 deletions src/ALE/MOM_ALE.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -23,21 +23,15 @@ module MOM_ALE
use MOM_EOS, only : calculate_density
use MOM_string_functions, only : uppercase, extractWord
use MOM_verticalGrid, only : verticalGrid_type

use regrid_edge_values, only : edgeValueArrays
use regrid_edge_values, only : edge_values_implicit_h4
use regrid_edge_values, only : triDiagEdgeWorkAllocate, triDiagEdgeWorkDeallocate
use regrid_edge_slopes, only : edgeSlopeArrays
use regrid_edge_slopes, only : triDiagSlopeWorkAllocate, triDiagSlopeWorkDeallocate

use PLM_functions, only : PLM_reconstruction, PLM_boundary_extrapolation
use PPM_functions, only : PPM_reconstruction, PPM_boundary_extrapolation

use P1M_functions, only : P1M_interpolation, P1M_boundary_extrapolation
use P3M_functions, only : P3M_interpolation, P3M_boundary_extrapolation
use MOM_regridding, only : initialize_regridding, regridding_main , end_regridding
use MOM_regridding, only : uniformResolution
use MOM_regridding, only : check_grid_integrity, setCoordinateResolution
use MOM_regridding, only : setcoordinateinterfaces
use MOM_regridding, only : regriddingCoordinateModeDoc, DEFAULT_COORDINATE_MODE
use MOM_regridding, only : regriddingInterpSchemeDoc, regriddingDefaultInterpScheme
use MOM_regridding, only : setRegriddingBoundaryExtrapolation
Expand Down Expand Up @@ -82,8 +76,6 @@ module MOM_ALE

type(regridding_CS) :: regridCS ! Regridding parameters and work arrays
type(remapping_CS) :: remapCS ! Remapping parameters and work arrays
type(edgeValueArrays) :: edgeValueWrk ! Work space for edge values
type(edgeSlopeArrays) :: edgeSlopeWrk ! Work space for edge slopes

! Used only for queries, not directly by this module
integer :: nk
Expand Down Expand Up @@ -424,7 +416,7 @@ subroutine pressure_gradient_ppm( CS, S_t, S_b, T_t, T_b, G, tv, h )
! Reconstruct salinity profile
ppoly_parab_E = 0.0
ppoly_parab_coefficients = 0.0
call edge_values_implicit_h4( G%ke, hTmp, tmp, CS%edgeValueWrk, ppoly_parab_E )
call edge_values_implicit_h4( G%ke, hTmp, tmp, ppoly_parab_E )
call PPM_reconstruction( G%ke, hTmp, tmp, ppoly_parab_E, ppoly_parab_coefficients )
if (CS%boundary_extrapolation_for_pressure) call &
PPM_boundary_extrapolation( G%ke, hTmp, tmp, ppoly_parab_E, ppoly_parab_coefficients )
Expand All @@ -438,7 +430,7 @@ subroutine pressure_gradient_ppm( CS, S_t, S_b, T_t, T_b, G, tv, h )
ppoly_parab_E = 0.0
ppoly_parab_coefficients = 0.0
tmp(:) = tv%T(i,j,:)
call edge_values_implicit_h4( G%ke, hTmp, tmp, CS%edgeValueWrk, ppoly_parab_E )
call edge_values_implicit_h4( G%ke, hTmp, tmp, ppoly_parab_E )
call PPM_reconstruction( G%ke, hTmp, tmp, ppoly_parab_E, ppoly_parab_coefficients )
if (CS%boundary_extrapolation_for_pressure) call &
PPM_boundary_extrapolation( G%ke, hTmp, tmp, ppoly_parab_E, ppoly_parab_coefficients )
Expand Down Expand Up @@ -476,10 +468,6 @@ subroutine ALE_memory_allocation( G, CS )

nz = G%ke

! Allocate memory for the tridiagonal system
call triDiagEdgeWorkAllocate( nz, CS%edgeValueWrk )
call triDiagSlopeWorkAllocate( nz, CS%edgeSlopeWrk )

! Work space
ALLOC_(CS%dzRegrid(G%isd:G%ied,G%jsd:G%jed,nz+1)); CS%dzRegrid(:,:,:) = 0.

Expand All @@ -496,10 +484,6 @@ subroutine ALE_memory_deallocation( CS )

type(ALE_CS), intent(inout) :: CS

! Reclaim memory for the tridiagonal system
call triDiagEdgeWorkDeallocate( CS%edgeValueWrk )
call triDiagSlopeWorkDeallocate( CS%edgeSlopeWrk )

! Work space
DEALLOC_(CS%dzRegrid)

Expand Down Expand Up @@ -623,6 +607,7 @@ subroutine ALE_initRegridding( G, param_file, mod, regridCS, dz )
endif
endif
call setCoordinateResolution( dz, regridCS )
call setCoordinateInterfaces( G, regridCS )

call get_param(param_file, mod, "MIN_THICKNESS", tmpReal, &
"When regridding, this is the minimum layer\n"//&
Expand Down
71 changes: 35 additions & 36 deletions src/ALE/MOM_regridding.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -19,13 +19,9 @@ module MOM_regridding
use MOM_EOS, only : calculate_density
use MOM_string_functions,only : uppercase

use regrid_edge_values, only : edgeValueArrays
use regrid_edge_values, only : edge_values_explicit_h2, edge_values_explicit_h4
use regrid_edge_values, only : edge_values_implicit_h4, edge_values_implicit_h6
use regrid_edge_values, only : triDiagEdgeWorkAllocate, triDiagEdgeWorkDeallocate
use regrid_edge_slopes, only : edgeSlopeArrays
use regrid_edge_slopes, only : edge_slopes_implicit_h3, edge_slopes_implicit_h5
use regrid_edge_slopes, only : triDiagSlopeWorkAllocate, triDiagSlopeWorkDeallocate

use PLM_functions, only : PLM_reconstruction, PLM_boundary_extrapolation
use PPM_functions, only : PPM_reconstruction, PPM_boundary_extrapolation
Expand Down Expand Up @@ -56,6 +52,7 @@ module MOM_regridding
! coorindate. It has the units of the target coordiante, e.g.
! meters for z*, non-dimensional for sigma, etc.
real, dimension(:), allocatable :: coordinateResolution
! This array is set by function setcoordinateInterfaces
! This array is nominal coordinate of interfaces and is the
! running sum of coordinateResolution. i.e.
! coordinateInterfaces(k) = coordinateResolution(k+1)-coordinateResolution(k)
Expand All @@ -82,8 +79,6 @@ module MOM_regridding
! Minimum thickness allowed when building the new grid through regridding
real :: min_thickness

type(edgeValueArrays) :: edgeValueWrk ! Work space for edge values
type(edgeSlopeArrays) :: edgeSlopeWrk ! Work space for edge slopes
end type

! -----------------------------------------------------------------------------
Expand All @@ -97,6 +92,7 @@ module MOM_regridding
public setRegriddingMinimumThickness
public uniformResolution
public setCoordinateResolution
public setCoordinateInterfaces
public getCoordinateResolution
public getCoordinateInterfaces
public getCoordinateUnits
Expand Down Expand Up @@ -250,8 +246,8 @@ subroutine regridding_main( remapCS, CS, G, h, tv, dzInterface )
!------------------------------------------------------------------------------

! Arguments
type(remapping_CS), intent(inout) :: remapCS ! Remapping parameters and options
type(regridding_CS), intent(inout) :: CS ! Regridding parameters and options
type(remapping_CS), intent(in) :: remapCS ! Remapping parameters and options
type(regridding_CS), intent(in) :: CS ! Regridding parameters and options
type(ocean_grid_type), intent(in) :: G ! Ocean grid informations
real, dimension(NIMEM_,NJMEM_, NKMEM_), intent(inout) :: h ! Current 3D grid obtained after the last time step
type(thermo_var_ptrs), intent(inout) :: tv ! Thermodynamical variables (T, S, ...)
Expand Down Expand Up @@ -571,8 +567,8 @@ subroutine buildGridRho( G, h, tv, dzInterface, remapCS, CS )
real, dimension(NIMEM_,NJMEM_, NKMEM_), intent(in) :: h
type(thermo_var_ptrs), intent(in) :: tv
real, dimension(NIMEM_,NJMEM_, NK_INTERFACE_), intent(inout) :: dzInterface
type(remapping_CS), intent(inout) :: remapCS
type(regridding_CS), intent(inout) :: CS
type(remapping_CS), intent(in) :: remapCS
type(regridding_CS), intent(in) :: CS

! Local variables
integer :: i, j, k, m
Expand Down Expand Up @@ -600,13 +596,6 @@ subroutine buildGridRho( G, h, tv, dzInterface, remapCS, CS )
threshold = CS%min_thickness
p_column(:) = 0.

! Prescribe target values
CS%coordinateInterfaces(1) = G%Rlay(1)+0.5*(G%Rlay(1)-G%Rlay(2))
CS%coordinateInterfaces(nz+1) = G%Rlay(nz)+0.5*(G%Rlay(nz)-G%Rlay(nz-1))
do k = 2,nz
CS%coordinateInterfaces(k) = CS%coordinateInterfaces(k-1) + CS%coordinateResolution(k)
end do

! Build grid based on target interface densities
do i = G%isc-1,G%iec+1
do j = G%jsc-1,G%jec+1
Expand Down Expand Up @@ -910,7 +899,7 @@ subroutine regridding_iteration( densities, target_values, CS, &
integer, intent(in) :: n1 ! Number of cells on target grid
real, dimension(:), intent(out) :: h1 ! cell widths on target grid
real, dimension(:), intent(out) :: x1 ! interface positions on target grid
type(regridding_CS), intent(inout) :: CS ! Parameters used for regridding
type(regridding_CS), intent(in) :: CS ! Parameters used for regridding

! Local variables
integer :: degree, k
Expand Down Expand Up @@ -946,7 +935,7 @@ subroutine regridding_iteration( densities, target_values, CS, &
case ( INTERPOLATION_P1M_IH4 )
degree = DEGREE_1
if ( n0 .ge. 4 ) then
call edge_values_implicit_h4( n0, h0, densities, CS%edgeValueWrk, ppoly0_E )
call edge_values_implicit_h4( n0, h0, densities, ppoly0_E )
else
call edge_values_explicit_h2( n0, h0, densities, ppoly0_E )
end if
Expand Down Expand Up @@ -983,7 +972,7 @@ subroutine regridding_iteration( densities, target_values, CS, &

if ( n0 .ge. 4 ) then
degree = DEGREE_2
call edge_values_implicit_h4( n0, h0, densities, CS%edgeValueWrk, ppoly0_E )
call edge_values_implicit_h4( n0, h0, densities, ppoly0_E )
call PPM_reconstruction( n0, h0, densities, ppoly0_E, ppoly0_coefficients )
if ( CS%boundary_extrapolation) then
call PPM_boundary_extrapolation( n0, h0, densities, ppoly0_E, ppoly0_coefficients )
Expand All @@ -1001,8 +990,8 @@ subroutine regridding_iteration( densities, target_values, CS, &

if ( n0 .ge. 4 ) then
degree = DEGREE_3
call edge_values_implicit_h4( n0, h0, densities, CS%edgeValueWrk, ppoly0_E )
call edge_slopes_implicit_h3( n0, h0, densities, CS%edgeSlopeWrk, ppoly0_S )
call edge_values_implicit_h4( n0, h0, densities, ppoly0_E )
call edge_slopes_implicit_h3( n0, h0, densities, ppoly0_S )
call P3M_interpolation( n0, h0, densities, ppoly0_E, ppoly0_S, ppoly0_coefficients )
if ( CS%boundary_extrapolation) then
call P3M_boundary_extrapolation( n0, h0, densities, ppoly0_E, ppoly0_S, ppoly0_coefficients )
Expand All @@ -1019,8 +1008,8 @@ subroutine regridding_iteration( densities, target_values, CS, &
case ( INTERPOLATION_P3M_IH6IH5 )
if ( n0 .ge. 6 ) then
degree = DEGREE_3
call edge_values_implicit_h6( n0, h0, densities, CS%edgeValueWrk, ppoly0_E )
call edge_slopes_implicit_h5( n0, h0, densities, CS%edgeSlopeWrk, ppoly0_S )
call edge_values_implicit_h6( n0, h0, densities, ppoly0_E )
call edge_slopes_implicit_h5( n0, h0, densities, ppoly0_S )
call P3M_interpolation( n0, h0, densities, ppoly0_E, ppoly0_S, ppoly0_coefficients )
if ( CS%boundary_extrapolation) then
call P3M_boundary_extrapolation( n0, h0, densities, ppoly0_E, ppoly0_S, ppoly0_coefficients )
Expand All @@ -1038,8 +1027,8 @@ subroutine regridding_iteration( densities, target_values, CS, &

if ( n0 .ge. 4 ) then
degree = DEGREE_4
call edge_values_implicit_h4( n0, h0, densities, CS%edgeValueWrk, ppoly0_E )
call edge_slopes_implicit_h3( n0, h0, densities, CS%edgeSlopeWrk, ppoly0_S )
call edge_values_implicit_h4( n0, h0, densities, ppoly0_E )
call edge_slopes_implicit_h3( n0, h0, densities, ppoly0_S )
call PQM_reconstruction( n0, h0, densities, ppoly0_E, ppoly0_S, ppoly0_coefficients )
if ( CS%boundary_extrapolation) then
call PQM_boundary_extrapolation_v1( n0, h0, densities, ppoly0_E, ppoly0_S, ppoly0_coefficients )
Expand All @@ -1056,8 +1045,8 @@ subroutine regridding_iteration( densities, target_values, CS, &
case ( INTERPOLATION_PQM_IH6IH5 )
if ( n0 .ge. 6 ) then
degree = DEGREE_4
call edge_values_implicit_h6( n0, h0, densities, CS%edgeValueWrk, ppoly0_E )
call edge_slopes_implicit_h5( n0, h0, densities, CS%edgeSlopeWrk, ppoly0_S )
call edge_values_implicit_h6( n0, h0, densities, ppoly0_E )
call edge_slopes_implicit_h5( n0, h0, densities, ppoly0_S )
call PQM_reconstruction( n0, h0, densities, ppoly0_E, ppoly0_S, ppoly0_coefficients )
if ( CS%boundary_extrapolation) then
call PQM_boundary_extrapolation_v1( n0, h0, densities, ppoly0_E, ppoly0_S, ppoly0_coefficients )
Expand Down Expand Up @@ -1527,6 +1516,23 @@ subroutine setCoordinateResolution( dz, CS )

end subroutine setCoordinateResolution

!-----------------------------------------------------------------------------
! Set the running sum of coordinateResolution
!-----------------------------------------------------------------------------
subroutine setCoordinateInterfaces( G, CS )
type(ocean_grid_type), intent(in) :: G ! Ocean grid informations
type(regridding_CS), intent(inout) :: CS
integer :: k, nz

nz = CS%nk
CS%coordinateInterfaces(1) = G%Rlay(1)+0.5*(G%Rlay(1)-G%Rlay(2))
CS%coordinateInterfaces(nz+1) = G%Rlay(nz)+0.5*(G%Rlay(nz)-G%Rlay(nz-1))
do k = 2,nz
CS%coordinateInterfaces(k) = CS%coordinateInterfaces(k-1) + CS%coordinateResolution(k)
end do

end subroutine setCoordinateInterfaces

!------------------------------------------------------------------------------
! Query the fixed resolution data
!------------------------------------------------------------------------------
Expand Down Expand Up @@ -1647,10 +1653,6 @@ subroutine allocate_regridding( CS )

! Local variables

! Allocate memory for the tridiagonal system
call triDiagEdgeWorkAllocate( CS%nk, CS%edgeValueWrk )
call triDiagSlopeWorkAllocate( CS%nk, CS%edgeSlopeWrk )

! Target values
allocate( CS%coordinateInterfaces(CS%nk+1) )

Expand All @@ -1670,12 +1672,9 @@ subroutine regridding_memory_deallocation( CS )

type(regridding_CS), intent(inout) :: CS

! Reclaim memory for the tridiagonal system
call triDiagEdgeWorkDeallocate( CS%edgeValueWrk )
call triDiagSlopeWorkDeallocate( CS%edgeSlopeWrk )

! Target values
deallocate( CS%coordinateInterfaces )
deallocate( CS%coordinateResolution )

end subroutine regridding_memory_deallocation

Expand Down
44 changes: 7 additions & 37 deletions src/ALE/MOM_remapping.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -13,13 +13,9 @@ module MOM_remapping
use MOM_string_functions, only : uppercase
use MOM_variables, only : ocean_grid_type, thermo_var_ptrs
use polynomial_functions, only : evaluation_polynomial, integration_polynomial
use regrid_edge_values, only : edgeValueArrays
use regrid_edge_values, only : edge_values_explicit_h4, edge_values_implicit_h4
use regrid_edge_values, only : edge_values_implicit_h4, edge_values_implicit_h6
use regrid_edge_values, only : triDiagEdgeWorkAllocate, triDiagEdgeWorkDeallocate
use regrid_edge_slopes, only : edgeSlopeArrays
use regrid_edge_slopes, only : edge_slopes_implicit_h3, edge_slopes_implicit_h5
use regrid_edge_slopes, only : triDiagSlopeWorkAllocate, triDiagSlopeWorkDeallocate
use PCM_functions, only : PCM_reconstruction
use PLM_functions, only : PLM_reconstruction, PLM_boundary_extrapolation
use PPM_functions, only : PPM_reconstruction, PPM_boundary_extrapolation
Expand All @@ -34,9 +30,6 @@ module MOM_remapping
! -----------------------------------------------------------------------------
type, public :: remapping_CS
private
! Work arrays
type(edgeValueArrays) :: edgeValueWrk ! Work space for edge values
type(edgeSlopeArrays) :: edgeSlopeWrk ! Work space for edge slopes
! Parameters
integer :: nk = 0 ! Number of layers/levels in vertical
integer :: remapping_scheme = -911 ! Determines which reconstruction to use
Expand Down Expand Up @@ -107,7 +100,7 @@ subroutine remapping_main( CS, G, h, dxInterface, tv, u, v )
!------------------------------------------------------------------------------

! Arguments
type(remapping_CS), intent(inout) :: CS
type(remapping_CS), intent(in) :: CS
type(ocean_grid_type), intent(in) :: G
real, dimension(NIMEM_,NJMEM_,NKMEM_), intent(in) :: h
real, dimension(NIMEM_,NJMEM_,NK_INTERFACE_), intent(in) :: dxInterface
Expand Down Expand Up @@ -373,7 +366,7 @@ subroutine remapping_core( CS, n0, h0, u0, n1, dx, u1 )
!------------------------------------------------------------------------------

! Arguments
type(remapping_CS), intent(inout) :: CS
type(remapping_CS), intent(in) :: CS
integer, intent(in) :: n0 ! Number of cells on source grid
real, dimension(:), intent(in) :: h0 ! cell widths on source grid
real, dimension(:), intent(in) :: u0 ! cell averages on source grid
Expand Down Expand Up @@ -460,23 +453,23 @@ subroutine remapping_core( CS, n0, h0, u0, n1, dx, u1 )
end if
iMethod = INTEGRATION_PPM
case ( REMAPPING_PPM_IH4 )
call edge_values_implicit_h4( n0, h0, u0, CS%edgeValueWrk, ppoly_r_E )
call edge_values_implicit_h4( n0, h0, u0, ppoly_r_E )
call PPM_reconstruction( n0, h0, u0, ppoly_r_E, ppoly_r_coefficients )
if ( CS%boundary_extrapolation) then
call PPM_boundary_extrapolation( n0, h0, u0, ppoly_r_E, ppoly_r_coefficients )
end if
iMethod = INTEGRATION_PPM
case ( REMAPPING_PQM_IH4IH3 )
call edge_values_implicit_h4( n0, h0, u0, CS%edgeValueWrk, ppoly_r_E )
call edge_slopes_implicit_h3( n0, h0, u0, CS%edgeSlopeWrk, ppoly_r_S )
call edge_values_implicit_h4( n0, h0, u0, ppoly_r_E )
call edge_slopes_implicit_h3( n0, h0, u0, ppoly_r_S )
call PQM_reconstruction( n0, h0, u0, ppoly_r_E, ppoly_r_S, ppoly_r_coefficients )
if ( CS%boundary_extrapolation) then
call PQM_boundary_extrapolation_v1( n0, h0, u0, ppoly_r_E, ppoly_r_S, ppoly_r_coefficients )
end if
iMethod = INTEGRATION_PQM
case ( REMAPPING_PQM_IH6IH5 )
call edge_values_implicit_h6( n0, h0, u0, CS%edgeValueWrk, ppoly_r_E )
call edge_slopes_implicit_h5( n0, h0, u0, CS%edgeSlopeWrk, ppoly_r_S )
call edge_values_implicit_h6( n0, h0, u0, ppoly_r_E )
call edge_slopes_implicit_h5( n0, h0, u0, ppoly_r_S )
call PQM_reconstruction( n0, h0, u0, ppoly_r_E, ppoly_r_S, ppoly_r_coefficients )
if ( CS%boundary_extrapolation) then
call PQM_boundary_extrapolation_v1( n0, h0, u0, ppoly_r_E, ppoly_r_S, ppoly_r_coefficients )
Expand Down Expand Up @@ -1115,13 +1108,6 @@ subroutine setReconstructionType(string,CS)
"Unrecognized choice for REMAPPING_SCHEME ("//trim(string)//").")
end select

if (CS%degree > 0 .and. degree/=CS%degree) then
! If the degree has changed then deallocate to force a re-allocation
call end_remapping(CS)
endif
if ( CS%degree == 0 ) then
call allocate_remapping( CS )
endif
CS%degree = degree


Expand All @@ -1145,29 +1131,13 @@ subroutine remapDisableBoundaryExtrapolation(CS)
CS%boundary_extrapolation = .false.
end subroutine remapDisableBoundaryExtrapolation

!------------------------------------------------------------------------------
! Memory allocation for remapping
!------------------------------------------------------------------------------
subroutine allocate_remapping( CS )
! Arguments
type(remapping_CS), intent(inout) :: CS

call triDiagEdgeWorkAllocate( CS%nk, CS%edgeValueWrk )
call triDiagSlopeWorkAllocate( CS%nk, CS%edgeSlopeWrk )

end subroutine allocate_remapping


!------------------------------------------------------------------------------
! Memory deallocation for remapping
!------------------------------------------------------------------------------
subroutine end_remapping(CS)
! Arguments
type(remapping_CS), intent(inout) :: CS

! Deallocate memory for grid
call triDiagEdgeWorkDeallocate( CS%edgeValueWrk )
call triDiagSlopeWorkDeallocate( CS%edgeSlopeWrk )
CS%degree = 0

end subroutine end_remapping
Expand Down
Loading