#include "cppdefs.h" MODULE mod_param ! !git $Id$ !svn $Id: mod_param.F 1202 2023-10-24 15:36:07Z arango $ !================================================== Hernan G. Arango === ! Copyright (c) 2002-2023 The ROMS/TOMS Group ! ! Licensed under a MIT/X style license ! ! See License_ROMS.md ! !======================================================================= ! ! ! Grid parameters: ! ! ! ! Im Number of global grid points in the XI-direction ! ! for each nested grid. ! ! Jm Number of global grid points in the ETA-direction ! ! for each nested grid. ! ! Lm Number of interior grid points in the XI-direction ! ! for each nested grid. ! ! Mm Number of internal grid points in the ETA-direction. ! ! for each nested grid. ! ! N Number of vertical levels for each nested grid. ! ! Ngrids Number of nested and/or connected grids to solve. ! ! NtileI Number of XI-direction tiles or domain partitions for ! ! each nested grid. Values used to compute tile ranges. ! ! NtileJ Number of ETA-direction tiles or domain partitions for ! ! each nested grid. Values used to compute tile ranges. ! ! NtileX Number of XI-direction tiles or domain partitions for ! ! each nested grid. Values used in parallel loops. ! ! NtileE Number of ETA-direction tiles or domain partitions for ! ! each nested grid. Values used in parallel loops. ! ! HaloBry Buffers halo size for exchanging boundary arrays. ! ! HaloSizeI Maximum halo size, in grid points, in XI-direction. ! ! HaloSizeJ Maximum halo size, in grid points, in ETA-direction. ! ! TileSide Maximun tile side length in XI- or ETA-directions. ! ! TileSize Maximum tile size. ! ! ! ! Configuration parameters: ! ! ! ! Nbico Number of balanced SSH elliptic equation iterations. ! ! Nfloats Number of floats trajectories. ! ! Nstation Number of output stations. ! ! MTC Maximum number of tidal components. ! ! ! ! State variables parameters: ! ! ! #ifdef PROPAGATOR ! Mstate Size of FULL state vector (water points only). ! ! Nstate Size of NODE partition state vector: Nstate=Mstate in ! ! serial applications. ! ! Nstr State vector node partition starting index. ! ! Nend State vector node partition ending index. ! ! Nsize Size of the eigenvalue problem: Nend-Nstr+1. ! # ifdef SO_SEMI ! Nsemi Number of time record saved in seminorm strochastic ! ! optimals adjoint state vector. ! # endif #endif ! NSA Number of state array for error covariance. ! ! NSV Number of model state variables. ! ! ! ! Tracer parameters: ! ! ! ! NAT Number of active tracer type variables (usually, ! ! NAT=2 for potential temperature and salinity). ! ! NBT Number of biological tracer type variables. ! ! NST Number of sediment tracer type variables (NCS+NNS). ! ! NPT Number of extra passive tracer type variables to ! ! advect and diffuse only (dyes, etc). ! ! MT Maximum number of tracer type variables. ! ! NT Total number of tracer type variables. ! ! NTCLM Total number of climatology tracer type variables to ! ! process. ! ! ! ! Nbed Number of sediment bed layers. ! ! NCS Number of cohesive (mud) sediment tracers. ! ! NNS Number of non-cohesive (sand) sediment tracers. ! ! ! ! Diagnostic fields parameters: ! ! ! # if defined ECOSIM && defined BIO_OPTIC ! NDbands Number of spectral irradiance bands. ! # endif ! NDbio2d Number of diagnostic 2D biology fields. ! ! NDbio3d Number of diagnostic 3D biology fields. ! ! NDbio4d Number of diagnostic 4D bio-optical fields. ! ! NDT Number of diagnostic tracer fields. ! ! NDM2d Number of diagnostic 2D momentum fields. ! ! NDM3d Number of diagnostic 3D momentum fields. ! ! NDrhs Number of diagnostic 3D right-hand-side fields. ! ! ! ! Estimated Dynamic and automatic memory parameters: ! ! ! ! BmemMax Maximum automatic memory for distributed-memory buffers ! ! (bytes). ! ! Dmem Dynamic memory requirements (array elements) ! ! ! !======================================================================= ! USE mod_kinds ! implicit none ! PUBLIC :: allocate_param PUBLIC :: deallocate_param PUBLIC :: initialize_param ! !----------------------------------------------------------------------- ! Grid nesting parameters. !----------------------------------------------------------------------- ! ! Number of nested and/or connected grids to solve. ! integer :: Ngrids ! ! Number of grid nesting layers. This parameter allows refinement and ! composite grid combinations. ! integer :: NestLayers = 1 ! ! Number of grids in each nested layer. ! integer, allocatable :: GridsInLayer(:) ! [NestLayers] ! ! Application grid number as a function of the maximum number of grids ! in each nested layer and number of nested layers. ! integer, allocatable :: GridNumber(:,:) ! [Ngrids,NestLayers] ! ! Maximum automatic memory (bytes) required in distributed-memory ! buffers due to parallel exchanges. ! real(r8), allocatable :: BmemMax(:) ! [Ngrids] ! ! Estimated dynamic memory requirement (number of array elements) for ! a particular application per nested grid. ! real(r8), allocatable :: Dmem(:) ! [Ngrids] ! !----------------------------------------------------------------------- ! Lower and upper bounds indices per domain partition for all grids. !----------------------------------------------------------------------- ! ! Notice that these indices have different values in periodic and ! nesting applications, and on tiles next to the boundaries. Special ! indices are required to process overlap regions (suffices P and T) ! lateral boundary conditions (suffices B and M) in nested grid ! applications. The halo indices are used in private computations ! which include ghost-points and are limited by MAX/MIN functions ! on the tiles next to the model boundaries. For more information ! about all these indices, see routine "var_bounds" in file ! "Utility/get_bounds.F". ! ! All the 1D array indices are of size -1:NtileI(ng)*NtileJ(ng)-1. The ! -1 index include the values for the full (no partitions) grid. ! ! Notice that the starting (Imin, Jmin) and ending (Imax, Jmax) indices ! for I/O processing are 3D arrays. The first dimension (1:4) is for ! 1=PSI, 2=RHO, 3=u, 4=v points; the second dimension (0:1) is number ! of ghost points (0: no ghost points, 1: Nghost points), and the ! the third dimension is for 0:NtileI(ng)*NtileJ(ng)-1. ! TYPE T_BOUNDS integer, pointer :: tile(:) ! tile partition integer, pointer :: LBi(:) ! lower bound I-dimension integer, pointer :: UBi(:) ! upper bound I-dimension integer, pointer :: LBj(:) ! lower bound J-dimension integer, pointer :: UBj(:) ! upper bound J-dimension integer :: LBij ! lower bound MIN(I,J)-dimension integer :: UBij ! upper bound MAX(I,J)-dimension integer :: edge(4,4) ! boundary edges I- or J-indices integer, pointer :: Istr(:) ! starting tile I-direction integer, pointer :: Iend(:) ! ending tile I-direction integer, pointer :: Jstr(:) ! starting tile J-direction integer, pointer :: Jend(:) ! ending tile J-direction integer, pointer :: IstrR(:) ! starting tile I-direction (RHO) integer, pointer :: IendR(:) ! ending tile I-direction (RHO) integer, pointer :: IstrU(:) ! starting tile I-direction (U) integer, pointer :: JstrR(:) ! starting tile J-direction (RHO) integer, pointer :: JendR(:) ! ending tile J-direction (RHO) integer, pointer :: JstrV(:) ! starting tile J-direction (V) integer, pointer :: IstrB(:) ! starting obc I-direction (RHO,V) integer, pointer :: IendB(:) ! ending obc I-direction (RHO,V) integer, pointer :: IstrM(:) ! starting obc I-direction (PSI,U) integer, pointer :: JstrB(:) ! starting obc J-direction (RHO,U) integer, pointer :: JendB(:) ! ending obc J-direction (RHO,U) integer, pointer :: JstrM(:) ! starting obc J-direction (PSI,V) integer, pointer :: IstrP(:) ! starting nest I-direction (PSI,U) integer, pointer :: IendP(:) ! ending nest I-direction (PSI) integer, pointer :: JstrP(:) ! starting nest J-direction (PSI,V) integer, pointer :: JendP(:) ! ending nest J-direction (PSI) integer, pointer :: IstrT(:) ! starting nest I-direction (RHO) integer, pointer :: IendT(:) ! ending nest I-direction (RHO) integer, pointer :: JstrT(:) ! starting nest J-direction (RHO) integer, pointer :: JendT(:) ! ending nest J-direction (RHO) integer, pointer :: Istrm3(:) ! starting I-halo, Istr-3 integer, pointer :: Istrm2(:) ! starting I-halo, Istr-2 integer, pointer :: Istrm1(:) ! starting I-halo, Istr-1 integer, pointer :: IstrUm2(:) ! starting I-halo, IstrU-2 integer, pointer :: IstrUm1(:) ! starting I-halo, IstrU-1 integer, pointer :: Iendp1(:) ! ending I-halo, Iend+1 integer, pointer :: Iendp2(:) ! ending I-halo, Iend+2 integer, pointer :: Iendp2i(:) ! ending I-halo, Iend+2 integer, pointer :: Iendp3(:) ! ending I-halo, Iend+3 integer, pointer :: Jstrm3(:) ! starting J-halo, Jstr-3 integer, pointer :: Jstrm2(:) ! starting J-halo, Jstr-2 integer, pointer :: Jstrm1(:) ! starting J-halo, Jstr-1 integer, pointer :: JstrVm2(:) ! starting J-halo, JstrV-2 integer, pointer :: JstrVm1(:) ! starting J-halo, JstrV-1 integer, pointer :: Jendp1(:) ! ending J-halo, Jend+1 integer, pointer :: Jendp2(:) ! ending J-halo, Jend+2 integer, pointer :: Jendp2i(:) ! ending J-halo, Jend+2 integer, pointer :: Jendp3(:) ! ending J-halo, Jend+3 integer, pointer :: Imin(:,:,:) ! starting ghost I-direction integer, pointer :: Imax(:,:,:) ! ending ghost I-direction integer, pointer :: Jmin(:,:,:) ! starting ghost J-direction integer, pointer :: Jmax(:,:,:) ! ending ghost J-direction END TYPE T_BOUNDS TYPE (T_BOUNDS), allocatable :: BOUNDS(:) ! !----------------------------------------------------------------------- ! Lower and upper bounds in NetCDF files. !----------------------------------------------------------------------- ! TYPE T_IOBOUNDS integer :: ILB_psi ! I-direction lower bound (PSI) integer :: IUB_psi ! I-direction upper bound (PSI) integer :: JLB_psi ! J-direction lower bound (PSI) integer :: JUB_psi ! J-direction upper bound (PSI) integer :: ILB_rho ! I-direction lower bound (RHO) integer :: IUB_rho ! I-direction upper bound (RHO) integer :: JLB_rho ! J-direction lower bound (RHO) integer :: JUB_rho ! J-direction upper bound (RHO) integer :: ILB_u ! I-direction lower bound (U) integer :: IUB_u ! I-direction upper bound (U) integer :: JLB_u ! J-direction lower bound (U) integer :: JUB_u ! J-direction upper bound (U) integer :: ILB_v ! I-direction lower bound (V) integer :: IUB_v ! I-direction upper bound (V) integer :: JLB_v ! J-direction lower bound (V) integer :: JUB_v ! J-direction upper bound (V) integer :: IorJ ! number of MAX(I,J)-direction points integer :: xi_psi ! number of I-direction points (PSI) integer :: xi_rho ! number of I-direction points (RHO) integer :: xi_u ! number of I-direction points (U) integer :: xi_v ! number of I-direction points (V) integer :: eta_psi ! number of J-direction points (PSI) integer :: eta_rho ! number of J-direction points (RHO) integer :: eta_u ! number of I-direction points (U) integer :: eta_v ! number of I-direction points (V) #if (defined READ_WATER || defined WRITE_WATER) && defined MASKING integer :: xy_psi ! number of IJ-water points (PSI) integer :: xy_rho ! number of IJ-water points (RHO) integer :: xy_u ! number of IJ-water points (U) integer :: xy_v ! number of IJ-water points (V) #endif END TYPE T_IOBOUNDS TYPE (T_IOBOUNDS), allocatable :: IOBOUNDS(:) ! !----------------------------------------------------------------------- ! Domain boundary edges switches and tiles minimum and maximum ! fractional grid coordinates. !----------------------------------------------------------------------- ! TYPE T_DOMAIN logical, pointer :: Eastern_Edge(:) logical, pointer :: Western_Edge(:) logical, pointer :: Northern_Edge(:) logical, pointer :: Southern_Edge(:) logical, pointer :: NorthEast_Corner(:) logical, pointer :: NorthWest_Corner(:) logical, pointer :: SouthEast_Corner(:) logical, pointer :: SouthWest_Corner(:) logical, pointer :: NorthEast_Test(:) logical, pointer :: NorthWest_Test(:) logical, pointer :: SouthEast_Test(:) logical, pointer :: SouthWest_Test(:) real(r8), pointer :: Xmin_psi(:) real(r8), pointer :: Xmax_psi(:) real(r8), pointer :: Ymin_psi(:) real(r8), pointer :: Ymax_psi(:) real(r8), pointer :: Xmin_rho(:) real(r8), pointer :: Xmax_rho(:) real(r8), pointer :: Ymin_rho(:) real(r8), pointer :: Ymax_rho(:) real(r8), pointer :: Xmin_u(:) real(r8), pointer :: Xmax_u(:) real(r8), pointer :: Ymin_u(:) real(r8), pointer :: Ymax_u(:) real(r8), pointer :: Xmin_v(:) real(r8), pointer :: Xmax_v(:) real(r8), pointer :: Ymin_v(:) real(r8), pointer :: Ymax_v(:) END TYPE T_DOMAIN TYPE (T_DOMAIN), allocatable :: DOMAIN(:) ! !----------------------------------------------------------------------- ! Lateral Boundary Conditions (LBC) switches structure. !----------------------------------------------------------------------- ! ! The lateral boundary conditions are specified by logical switches ! to facilitate applications with nested grids. It also allows to ! set different boundary conditions between the nonlinear model and ! the adjoint/tangent models. The LBC structure is allocated as: ! ! LBC(1:4, nLBCvar, Ngrids) ! ! where 1:4 are the number boundary edges, nLBCvar are the number ! LBC state variables, and Ngrids is the number of nested grids. ! For example, for free-surface gradient boundary conditions we ! have: ! ! LBC(iwest, isFsur, ng) % gradient ! LBC(ieast, isFsur, ng) % gradient ! LBC(isouth, isFsur, ng) % gradient ! LBC(inorth, isFsur, ng) % gradient ! integer :: nLBCvar ! TYPE T_LBC logical :: acquire ! process lateral boundary data logical :: Chapman_explicit logical :: Chapman_implicit logical :: clamped logical :: closed logical :: Flather logical :: gradient logical :: mixed logical :: nested logical :: nudging logical :: periodic logical :: radiation logical :: reduced logical :: Shchepetkin END TYPE T_LBC TYPE (T_LBC), allocatable :: LBC(:,:,:) #if defined ADJOINT || defined TANGENT || defined TL_IOMS TYPE (T_LBC), allocatable :: ad_LBC(:,:,:) TYPE (T_LBC), allocatable :: tl_LBC(:,:,:) #endif #ifdef SOLVE3D ! !----------------------------------------------------------------------- ! Tracer horizontal and vertical advection switches structure. !----------------------------------------------------------------------- ! ! A different advection scheme is allowed for each tracer. For ! example, a positive-definite algorithm can be activated for salinity, ! biological, and sediment tracers, while a different one is set for ! temperature. ! TYPE T_ADV logical :: AKIMA4 ! fourth-order Akima logical :: CENTERED2 ! second-order centered differences logical :: CENTERED4 ! fourth-order centered differences logical :: HSIMT ! third-order HSIMT-TVD logical :: MPDATA ! recursive flux corrected MPDATA logical :: SPLINES ! parabolic splines (vertical only) logical :: SPLIT_U3 ! split third-order upstream logical :: UPSTREAM3 ! third-order upstream bias END TYPE T_ADV TYPE (T_ADV), allocatable :: Hadvection(:,:) ! horizontal TYPE (T_ADV), allocatable :: Vadvection(:,:) ! vertical # if defined ADJOINT || defined TANGENT || defined TL_IOMS TYPE (T_ADV), allocatable :: ad_Hadvection(:,:) ! horizontal TYPE (T_ADV), allocatable :: ad_Vadvection(:,:) ! vertical # endif # if defined TANGENT || defined TL_IOMS TYPE (T_ADV), allocatable :: tl_Hadvection(:,:) ! horizontal TYPE (T_ADV), allocatable :: tl_Vadvection(:,:) ! vertical # endif #endif ! !----------------------------------------------------------------------- ! Field statistics STATS structure. !----------------------------------------------------------------------- ! TYPE T_STATS integer(i8b) :: checksum ! bit sum hash integer :: count ! processed values count real(r8) :: min ! minimum value real(r8) :: max ! maximum value real(r8) :: avg ! arithmetic mean real(r8) :: rms ! root mean square END TYPE T_STATS #ifdef PROPAGATOR ! !----------------------------------------------------------------------- ! Generalized Stability Theory (GST) nested grid pointers for the ! state vector. Mostly used for memory management with the ARPACK ! library. !----------------------------------------------------------------------- ! TYPE T_GST real(r8), pointer :: vector(:) END TYPE #endif ! !----------------------------------------------------------------------- ! Model grid(s) parameters. !----------------------------------------------------------------------- ! ! Number of interior RHO-points in the XI- and ETA-directions. The ! size of models state variables (C-grid) at input and output are: ! ! RH0-type variables: [0:Lm+1, 0:Mm+1] ----v(i,j+1)---- ! PSI-type variables: [1:Lm+1, 1:Mm+1] | | ! U-type variables: [1:Lm+1, 0:Mm+1] u(i,j) r(i,j) u(i+1,j) ! V-type variables: [0:Lm+1, 1:Mm+1] | | ! -----v(i,j)----- integer, allocatable :: Lm(:) integer, allocatable :: Mm(:) ! ! Global horizontal size of model arrays including padding. All the ! model state arrays are of same size to facilitate parallelization. ! integer, allocatable :: Im(:) integer, allocatable :: Jm(:) ! ! Number of vertical levels. The vertical ranges of model state ! variables are: ! -----W(i,j,k)----- ! RHO-, U-, V-type variables: [1:N] | | ! W-type variables: [0:N] | r(i,j,k) | ! | | ! ----W(i,j,k-1)---- integer, allocatable :: N(:) ! !----------------------------------------------------------------------- ! Tracers parameters. !----------------------------------------------------------------------- ! ! Total number of tracer type variables, NT(:) = NAT + NBT + NPT + NST. ! The MT corresponds to the maximum number of tracers between all ! nested grids. ! integer, allocatable :: NT(:) integer :: MT ! ! Total number of climatology tracer type variables to process. ! integer, allocatable :: NTCLM(:) ! ! Number of active tracers. Usually, NAT=2 for potential temperature ! and salinity. ! integer :: NAT = 0 ! ! Total number of inert passive tracers to advect and diffuse only ! (like dyes, etc). This parameter is independent of the number of ! biological and/or sediment tracers. ! integer :: NPT = 0 ! ! Number of biological tracers. ! integer :: NBT = 0 ! !----------------------------------------------------------------------- ! Sediment tracers parameters. !----------------------------------------------------------------------- ! ! Number of sediment bed layes. ! integer :: Nbed = 0 ! ! Total number of sediment tracers, NST = NCS + NNS. ! integer :: NST = 0 ! ! Number of cohesive (mud) sediments. ! integer :: NCS = 0 ! ! Number of non-cohesive (sand) sediments. ! integer :: NNS = 0 #ifdef FLOATS ! !----------------------------------------------------------------------- ! Floats trajectories parameters. !----------------------------------------------------------------------- ! ! Number of trajectory time-stepping levels, [0:NFT]. ! integer, parameter :: NFT = 4 ! ! Total number of floats to track. ! integer, allocatable :: Nfloats(:) ! ! Total number of float variables to process and output. ! integer, allocatable :: NFV(:) #endif #ifdef STATIONS ! !----------------------------------------------------------------------- ! Stations parameters. !----------------------------------------------------------------------- ! ! Number of output stations. ! integer, allocatable :: Nstation(:) #endif ! !----------------------------------------------------------------------- ! Maximum number of tidal constituents to process. !----------------------------------------------------------------------- ! integer :: MTC #ifdef DIAGNOSTICS ! !----------------------------------------------------------------------- ! Diagnostic fields parameters. !----------------------------------------------------------------------- ! ! Number of diagnostic tracer fields. ! integer :: NDT ! ! Number of diagnostic momentum fields. ! integer :: NDM2d ! 2D momentum integer :: NDM3d ! 3D momentum ! ! Number of diagnostic biology/bio-optical fields. Currently, only ! available for the Fennel and EcoSim models. ! integer :: NDbio2d ! 2D fields integer :: NDbio3d ! 3D fields integer :: NDbio4d ! 4D fields ! ! Number of diagnostic 3D right-hand-side fields. ! integer :: NDrhs # if defined ECOSIM && defined BIO_OPTIC ! ! Number of light spectral irradiance bands to process. ! integer :: NDbands # endif #endif #if defined FOUR_DVAR || defined IMPULSE ! !----------------------------------------------------------------------- ! 4DVAR parameters. !----------------------------------------------------------------------- ! ! Number of state arrays (size of additional dimension) for error ! covariance normalization and standard deviation factors. This ! paramenter is NSA=1 for strong constraint and NSA=2 for weak ! constraint. In weak constraint, it contains fields for initial ! conditions and model error covariance states. ! integer :: NSA #endif #if defined BALANCE_OPERATOR && defined ZETA_ELLIPTIC ! ! Number of biconjugate gradient algorithm iteractions. ! integer, allocatable :: Nbico(:) #endif #if defined TIME_CONV && defined WEAK_CONSTRAINT ! ! Number of full state records needed for time convolutions. It is ! assign to (ntimes(ng)/nADJ(ng))+1 in "read_phypar.F" ! integer, allocatable :: NrecTC(:) #endif ! !----------------------------------------------------------------------- ! Model state parameters. !----------------------------------------------------------------------- ! ! Number of model state variables. ! integer, allocatable :: NSV(:) #ifdef PROPAGATOR ! ! Size of the full packed state vector (water points only), its ! parallel node partition, and associated node starting and ending ! indices. ! integer, allocatable :: Mstate(:) ! full state integer, allocatable :: Nstate(:) ! partioned state integer, allocatable :: Nstr(:) ! starting index integer, allocatable :: Nend(:) ! ending index integer, allocatable :: Nsize(:) ! Nend-Nstr+1 # ifdef SO_SEMI ! ! Number of time record saved in seminorm strochastic optimals adjoint ! state vector. ! integer, allocatable :: Nsemi(:) # endif #endif ! ! Set nonlinear, tangent linear, representer (finite-amplitude tangent ! linear) and adjoint models identifiers. ! integer, parameter :: iNLM = 1 ! nonlinear integer, parameter :: iTLM = 2 ! perturbation tangent linear integer, parameter :: iRPM = 3 ! finite-amplitude tangent linear integer, parameter :: iADM = 4 ! adjoint ! character (len=3), dimension(4) :: KernelString = & & (/ 'NLM','TLM','RPM','ADM' /) ! !----------------------------------------------------------------------- ! Domain partition parameters. !----------------------------------------------------------------------- ! ! Number of tiles or domain partitions in the XI- and ETA-directions. ! These values are used to compute tile ranges [Istr:Iend, Jstr:Jend]. ! integer, allocatable :: NtileI(:) integer, allocatable :: NtileJ(:) ! ! Number of tiles or domain partitions in the XI- and ETA-directions. ! These values are used to parallel loops to differentiate between ! shared-memory and distributed-memory. Notice that in distributed ! memory both values are set to one. ! integer, allocatable :: NtileX(:) integer, allocatable :: NtileE(:) ! ! Maximun number of points in the halo region for exchanging state ! boundary arrays during convolutions. ! integer, allocatable :: HaloBry(:) ! ! Maximum number of points in the halo region in the XI- and ! ETA-directions. ! integer, allocatable :: HaloSizeI(:) integer, allocatable :: HaloSizeJ(:) ! ! Maximum tile side length in XI- or ETA-directions. ! integer, allocatable :: TileSide(:) ! ! Maximum number of points in a tile partition. ! integer, allocatable :: TileSize(:) ! ! Set number of ghost-points in the halo region. It is used in ! periodic, nested, and distributed-memory applications. ! integer :: NghostPoints ! !----------------------------------------------------------------------- ! Staggered C-grid location identifiers. !----------------------------------------------------------------------- ! integer, parameter :: p2dvar = 1 ! 2D PSI-variable integer, parameter :: r2dvar = 2 ! 2D RHO-variable integer, parameter :: u2dvar = 3 ! 2D U-variable integer, parameter :: v2dvar = 4 ! 2D V-variable integer, parameter :: p3dvar = 5 ! 3D PSI-variable integer, parameter :: r3dvar = 6 ! 3D RHO-variable integer, parameter :: u3dvar = 7 ! 3D U-variable integer, parameter :: v3dvar = 8 ! 3D V-variable integer, parameter :: w3dvar = 9 ! 3D W-variable integer, parameter :: b3dvar = 10 ! 3D BED-sediment integer, parameter :: l3dvar = 11 ! 3D spectral light integer, parameter :: l4dvar = 12 ! 4D spectral light integer, parameter :: r2dobc = 13 ! 2D OBC RHO-variable integer, parameter :: u2dobc = 14 ! 2D OBC U-variable integer, parameter :: v2dobc = 15 ! 2D OBC V-variable integer, parameter :: r3dobc = 16 ! 3D OBC RHO-variable integer, parameter :: u3dobc = 17 ! 3D OBC U-variable integer, parameter :: v3dobc = 18 ! 3D OBC V-variable ! CONTAINS ! SUBROUTINE allocate_param ! !======================================================================= ! ! ! This routine allocates several variables in the module that depend ! ! on the number of nested grids. ! ! ! !======================================================================= ! !----------------------------------------------------------------------- ! Allocate dimension parameters. !----------------------------------------------------------------------- ! IF (.not.allocated(Lm)) THEN allocate ( Lm(Ngrids) ) END IF IF (.not.allocated(Mm)) THEN allocate ( Mm(Ngrids) ) END IF IF (.not.allocated(Im)) THEN allocate ( Im(Ngrids) ) END IF IF (.not.allocated(Jm)) THEN allocate ( Jm(Ngrids) ) END IF IF (.not.allocated(N)) THEN allocate ( N(Ngrids) ) END IF IF (.not.allocated(NT)) THEN allocate ( NT(Ngrids) ) END IF #if defined TIME_CONV && defined WEAK_CONSTRAINT IF (.not.allocated(NrecTC)) THEN allocate ( NrecTC(Ngrids) ) END IF #endif IF (.not.allocated(NTCLM)) THEN allocate ( NTCLM(Ngrids) ) END IF IF (.not.allocated(NSV)) THEN allocate ( NSV(Ngrids) ) END IF IF (.not.allocated(NtileI)) THEN allocate ( NtileI(Ngrids) ) END IF IF (.not.allocated(NtileJ)) THEN allocate ( NtileJ(Ngrids) ) END IF IF (.not.allocated(NtileX)) THEN allocate ( NtileX(Ngrids) ) END IF IF (.not.allocated(NtileE)) THEN allocate ( NtileE(Ngrids) ) END IF IF (.not.allocated(HaloBry)) THEN allocate ( HaloBry(Ngrids) ) END IF IF (.not.allocated(HaloSizeI)) THEN allocate ( HaloSizeI(Ngrids) ) END IF IF (.not.allocated(HaloSizeJ)) THEN allocate ( HaloSizeJ(Ngrids) ) END IF IF (.not.allocated(TileSide)) THEN allocate ( TileSide(Ngrids) ) END IF IF (.not.allocated(TileSize)) THEN allocate ( TileSize(Ngrids) ) END IF IF (.not.allocated(BmemMax)) THEN allocate ( BmemMax(Ngrids) ) BmemMax=0.0_r8 END IF IF (.not.allocated(Dmem)) THEN allocate ( Dmem(Ngrids) ) Dmem=0.0_r8 END IF #ifdef PROPAGATOR IF (.not.allocated(Mstate)) THEN allocate ( Mstate(Ngrids) ) END IF IF (.not.allocated(Nstate)) THEN allocate ( Nstate(Ngrids) ) END IF IF (.not.allocated(Nstr)) THEN allocate ( Nstr(Ngrids) ) END IF IF (.not.allocated(Nend)) THEN allocate ( Nend(Ngrids) ) END IF IF (.not.allocated(Nsize)) THEN allocate ( Nsize(Ngrids) ) END IF # ifdef SO_SEMI IF (.not.allocated(Nsemi)) THEN allocate ( Nsemi(Ngrids) ) END IF # endif #endif #ifdef FLOATS IF (.not.allocated(Nfloats)) THEN allocate ( Nfloats(Ngrids) ) END IF IF (.not.allocated(NFV)) THEN allocate ( NFV(Ngrids) ) END IF #endif #ifdef STATIONS IF (.not.allocated(Nstation)) THEN allocate ( Nstation(Ngrids) ) END IF #endif #if defined BALANCE_OPERATOR && defined ZETA_ELLIPTIC IF (.not.allocated(Nbico)) THEN allocate ( Nbico(Ngrids) ) END IF #endif ! RETURN END SUBROUTINE allocate_param ! SUBROUTINE deallocate_param ! !======================================================================= ! ! ! This routine deallocates som of the variables in the module. ! ! Notice that "destroy" cannot be use to deallocate the pointer ! ! variables because of cyclic dependencies. ! ! ! !======================================================================= ! ! Local variable declarations. ! integer :: ng #ifdef SUBOBJECT_DEALLOCATION ! !----------------------------------------------------------------------- ! Deallocate each variable in the derived-type T_BOUNDS structure ! separately. !----------------------------------------------------------------------- ! DO ng=1,Ngrids IF (associated(BOUNDS(ng)%tile)) THEN deallocate ( BOUNDS(ng)%tile ) END IF ! IF (associated(BOUNDS(ng)%LBi)) THEN deallocate ( BOUNDS(ng)%LBi ) END IF IF (associated(BOUNDS(ng)%UBi)) THEN deallocate ( BOUNDS(ng)%UBi ) END IF IF (associated(BOUNDS(ng)%LBj)) THEN deallocate ( BOUNDS(ng)%LBj ) END IF IF (associated(BOUNDS(ng)%UBj)) THEN deallocate ( BOUNDS(ng)%UBj ) END IF ! IF (associated(BOUNDS(ng)%Istr)) THEN deallocate ( BOUNDS(ng)%Istr ) END IF IF (associated(BOUNDS(ng)%Iend)) THEN deallocate ( BOUNDS(ng)%Iend ) END IF IF (associated(BOUNDS(ng)%Jstr)) THEN deallocate ( BOUNDS(ng)%Jstr ) END IF IF (associated(BOUNDS(ng)%Jend)) THEN deallocate ( BOUNDS(ng)%Jend ) END IF ! IF (associated(BOUNDS(ng)%IstrR)) THEN deallocate ( BOUNDS(ng)%IstrR ) END IF IF (associated(BOUNDS(ng)%IendR)) THEN deallocate ( BOUNDS(ng)%IendR ) END IF IF (associated(BOUNDS(ng)%IstrU)) THEN deallocate ( BOUNDS(ng)%IstrU ) END IF ! IF (associated(BOUNDS(ng)%JstrR)) THEN deallocate ( BOUNDS(ng)%JstrR ) END IF IF (associated(BOUNDS(ng)%JendR)) THEN deallocate ( BOUNDS(ng)%JendR ) END IF IF (associated(BOUNDS(ng)%JstrV)) THEN deallocate ( BOUNDS(ng)%JstrV ) END IF ! IF (associated(BOUNDS(ng)%IstrB)) THEN deallocate ( BOUNDS(ng)%IstrB ) END IF IF (associated(BOUNDS(ng)%IendB)) THEN deallocate ( BOUNDS(ng)%IendB ) END IF IF (associated(BOUNDS(ng)%IstrM)) THEN deallocate ( BOUNDS(ng)%IstrM ) END IF ! IF (associated(BOUNDS(ng)%JstrB)) THEN deallocate ( BOUNDS(ng)%JstrB ) END IF IF (associated(BOUNDS(ng)%JendB)) THEN deallocate ( BOUNDS(ng)%JendB ) END IF IF (associated(BOUNDS(ng)%JstrM)) THEN deallocate ( BOUNDS(ng)%JstrM ) END IF ! IF (associated(BOUNDS(ng)%IstrP)) THEN deallocate ( BOUNDS(ng)%IstrP ) END IF IF (associated(BOUNDS(ng)%IendP)) THEN deallocate ( BOUNDS(ng)%IendP ) END IF IF (associated(BOUNDS(ng)%JstrP)) THEN deallocate ( BOUNDS(ng)%JstrP ) END IF IF (associated(BOUNDS(ng)%JendP)) THEN deallocate ( BOUNDS(ng)%JendP ) END IF ! IF (associated(BOUNDS(ng)%IstrT)) THEN deallocate ( BOUNDS(ng)%IstrT ) END IF IF (associated(BOUNDS(ng)%IendT)) THEN deallocate ( BOUNDS(ng)%IendT ) END IF IF (associated(BOUNDS(ng)%JstrT)) THEN deallocate ( BOUNDS(ng)%JstrT ) END IF IF (associated(BOUNDS(ng)%JendT)) THEN deallocate ( BOUNDS(ng)%JendT ) END IF ! IF (associated(BOUNDS(ng)%Istrm3)) THEN deallocate ( BOUNDS(ng)%Istrm3 ) END IF IF (associated(BOUNDS(ng)%Istrm2)) THEN deallocate ( BOUNDS(ng)%Istrm2 ) END IF IF (associated(BOUNDS(ng)%Istrm1)) THEN deallocate ( BOUNDS(ng)%Istrm1 ) END IF IF (associated(BOUNDS(ng)%IstrUm2)) THEN deallocate ( BOUNDS(ng)%IstrUm2 ) END IF IF (associated(BOUNDS(ng)%IstrUm1)) THEN deallocate ( BOUNDS(ng)%IstrUm1 ) END IF IF (associated(BOUNDS(ng)%Iendp1)) THEN deallocate ( BOUNDS(ng)%Iendp1 ) END IF IF (associated(BOUNDS(ng)%Iendp2)) THEN deallocate ( BOUNDS(ng)%Iendp2 ) END IF IF (associated(BOUNDS(ng)%Iendp2i)) THEN deallocate ( BOUNDS(ng)%Iendp2i ) END IF IF (associated(BOUNDS(ng)%Iendp3)) THEN deallocate ( BOUNDS(ng)%Iendp3 ) END IF ! IF (associated(BOUNDS(ng)%Jstrm3)) THEN deallocate ( BOUNDS(ng)%Jstrm3 ) END IF IF (associated(BOUNDS(ng)%Jstrm2)) THEN deallocate ( BOUNDS(ng)%Jstrm2 ) END IF IF (associated(BOUNDS(ng)%Jstrm1)) THEN deallocate ( BOUNDS(ng)%Jstrm1 ) END IF IF (associated(BOUNDS(ng)%JstrVm2)) THEN deallocate ( BOUNDS(ng)%JstrVm2 ) END IF IF (associated(BOUNDS(ng)%JstrVm1)) THEN deallocate ( BOUNDS(ng)%JstrVm1 ) END IF IF (associated(BOUNDS(ng)%Jendp1)) THEN deallocate ( BOUNDS(ng)%Jendp1 ) END IF IF (associated(BOUNDS(ng)%Jendp2)) THEN deallocate ( BOUNDS(ng)%Jendp2 ) END IF IF (associated(BOUNDS(ng)%Jendp2i)) THEN deallocate ( BOUNDS(ng)%Jendp2i ) END IF IF (associated(BOUNDS(ng)%Jendp3)) THEN deallocate ( BOUNDS(ng)%Jendp3 ) END IF ! IF (associated(BOUNDS(ng)%Imin)) THEN deallocate ( BOUNDS(ng)%Imin ) END IF IF (associated(BOUNDS(ng)%Imax)) THEN deallocate ( BOUNDS(ng)%Imax ) END IF IF (associated(BOUNDS(ng)%Jmin)) THEN deallocate ( BOUNDS(ng)%Jmin ) END IF IF (associated(BOUNDS(ng)%Jmax)) THEN deallocate ( BOUNDS(ng)%Jmax ) END IF END DO ! !----------------------------------------------------------------------- ! Deallocate each variable in the derived-type T_DOMAIN structure ! separately. !----------------------------------------------------------------------- ! DO ng=1,Ngrids IF (associated(DOMAIN(ng)%Eastern_Edge)) THEN deallocate ( DOMAIN(ng)%Eastern_Edge ) END IF IF (associated(DOMAIN(ng)%Western_Edge)) THEN deallocate ( DOMAIN(ng)%Western_Edge ) END IF IF (associated(DOMAIN(ng)%Northern_Edge)) THEN deallocate ( DOMAIN(ng)%Northern_Edge ) END IF IF (associated(DOMAIN(ng)%Southern_Edge)) THEN deallocate ( DOMAIN(ng)%Southern_Edge ) END IF ! IF (associated(DOMAIN(ng)%NorthEast_Corner)) THEN deallocate ( DOMAIN(ng)%NorthEast_Corner ) END IF IF (associated(DOMAIN(ng)%NorthWest_Corner)) THEN deallocate ( DOMAIN(ng)%NorthWest_Corner ) END IF IF (associated(DOMAIN(ng)%SouthEast_Corner)) THEN deallocate ( DOMAIN(ng)%SouthEast_Corner ) END IF IF (associated(DOMAIN(ng)%SouthWest_Corner)) THEN deallocate ( DOMAIN(ng)%SouthWest_Corner ) END IF ! IF (associated(DOMAIN(ng)%NorthEast_Test)) THEN deallocate ( DOMAIN(ng)%NorthEast_Test ) END IF IF (associated(DOMAIN(ng)%NorthWest_Test)) THEN deallocate ( DOMAIN(ng)%NorthWest_Test ) END IF IF (associated(DOMAIN(ng)%SouthEast_Test)) THEN deallocate ( DOMAIN(ng)%SouthEast_Test ) END IF IF (associated(DOMAIN(ng)%SouthWest_Test)) THEN deallocate ( DOMAIN(ng)%SouthWest_Test ) END IF ! IF (associated(DOMAIN(ng)%Xmin_psi)) THEN deallocate ( DOMAIN(ng)%Xmin_psi ) END IF IF (associated(DOMAIN(ng)%Xmax_psi)) THEN deallocate ( DOMAIN(ng)%Xmax_psi ) END IF IF (associated(DOMAIN(ng)%Ymin_psi)) THEN deallocate ( DOMAIN(ng)%Ymin_psi ) END IF IF (associated(DOMAIN(ng)%Ymax_psi)) THEN deallocate ( DOMAIN(ng)%Ymax_psi ) END IF ! IF (associated(DOMAIN(ng)%Xmin_rho)) THEN deallocate ( DOMAIN(ng)%Xmin_rho ) END IF IF (associated(DOMAIN(ng)%Xmax_rho)) THEN deallocate ( DOMAIN(ng)%Xmax_rho ) END IF IF (associated(DOMAIN(ng)%Ymin_rho)) THEN deallocate ( DOMAIN(ng)%Ymin_rho ) END IF IF (associated(DOMAIN(ng)%Ymax_rho)) THEN deallocate ( DOMAIN(ng)%Ymax_rho ) END IF ! IF (associated(DOMAIN(ng)%Xmin_u)) THEN deallocate ( DOMAIN(ng)%Xmin_u ) END IF IF (associated(DOMAIN(ng)%Xmax_u)) THEN deallocate ( DOMAIN(ng)%Xmax_u ) END IF IF (associated(DOMAIN(ng)%Ymin_u)) THEN deallocate ( DOMAIN(ng)%Ymin_u ) END IF IF (associated(DOMAIN(ng)%Ymax_u)) THEN deallocate ( DOMAIN(ng)%Ymax_u ) END IF ! IF (associated(DOMAIN(ng)%Xmin_v)) THEN deallocate ( DOMAIN(ng)%Xmin_v ) END IF IF (associated(DOMAIN(ng)%Xmax_v)) THEN deallocate ( DOMAIN(ng)%Xmax_v ) END IF IF (associated(DOMAIN(ng)%Ymin_v)) THEN deallocate ( DOMAIN(ng)%Ymin_v ) END IF IF (associated(DOMAIN(ng)%Ymax_v)) THEN deallocate ( DOMAIN(ng)%Ymax_v ) END IF END DO #endif ! !----------------------------------------------------------------------- ! Deallocate derived-type module structures: !----------------------------------------------------------------------- ! IF (allocated(BOUNDS)) deallocate ( BOUNDS ) IF (allocated(IOBOUNDS)) deallocate ( IOBOUNDS ) IF (allocated(DOMAIN)) deallocate ( DOMAIN ) #ifdef SOLVE3D IF (allocated(Hadvection)) deallocate ( Hadvection ) IF (allocated(Vadvection)) deallocate ( Vadvection ) # if defined ADJOINT || defined TANGENT || defined TL_IOMS IF (allocated(ad_Hadvection)) deallocate ( ad_Hadvection ) IF (allocated(ad_Vadvection)) deallocate ( ad_Vadvection ) # endif # if defined TANGENT || defined TL_IOMS IF (allocated(tl_Hadvection)) deallocate ( tl_Hadvection ) IF (allocated(tl_Vadvection)) deallocate ( tl_Vadvection ) # endif #endif ! !----------------------------------------------------------------------- ! Deallocate dimension parameters. !----------------------------------------------------------------------- ! IF (allocated(Lm)) deallocate ( Lm ) IF (allocated(Mm)) deallocate ( Mm ) IF (allocated(Im)) deallocate ( Im ) IF (allocated(Jm)) deallocate ( Jm ) IF (allocated(N)) deallocate ( N ) IF (allocated(NT)) deallocate ( NT ) #if defined TIME_CONV && defined WEAK_CONSTRAINT IF (allocated(NrecTC)) deallocate ( NrecTC ) #endif IF (allocated(NTCLM)) deallocate ( NTCLM ) IF (allocated(NSV)) deallocate ( NSV ) IF (allocated(NtileI)) deallocate ( NtileI ) IF (allocated(NtileJ)) deallocate ( NtileJ ) IF (allocated(NtileX)) deallocate ( NtileX ) IF (allocated(NtileE)) deallocate ( NtileE ) IF (allocated(HaloBry)) deallocate ( HaloBry ) IF (allocated(HaloSizeI)) deallocate ( HaloSizeI ) IF (allocated(HaloSizeJ)) deallocate ( HaloSizeJ ) IF (allocated(TileSide)) deallocate ( TileSide ) IF (allocated(TileSize)) deallocate ( TileSize ) IF (allocated(BmemMax)) deallocate ( BmemMax ) IF (allocated(Dmem)) deallocate ( Dmem ) #ifdef PROPAGATOR IF (allocated(Mstate)) deallocate ( Mstate ) IF (allocated(Nstate)) deallocate ( Nstate ) IF (allocated(Nstr)) deallocate ( Nstr ) IF (allocated(Nend)) deallocate ( Nend ) IF (allocated(Nsize)) deallocate ( Nsize ) # ifdef SO_SEMI IF (allocated(Nsemi)) deallocate ( Nsemi ) # endif #endif IF (allocated(GridsInLayer)) deallocate ( GridsInLayer ) IF (allocated(GridNumber)) deallocate ( GridNumber ) #ifdef FLOATS IF (allocated(Nfloats)) deallocate ( Nfloats ) IF (allocated(NFV)) deallocate ( NFV ) #endif #ifdef STATIONS IF (allocated(Nstation)) deallocate ( Nstation ) #endif #if defined BALANCE_OPERATOR && defined ZETA_ELLIPTIC IF (allocated(Nbico)) deallocate ( Nbico ) #endif ! RETURN END SUBROUTINE deallocate_param ! SUBROUTINE initialize_param ! !======================================================================= ! ! ! This routine initializes several parameters in module "mod_param" ! ! for all nested grids. ! ! ! !======================================================================= ! ! Local variable declarations ! integer :: I_padd, J_padd, Ntiles integer :: ibry, itrc, ivar, ng ! !----------------------------------------------------------------------- ! Now that we know the values for the tile partitions (NtileI,NtileJ), ! allocate module structures. !----------------------------------------------------------------------- ! ! Allocate lower and upper bounds indices structure. ! IF (.not.allocated(BOUNDS)) THEN allocate ( BOUNDS(Ngrids) ) DO ng=1,Ngrids Ntiles=NtileI(ng)*NtileJ(ng)-1 allocate ( BOUNDS(ng) % tile(-1:Ntiles) ) allocate ( BOUNDS(ng) % LBi(-1:Ntiles) ) allocate ( BOUNDS(ng) % UBi(-1:Ntiles) ) allocate ( BOUNDS(ng) % LBj(-1:Ntiles) ) allocate ( BOUNDS(ng) % UBj(-1:Ntiles) ) allocate ( BOUNDS(ng) % Istr(-1:Ntiles) ) allocate ( BOUNDS(ng) % Iend(-1:Ntiles) ) allocate ( BOUNDS(ng) % Jstr(-1:Ntiles) ) allocate ( BOUNDS(ng) % Jend(-1:Ntiles) ) allocate ( BOUNDS(ng) % IstrR(-1:Ntiles) ) allocate ( BOUNDS(ng) % IendR(-1:Ntiles) ) allocate ( BOUNDS(ng) % IstrU(-1:Ntiles) ) allocate ( BOUNDS(ng) % JstrR(-1:Ntiles) ) allocate ( BOUNDS(ng) % JendR(-1:Ntiles) ) allocate ( BOUNDS(ng) % JstrV(-1:Ntiles) ) allocate ( BOUNDS(ng) % IstrB(-1:Ntiles) ) allocate ( BOUNDS(ng) % IendB(-1:Ntiles) ) allocate ( BOUNDS(ng) % IstrM(-1:Ntiles) ) allocate ( BOUNDS(ng) % JstrB(-1:Ntiles) ) allocate ( BOUNDS(ng) % JendB(-1:Ntiles) ) allocate ( BOUNDS(ng) % JstrM(-1:Ntiles) ) allocate ( BOUNDS(ng) % IstrP(-1:Ntiles) ) allocate ( BOUNDS(ng) % IendP(-1:Ntiles) ) allocate ( BOUNDS(ng) % IstrT(-1:Ntiles) ) allocate ( BOUNDS(ng) % IendT(-1:Ntiles) ) allocate ( BOUNDS(ng) % JstrP(-1:Ntiles) ) allocate ( BOUNDS(ng) % JendP(-1:Ntiles) ) allocate ( BOUNDS(ng) % JstrT(-1:Ntiles) ) allocate ( BOUNDS(ng) % JendT(-1:Ntiles) ) allocate ( BOUNDS(ng) % Istrm3(-1:Ntiles) ) allocate ( BOUNDS(ng) % Istrm2(-1:Ntiles) ) allocate ( BOUNDS(ng) % Istrm1(-1:Ntiles) ) allocate ( BOUNDS(ng) % IstrUm2(-1:Ntiles) ) allocate ( BOUNDS(ng) % IstrUm1(-1:Ntiles) ) allocate ( BOUNDS(ng) % Iendp1(-1:Ntiles) ) allocate ( BOUNDS(ng) % Iendp2(-1:Ntiles) ) allocate ( BOUNDS(ng) % Iendp2i(-1:Ntiles) ) allocate ( BOUNDS(ng) % Iendp3(-1:Ntiles) ) allocate ( BOUNDS(ng) % Jstrm3(-1:Ntiles) ) allocate ( BOUNDS(ng) % Jstrm2(-1:Ntiles) ) allocate ( BOUNDS(ng) % Jstrm1(-1:Ntiles) ) allocate ( BOUNDS(ng) % JstrVm2(-1:Ntiles) ) allocate ( BOUNDS(ng) % JstrVm1(-1:Ntiles) ) allocate ( BOUNDS(ng) % Jendp1(-1:Ntiles) ) allocate ( BOUNDS(ng) % Jendp2(-1:Ntiles) ) allocate ( BOUNDS(ng) % Jendp2i(-1:Ntiles) ) allocate ( BOUNDS(ng) % Jendp3(-1:Ntiles) ) allocate ( BOUNDS(ng) % Imin(4,0:1,0:Ntiles) ) allocate ( BOUNDS(ng) % Imax(4,0:1,0:Ntiles) ) allocate ( BOUNDS(ng) % Jmin(4,0:1,0:Ntiles) ) allocate ( BOUNDS(ng) % Jmax(4,0:1,0:Ntiles) ) END DO END IF ! ! DOMAIN structure containing boundary edges switches and fractional ! grid lower/upper bounds for each tile. ! IF (.not.allocated(DOMAIN)) THEN allocate ( DOMAIN(Ngrids) ) DO ng=1,Ngrids Ntiles=NtileI(ng)*NtileJ(ng)-1 allocate ( DOMAIN(ng) % Eastern_Edge(-1:Ntiles) ) allocate ( DOMAIN(ng) % Western_Edge(-1:Ntiles) ) allocate ( DOMAIN(ng) % Northern_Edge(-1:Ntiles) ) allocate ( DOMAIN(ng) % Southern_Edge(-1:Ntiles) ) allocate ( DOMAIN(ng) % NorthEast_Corner(-1:Ntiles) ) allocate ( DOMAIN(ng) % NorthWest_Corner(-1:Ntiles) ) allocate ( DOMAIN(ng) % SouthEast_Corner(-1:Ntiles) ) allocate ( DOMAIN(ng) % SouthWest_Corner(-1:Ntiles) ) allocate ( DOMAIN(ng) % NorthEast_Test(-1:Ntiles) ) allocate ( DOMAIN(ng) % NorthWest_Test(-1:Ntiles) ) allocate ( DOMAIN(ng) % SouthEast_Test(-1:Ntiles) ) allocate ( DOMAIN(ng) % SouthWest_Test(-1:Ntiles) ) allocate ( DOMAIN(ng) % Xmin_psi(0:Ntiles) ) allocate ( DOMAIN(ng) % Xmax_psi(0:Ntiles) ) allocate ( DOMAIN(ng) % Ymin_psi(0:Ntiles) ) allocate ( DOMAIN(ng) % Ymax_psi(0:Ntiles) ) allocate ( DOMAIN(ng) % Xmin_rho(0:Ntiles) ) allocate ( DOMAIN(ng) % Xmax_rho(0:Ntiles) ) allocate ( DOMAIN(ng) % Ymin_rho(0:Ntiles) ) allocate ( DOMAIN(ng) % Ymax_rho(0:Ntiles) ) allocate ( DOMAIN(ng) % Xmin_u(0:Ntiles) ) allocate ( DOMAIN(ng) % Xmax_u(0:Ntiles) ) allocate ( DOMAIN(ng) % Ymin_u(0:Ntiles) ) allocate ( DOMAIN(ng) % Ymax_u(0:Ntiles) ) allocate ( DOMAIN(ng) % Xmin_v(0:Ntiles) ) allocate ( DOMAIN(ng) % Xmax_v(0:Ntiles) ) allocate ( DOMAIN(ng) % Ymin_v(0:Ntiles) ) allocate ( DOMAIN(ng) % Ymax_v(0:Ntiles) ) END DO END IF ! ! Allocate lower and upper bounds structure for I/O NetCDF files. ! IF (.not.allocated(IOBOUNDS)) THEN allocate ( IOBOUNDS(Ngrids) ) END IF #ifdef DIAGNOSTICS ! !----------------------------------------------------------------------- ! Determine number of diagnostic variables. !----------------------------------------------------------------------- # ifdef DIAGNOSTICS_TS ! ! Tracer diagnostics. ! NDT=6 ! Acceleration, advection, vertical diffusion # if defined TS_DIF2 || defined TS_DIF4 NDT=NDT+3 ! Horizontal (total, X-, Y-) diffusion # if defined MIX_GEO_TS || defined MIX_ISO_TS NDT=NDT+1 ! Horizontal S-diffusion due to rotated tensor # endif # endif # else NDT=0 ! No tracer diagnostics # endif # ifdef DIAGNOSTICS_UV ! ! 2D Momentum diagnostics. ! NDM2d=4 ! Acceleration, 2D P-Gradient, stresses # ifdef UV_ADV NDM2d=NDM2d+3 ! Horizontal total-, X-, and Y-advection # endif # ifdef WEC_VF # ifdef UV_COR NDM2d=NDM2d+1 ! Coriolis # endif # ifdef BOTTOM_STREAMING NDM2d=NDM2d+1 ! Bottom streaming # endif # ifdef SURFACE_STREAMING NDM2d=NDM2d+1 ! Surface streaming # endif NDM2d=NDM2d+8 ! zeta, zetaw, zqsp, zbeh (all 4 = prsgrd), ! K- and Horiz- VF, roller, break # endif # if defined VEGETATION && defined VEG_DRAG NDM2d=NDM2d+1 ! Resistance imposed on the flow by vegetation # endif # ifdef UV_COR NDM2d=NDM2d+1 ! Coriolis # endif # if defined UV_VIS2 || defined UV_VIS4 NDM2d=NDM2d+3 ! Horizontal total-, X-, and Y-viscosity # endif # ifdef SOLVE3D ! ! 3D Momentum diagnostics and right-hand-side terms. ! NDM3d=3 ! Acceleration, 3D P-Gradient, vertical viscosity NDrhs=1 ! 3D P-Gradient # ifdef UV_ADV NDM3d=NDM3d+4 ! Horizontal (total, X, Y) and vertical advection NDrhs=NDrhs+4 # endif # ifdef WEC_VF # ifdef UV_COR NDM3d=NDM3d+1 ! Coriolis NDrhs=NDrhs+1 # endif # ifdef BOTTOM_STREAMING NDM3d=NDM3d+1 ! Bottom streaming NDrhs=NDrhs+1 # endif # ifdef SURFACE_STREAMING NDM3d=NDM3d+1 ! Surface streaming NDrhs=NDrhs+1 # endif NDM3d=NDM3d+5 ! wave break, roller, K- Horiz and Jvert- VF NDrhs=NDrhs+5 # endif # ifdef UV_COR NDM3d=NDM3d+1 ! Coriolis NDrhs=NDrhs+1 # endif # if defined UV_VIS2 || defined UV_VIS4 NDM3d=NDM3d+3 ! Horizontal (total, X, Y) viscosity # endif # ifdef BODYFORCE NDrhs=NDrhs+1 ! Vertical viscosity # endif # else NDM3d=0 ! No 3D momentum diagnostics NDrhs=0 # endif # endif #endif ! !----------------------------------------------------------------------- ! Derived dimension parameters. !----------------------------------------------------------------------- ! DO ng=1,Ngrids I_padd=(Lm(ng)+2)/2-(Lm(ng)+1)/2 J_padd=(Mm(ng)+2)/2-(Mm(ng)+1)/2 Im(ng)=Lm(ng)+I_padd Jm(ng)=Mm(ng)+J_padd NT(ng)=NAT+NBT+NST+NPT #ifdef FLOATS ! ! Number of floats variables. ! # ifdef FLOAT_OYSTER # ifdef FLOAT_VWALK NFV(ng)=NT(ng)+18 # else NFV(ng)=NT(ng)+16 # endif # else # ifdef FLOAT_VWALK NFV(ng)=NT(ng)+12 # else NFV(ng)=NT(ng)+10 # endif # endif #endif ! ! Number of state variables. The order is irrelevant for determining ! the "NSV" dimension. ! #ifdef SOLVE3D NSV(ng)=5+NT(ng) ! zeta, ubar, vbar, u, v, Tvar(1:MT) # ifdef ADJUST_STFLUX NSV(ng)=NSV(ng)+NT(ng) ! Tsur(1:MT) # endif #else NSV(ng)=3 ! zeta, ubar, vbar #endif #ifdef ADJUST_WSTRESS NSV(ng)=NSV(ng)+2 ! Ustr, Vstr #endif #ifdef WEC NSV(ng)=NSV(ng)+2 ! ubar_stokes, vbar_stokes # if defined SOLVE3D NSV(ng)=NSV(ng)+2 ! u_stokes, v_stokes # endif #endif #ifdef SOLVE3D # if defined GLS_MIXING || defined MY25_MIXING NSV(ng)=NSV(ng)+1 ! TKE # endif NSV(ng)=NSV(ng)+1 ! W #endif END DO ! ! Set the maximum number of tracers between all nested grids. It cannot ! be zero. Usually, NAT=2 (temperature and salinity), but some Test ! Cases do not include salinity. ! MT=MAX(1,MAX(NAT,MAXVAL(NT))) ! !----------------------------------------------------------------------- ! Allocate Lateral Boundary Conditions switches structure. !----------------------------------------------------------------------- ! ! Here, "nLBCvar" is larger than needed because it includes unused ! state variables that need to be included because of their indices ! are in a particular order due to data assimilation variables ! associated with the state vector. ! #ifdef SOLVE3D nLBCvar=5+MT ! zeta, ubar, vbar, u, v, Tvar(1:MT) # ifdef ADJUST_STFLUX nLBCvar=nLBCvar+MT ! Tsur(1:MT); unused # endif #else nLBCvar=3 ! zeta, ubar, vbar #endif #ifdef ADJUST_WSTRESS nLBCvar=nLBCvar+2 ! Ustr, Vstr; unused #endif #ifdef WEC nLBCvar=nLBCvar+2 ! ubar_stokes, vbar_stokes # ifdef SOLVE3D nLBCvar=nLBCvar+2 ! u_stokes, v_stokes # endif #endif #ifdef SOLVE3D # if defined GLS_MIXING || defined MY25_MIXING nLBCvar=nLBCvar+1 ! TKE # endif #endif ! IF (.not.allocated(LBC)) THEN allocate ( LBC(4,nLBCvar,Ngrids) ) DO ng=1,Ngrids DO ivar=1,nLBCvar DO ibry=1,4 LBC(ibry,ivar,ng)%acquire = .FALSE. LBC(ibry,ivar,ng)%Chapman_explicit = .FALSE. LBC(ibry,ivar,ng)%Chapman_implicit = .FALSE. LBC(ibry,ivar,ng)%clamped = .FALSE. LBC(ibry,ivar,ng)%closed = .FALSE. LBC(ibry,ivar,ng)%Flather = .FALSE. LBC(ibry,ivar,ng)%gradient = .FALSE. LBC(ibry,ivar,ng)%mixed = .FALSE. LBC(ibry,ivar,ng)%nested = .FALSE. LBC(ibry,ivar,ng)%nudging = .FALSE. LBC(ibry,ivar,ng)%periodic = .FALSE. LBC(ibry,ivar,ng)%radiation = .FALSE. LBC(ibry,ivar,ng)%reduced = .FALSE. LBC(ibry,ivar,ng)%Shchepetkin = .FALSE. END DO END DO END DO END IF #if defined ADJOINT || defined TANGENT || defined TL_IOMS IF (.not.allocated(ad_LBC)) THEN allocate ( ad_LBC(4,nLBCvar,Ngrids) ) DO ng=1,Ngrids DO ivar=1,nLBCvar DO ibry=1,4 ad_LBC(ibry,ivar,ng)%acquire = .FALSE. ad_LBC(ibry,ivar,ng)%Chapman_implicit = .FALSE. ad_LBC(ibry,ivar,ng)%Chapman_explicit = .FALSE. ad_LBC(ibry,ivar,ng)%clamped = .FALSE. ad_LBC(ibry,ivar,ng)%closed = .FALSE. ad_LBC(ibry,ivar,ng)%Flather = .FALSE. ad_LBC(ibry,ivar,ng)%gradient = .FALSE. ad_LBC(ibry,ivar,ng)%mixed = .FALSE. ad_LBC(ibry,ivar,ng)%nested = .FALSE. ad_LBC(ibry,ivar,ng)%nudging = .FALSE. ad_LBC(ibry,ivar,ng)%periodic = .FALSE. ad_LBC(ibry,ivar,ng)%radiation = .FALSE. ad_LBC(ibry,ivar,ng)%reduced = .FALSE. ad_LBC(ibry,ivar,ng)%Shchepetkin = .FALSE. END DO END DO END DO END IF #endif #if defined TANGENT || defined TL_IOMS IF (.not.allocated(tl_LBC)) THEN allocate ( tl_LBC(4,nLBCvar,Ngrids) ) DO ng=1,Ngrids DO ivar=1,nLBCvar DO ibry=1,4 tl_LBC(ibry,ivar,ng)%acquire = .FALSE. tl_LBC(ibry,ivar,ng)%Chapman_implicit = .FALSE. tl_LBC(ibry,ivar,ng)%Chapman_explicit = .FALSE. tl_LBC(ibry,ivar,ng)%clamped = .FALSE. tl_LBC(ibry,ivar,ng)%closed = .FALSE. tl_LBC(ibry,ivar,ng)%Flather = .FALSE. tl_LBC(ibry,ivar,ng)%gradient = .FALSE. tl_LBC(ibry,ivar,ng)%mixed = .FALSE. tl_LBC(ibry,ivar,ng)%nested = .FALSE. tl_LBC(ibry,ivar,ng)%nudging = .FALSE. tl_LBC(ibry,ivar,ng)%periodic =.FALSE. tl_LBC(ibry,ivar,ng)%radiation = .FALSE. tl_LBC(ibry,ivar,ng)%reduced = .FALSE. tl_LBC(ibry,ivar,ng)%Shchepetkin = .FALSE. END DO END DO END DO END IF #endif #ifdef SOLVE3D ! !----------------------------------------------------------------------- ! Allocate rracer horizontal and vertical advection switches structure. !----------------------------------------------------------------------- ! ! Nonlinear model kernel. ! IF (.not.allocated(Hadvection)) THEN allocate ( Hadvection(MAXVAL(NT),Ngrids) ) DO ng=1,Ngrids DO itrc=1,NT(ng) Hadvection(itrc,ng)%AKIMA4 = .FALSE. Hadvection(itrc,ng)%CENTERED2 = .FALSE. Hadvection(itrc,ng)%CENTERED4 = .FALSE. Hadvection(itrc,ng)%HSIMT = .FALSE. Hadvection(itrc,ng)%MPDATA = .FALSE. Hadvection(itrc,ng)%SPLINES = .FALSE. Hadvection(itrc,ng)%SPLIT_U3 = .FALSE. Hadvection(itrc,ng)%UPSTREAM3 = .FALSE. END DO END DO END IF ! IF (.not.allocated(Vadvection)) THEN allocate ( Vadvection(MAXVAL(NT),Ngrids) ) DO ng=1,Ngrids DO itrc=1,NT(ng) Vadvection(itrc,ng)%AKIMA4 = .FALSE. Vadvection(itrc,ng)%CENTERED2 = .FALSE. Vadvection(itrc,ng)%CENTERED4 = .FALSE. Vadvection(itrc,ng)%HSIMT = .FALSE. Vadvection(itrc,ng)%MPDATA = .FALSE. Vadvection(itrc,ng)%SPLINES = .FALSE. Vadvection(itrc,ng)%SPLIT_U3 = .FALSE. Vadvection(itrc,ng)%UPSTREAM3 = .FALSE. END DO END DO END IF # if defined ADJOINT || defined TANGENT || defined TL_IOMS ! ! Adjoint model kernel. ! IF (.not.allocated(ad_Hadvection)) THEN allocate ( ad_Hadvection(MAXVAL(NT),Ngrids) ) DO ng=1,Ngrids DO itrc=1,NT(ng) ad_Hadvection(itrc,ng)%AKIMA4 = .FALSE. ad_Hadvection(itrc,ng)%CENTERED2 = .FALSE. ad_Hadvection(itrc,ng)%CENTERED4 = .FALSE. ad_Hadvection(itrc,ng)%HSIMT = .FALSE. ad_Hadvection(itrc,ng)%MPDATA = .FALSE. ad_Hadvection(itrc,ng)%SPLINES = .FALSE. ad_Hadvection(itrc,ng)%SPLIT_U3 = .FALSE. ad_Hadvection(itrc,ng)%UPSTREAM3 = .FALSE. END DO END DO END IF ! IF (.not.allocated(ad_Vadvection)) THEN allocate ( ad_Vadvection(MAXVAL(NT),Ngrids) ) DO ng=1,Ngrids DO itrc=1,NT(ng) ad_Vadvection(itrc,ng)%AKIMA4 = .FALSE. ad_Vadvection(itrc,ng)%CENTERED2 = .FALSE. ad_Vadvection(itrc,ng)%CENTERED4 = .FALSE. ad_Vadvection(itrc,ng)%HSIMT = .FALSE. ad_Vadvection(itrc,ng)%MPDATA = .FALSE. ad_Vadvection(itrc,ng)%SPLINES = .FALSE. ad_Vadvection(itrc,ng)%SPLIT_U3 = .FALSE. ad_Vadvection(itrc,ng)%UPSTREAM3 = .FALSE. END DO END DO END IF # endif # if defined TANGENT || defined TL_IOMS ! ! Perturbation or finite-amplitude tangent linear model kernels. ! IF (.not.allocated(tl_Hadvection)) THEN allocate ( tl_Hadvection(MAXVAL(NT),Ngrids) ) DO ng=1,Ngrids DO itrc=1,NT(ng) tl_Hadvection(itrc,ng)%AKIMA4 = .FALSE. tl_Hadvection(itrc,ng)%CENTERED2 = .FALSE. tl_Hadvection(itrc,ng)%CENTERED4 = .FALSE. tl_Hadvection(itrc,ng)%HSIMT = .FALSE. tl_Hadvection(itrc,ng)%MPDATA = .FALSE. tl_Hadvection(itrc,ng)%SPLINES = .FALSE. tl_Hadvection(itrc,ng)%SPLIT_U3 = .FALSE. tl_Hadvection(itrc,ng)%UPSTREAM3 = .FALSE. END DO END DO END IF ! IF (.not.allocated(tl_Vadvection)) THEN allocate ( tl_Vadvection(MAXVAL(NT),Ngrids) ) DO ng=1,Ngrids DO itrc=1,NT(ng) tl_Vadvection(itrc,ng)%AKIMA4 = .FALSE. tl_Vadvection(itrc,ng)%CENTERED2 = .FALSE. tl_Vadvection(itrc,ng)%CENTERED4 = .FALSE. tl_Vadvection(itrc,ng)%HSIMT = .FALSE. tl_Vadvection(itrc,ng)%MPDATA = .FALSE. tl_Vadvection(itrc,ng)%SPLINES = .FALSE. tl_Vadvection(itrc,ng)%SPLIT_U3 = .FALSE. tl_Vadvection(itrc,ng)%UPSTREAM3 = .FALSE. END DO END DO END IF # endif #endif ! RETURN END SUBROUTINE initialize_param END MODULE mod_param