MODULE mod_forces ! !git $Id$ !svn $Id: mod_forces.F 1178 2023-07-11 17:50:57Z arango $ !================================================== Hernan G. Arango === ! Copyright (c) 2002-2023 The ROMS/TOMS Group ! ! Licensed under a MIT/X style license ! ! See License_ROMS.md ! !======================================================================= ! ! ! Surface momentum stresses. ! ! ! ! sustr Surface momentum flux (wind stress) in the ! ! XI-direction (m2/s2) at horizontal U-points. ! ! sustrG Latest two-time snapshots of input "sustr" grided ! ! data used for interpolation. ! ! svstr Surface momentum flux (wind stress) in the ! ! ETA-direction (m2/s2) at horizontal V-points. ! ! svstrG Latest two-time snapshots of input "svstr" grided ! ! data used for interpolation. ! ! ! ! Bottom momentum stresses. ! ! ! ! bustr Bottom momentum flux (bottom stress) in the ! ! XI-direction (m2/s2) at horizontal U-points. ! ! bvstr Bottom momentum flux (bottom stress) in the ! ! ETA-direction (m2/s2) at horizontal V-points. ! ! ! ! Surface wind induced waves. ! ! ! ! Hwave Surface wind induced wave height (m). ! ! HwaveG Latest two-time snapshots of input "Hwave" grided ! ! data used for interpolation. ! ! Dwave Surface wind induced mean wave direction (radians). ! ! DwaveG Latest two-time snapshots of input "Dwave" grided ! ! data used for interpolation. ! ! Dwavep Surface wind induced peak wave direction (radians). ! ! DwavepG Latest two-time snapshots of input "Dwavep" grided ! ! data used for interpolation. ! ! Lwave Mean surface wavelength read in from swan output ! ! LwaveG Latest two-time snapshots of input "Lwave" grided ! ! data used for interpolation. ! ! Lwavep Peak surface wavelength read in from swan output ! ! LwavepG Latest two-time snapshots of input "Lwavep" grided ! ! data used for interpolation. ! ! Pwave_top Wind induced surface wave period (s). ! ! Pwave_topG Latest two-time snapshots of input "Pwave_top" grided ! ! data used for interpolation. ! ! Pwave_bot Wind induced bottom wave period (s). ! ! Pwave_botG Latest two-time snapshots of input "Pwave_bot" grided ! ! data used for interpolation. ! ! Uwave_rms Bottom orbital velocity read in from swan output ! ! Uwave_rmsG Latest two-time snapshots of input "Uwave_rms" grided ! ! data used for interpolation. ! ! Wave_break Percent of wave breaking for use with roller model. ! ! Wave_breakG Latest two-time snapshots of input "wave_break" ! ! gridded data used for interpolation. ! ! Wave_ds Wave directional spreading. ! ! Wave_qp Wave spectrum peakedness. ! ! ! ! Solar shortwave radiation flux. ! ! ! ! srflx Surface shortwave solar radiation flux (degC m/s) ! ! at horizontal RHO-points ! ! srflxG Latest two-time snapshots of input "srflx" grided ! ! data used for interpolation. ! ! ! ! Cloud fraction. ! ! ! ! cloud Cloud fraction (percentage/100). ! ! cloudG Latest two-time snapshots of input "cloud" grided ! ! data used for interpolation. ! ! ! ! Surface heat fluxes, Atmosphere-Ocean bulk parameterization. ! ! ! ! lhflx Latent heat flux (degC m/s). ! ! lrflx Longwave radiation (degC m/s). ! ! shflx Sensible heat flux (degC m/s). ! ! ! ! Surface air humidity. ! ! ! ! Hair Surface air specific (g/kg) or relative humidity ! ! (percentage). ! ! HairG Latest two-time snapshots of input "Hair" grided ! ! data used for interpolation. ! ! ! ! Surface air pressure. ! ! ! ! Pair Surface air pressure (mb). ! ! PairG Latest two-time snapshots of input "Pair" grided ! ! data used for interpolation. ! ! ! ! Surface air temperature. ! ! ! ! Tair Surface air temperature (Celsius) ! ! TairG Latest two-time snapshots of input "Tair" grided ! ! data used for interpolation. ! ! Surface Winds. ! ! ! ! Uwind Surface wind in the XI-direction (m/s) at ! ! horizontal RHO-points. ! ! UwindG Latest two-time snapshots of input "Uwind" grided ! ! data used for interpolation. ! ! Vwind Surface wind in the ETA-direction (m/s) at ! ! horizontal RHO-points. ! ! VwindG Latest two-time snapshots of input "Vwind" grided ! ! data used for interpolation. ! ! ! ! Rain fall rate. ! ! ! ! evap Evaporation rate (kg/m2/s). ! ! rain Rain fall rate (kg/m2/s). ! ! rainG Latest two-time snapshots of input "rain" grided ! ! data used for interpolation. ! ! ! ! Surface tracer fluxes. ! ! ! ! stflux Forcing surface flux of tracer type variables from ! ! data, coupling, bulk flux parameterization, or ! ! analytical formulas. ! ! ! ! stflux(:,:,itemp) surface net heat flux ! ! stflux(:,:,isalt) surface net freshwater flux (E-P) ! ! ! ! stfluxG Latest two-time snapshots of input "stflux" grided ! ! data used for interpolation. ! ! ! ! stflx ROMS state surface flux of tracer type variables ! ! (TracerUnits m/s) at horizontal RHO-points, as used ! ! in the governing equations. ! ! ! ! Bottom tracer fluxes. ! ! ! ! btflux Forcing bottom flux of tracer type variables from ! ! data or analytical formulas. Usually, the bottom ! ! flux of tracer is zero. ! ! ! ! btflux(:,:,itemp) bottom heat flux ! ! btflux(:,:,isalt) bottom freshwater flux ! ! ! ! btfluxG Latest two-time snapshots of input "vtflux" grided ! ! data used for interpolation. ! ! ! ! btflx ROMS state bottom flux of tracer type variables ! ! (TracerUnits m/s) at horizontal RHO-points, as used ! ! in the governing equations. ! ! ! ! Surface heat flux correction. ! ! ! ! dqdt Surface net heat flux sensitivity to SST, ! ! d(Q)/d(SST), (m/s). ! ! dqdtG Latest two-time snapshots of input "dqdt" grided ! ! data used for interpolation. ! ! sst Sea surface temperature (Celsius). ! ! sstG Latest two-time snapshots of input "sst" grided ! ! data used for interpolation. ! ! ! ! Surface freshwater flux correction. ! ! ! ! sss Sea surface salinity (PSU). ! ! sssG Latest two-time snapshots of input "sss" grided ! ! data used for interpolation. ! ! ! ! Surface spectral downwelling irradiance. ! ! ! ! SpecIr Spectral irradiance (NBands) from 400-700 nm at ! ! 5 nm bandwidth. ! ! avcos Cosine of average zenith angle of downwelling ! ! spectral photons. ! ! ! !======================================================================= ! USE mod_kinds ! implicit none ! PUBLIC :: allocate_forces PUBLIC :: deallocate_forces PUBLIC :: initialize_forces ! !----------------------------------------------------------------------- ! Define T_FORCES structure. !----------------------------------------------------------------------- ! TYPE T_FORCES ! ! Nonlinear model state. ! real(r8), pointer :: sustr(:,:) real(r8), pointer :: svstr(:,:) real(r8), pointer :: sustrG(:,:,:) real(r8), pointer :: svstrG(:,:,:) real(r8), pointer :: bustr(:,:) real(r8), pointer :: bvstr(:,:) real(r8), pointer :: Pair(:,:) real(r8), pointer :: PairG(:,:,:) real(r8), pointer :: srflx(:,:) real(r8), pointer :: srflxG(:,:,:) real(r8), pointer :: lhflx(:,:) real(r8), pointer :: lrflx(:,:) real(r8), pointer :: shflx(:,:) real(r8), pointer :: lrflxG(:,:,:) real(r8), pointer :: Hair(:,:) real(r8), pointer :: HairG(:,:,:) real(r8), pointer :: Tair(:,:) real(r8), pointer :: TairG(:,:,:) real(r8), pointer :: Uwind(:,:) real(r8), pointer :: Vwind(:,:) real(r8), pointer :: UwindG(:,:,:) real(r8), pointer :: VwindG(:,:,:) real(r8), pointer :: rain(:,:) real(r8), pointer :: rainG(:,:,:) real(r8), pointer :: evap(:,:) real(r8), pointer :: stflux(:,:,:) real(r8), pointer :: stfluxG(:,:,:,:) real(r8), pointer :: stflx(:,:,:) real(r8), pointer :: btflux(:,:,:) real(r8), pointer :: btfluxG(:,:,:,:) real(r8), pointer :: btflx(:,:,:) ! ! Tangent linear model state. ! real(r8), pointer :: tl_sustr(:,:) real(r8), pointer :: tl_svstr(:,:) real(r8), pointer :: tl_bustr(:,:) real(r8), pointer :: tl_bvstr(:,:) real(r8), pointer :: tl_stflux(:,:,:) real(r8), pointer :: tl_stflx (:,:,:) real(r8), pointer :: tl_btflux(:,:,:) real(r8), pointer :: tl_btflx (:,:,:) real(r8), pointer :: tl_srflx(:,:) real(r8), pointer :: tl_lhflx(:,:) real(r8), pointer :: tl_lrflx(:,:) real(r8), pointer :: tl_shflx(:,:) real(r8), pointer :: tl_evap(:,:) ! ! Adjoint model state. ! real(r8), pointer :: ad_sustr(:,:) real(r8), pointer :: ad_svstr(:,:) real(r8), pointer :: ad_bustr(:,:) real(r8), pointer :: ad_bvstr(:,:) real(r8), pointer :: ad_bustr_sol(:,:) real(r8), pointer :: ad_bvstr_sol(:,:) real(r8), pointer :: ad_stflx(:,:,:) real(r8), pointer :: ad_btflx(:,:,:) real(r8), pointer :: ad_srflx(:,:) real(r8), pointer :: ad_lhflx(:,:) real(r8), pointer :: ad_lrflx(:,:) real(r8), pointer :: ad_shflx(:,:) real(r8), pointer :: ad_evap(:,:) END TYPE T_FORCES ! TYPE (T_FORCES), allocatable :: FORCES(:) ! CONTAINS ! SUBROUTINE allocate_forces (ng, LBi, UBi, LBj, UBj) ! !======================================================================= ! ! ! This routine allocates all variables in the module for all nested ! ! grids. ! ! ! !======================================================================= ! USE mod_param ! ! Local variable declarations. ! integer, intent(in) :: ng, LBi, UBi, LBj, UBj ! ! Local variable declarations. ! real(r8) :: size2d ! !----------------------------------------------------------------------- ! Allocate module variables. !----------------------------------------------------------------------- ! IF (ng.eq.1) allocate ( FORCES(Ngrids) ) ! ! Set horizontal array size. ! size2d=REAL((UBi-LBi+1)*(UBj-LBj+1),r8) ! ! Nonlinear model state ! allocate ( FORCES(ng) % sustr(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % svstr(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % sustrG(LBi:UBi,LBj:UBj,2) ) Dmem(ng)=Dmem(ng)+2.0_r8*size2d allocate ( FORCES(ng) % svstrG(LBi:UBi,LBj:UBj,2) ) Dmem(ng)=Dmem(ng)+2.0_r8*size2d allocate ( FORCES(ng) % bustr(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % bvstr(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % Pair(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % PairG(LBi:UBi,LBj:UBj,2) ) Dmem(ng)=Dmem(ng)+2.0_r8*size2d allocate ( FORCES(ng) % srflx(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % srflxG(LBi:UBi,LBj:UBj,2) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % lhflx(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % lrflx(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % shflx(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % lrflxG(LBi:UBi,LBj:UBj,2) ) Dmem(ng)=Dmem(ng)+2.0_r8*size2d allocate ( FORCES(ng) % Hair(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % HairG(LBi:UBi,LBj:UBj,2) ) Dmem(ng)=Dmem(ng)+2.0_r8*size2d allocate ( FORCES(ng) % Tair(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % TairG(LBi:UBi,LBj:UBj,2) ) Dmem(ng)=Dmem(ng)+2.0_r8*size2d allocate ( FORCES(ng) % Uwind(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % Vwind(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % UwindG(LBi:UBi,LBj:UBj,2) ) Dmem(ng)=Dmem(ng)+2.0_r8*size2d allocate ( FORCES(ng) % VwindG(LBi:UBi,LBj:UBj,2) ) Dmem(ng)=Dmem(ng)+2.0_r8*size2d allocate ( FORCES(ng) % rain(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % rainG(LBi:UBi,LBj:UBj,2) ) Dmem(ng)=Dmem(ng)+2.0_r8*size2d allocate ( FORCES(ng) % evap(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % stflux(LBi:UBi,LBj:UBj,NT(ng)) ) Dmem(ng)=Dmem(ng)+REAL(NT(ng),r8)*size2d allocate ( FORCES(ng) % stfluxG(LBi:UBi,LBj:UBj,2,NT(ng)) ) Dmem(ng)=Dmem(ng)+2.0_r8*REAL(NT(ng),r8)*size2d allocate ( FORCES(ng) % stflx(LBi:UBi,LBj:UBj,NT(ng)) ) Dmem(ng)=Dmem(ng)+REAL(NT(ng),r8)*size2d allocate ( FORCES(ng) % btflux(LBi:UBi,LBj:UBj,NT(ng)) ) Dmem(ng)=Dmem(ng)+REAL(NT(ng),r8)*size2d allocate ( FORCES(ng) % btfluxG(LBi:UBi,LBj:UBj,2,NT(ng)) ) Dmem(ng)=Dmem(ng)+2.0_r8*REAL(NT(ng),r8)*size2d allocate ( FORCES(ng) % btflx(LBi:UBi,LBj:UBj,NT(ng)) ) Dmem(ng)=Dmem(ng)+REAL(NT(ng),r8)*size2d ! ! Tangent linear model state ! allocate ( FORCES(ng) % tl_sustr(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % tl_svstr(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % tl_bustr(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % tl_bvstr(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % tl_stflux(LBi:UBi,LBj:UBj,NT(ng)) ) Dmem(ng)=Dmem(ng)+REAL(NT(ng),r8)*size2d allocate ( FORCES(ng) % tl_stflx(LBi:UBi,LBj:UBj,NT(ng)) ) Dmem(ng)=Dmem(ng)+REAL(NT(ng),r8)*size2d allocate ( FORCES(ng) % tl_btflux(LBi:UBi,LBj:UBj,NT(ng)) ) Dmem(ng)=Dmem(ng)+REAL(NT(ng),r8)*size2d allocate ( FORCES(ng) % tl_btflx(LBi:UBi,LBj:UBj,NT(ng)) ) Dmem(ng)=Dmem(ng)+REAL(NT(ng),r8)*size2d allocate ( FORCES(ng) % tl_srflx(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % tl_lhflx(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % tl_lrflx(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % tl_shflx(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % tl_evap(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d ! ! Adjoint model state ! allocate ( FORCES(ng) % ad_sustr(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % ad_svstr(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % ad_bustr(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % ad_bvstr(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % ad_bustr_sol(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % ad_bvstr_sol(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % ad_stflx(LBi:UBi,LBj:UBj,NT(ng)) ) Dmem(ng)=Dmem(ng)+REAL(NT(ng),r8)*size2d allocate ( FORCES(ng) % ad_btflx(LBi:UBi,LBj:UBj,NT(ng)) ) Dmem(ng)=Dmem(ng)+REAL(NT(ng),r8)*size2d allocate ( FORCES(ng) % ad_srflx(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % ad_lhflx(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % ad_lrflx(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % ad_shflx(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d allocate ( FORCES(ng) % ad_evap(LBi:UBi,LBj:UBj) ) Dmem(ng)=Dmem(ng)+size2d ! RETURN END SUBROUTINE allocate_forces ! SUBROUTINE deallocate_forces (ng) ! !======================================================================= ! ! ! This routine deallocates all variables in the module for all nested ! ! grids. ! ! ! !======================================================================= ! USE mod_param, ONLY : Ngrids ! ! Imported variable declarations. ! integer, intent(in) :: ng ! ! Local variable declarations. ! character (len=*), parameter :: MyFile = & & "ROMS/Modules/mod_forces.F"//", deallocate_forces" ! !----------------------------------------------------------------------- ! Deallocate derived-type FORCES structure. !----------------------------------------------------------------------- ! IF (ng.eq.Ngrids) THEN IF (allocated(FORCES)) deallocate ( FORCES ) END IF ! RETURN END SUBROUTINE deallocate_forces ! SUBROUTINE initialize_forces (ng, tile, model) ! !======================================================================= ! ! ! This routine initialize all variables in the module using first ! ! touch distribution policy. In shared-memory configuration, this ! ! operation actually performs propagation of the "shared arrays" ! ! across the cluster, unless another policy is specified to ! ! override the default. ! ! ! !======================================================================= ! USE mod_param ! ! Imported variable declarations. ! integer, intent(in) :: ng, tile, model ! ! Local variable declarations. ! integer :: Imin, Imax, Jmin, Jmax integer :: i, j, k integer :: itrc real(r8), parameter :: IniVal = 0.0_r8 ! !----------------------------------------------------------------------- ! Set lower and upper tile bounds and staggered variables bounds for ! this horizontal domain partition. Notice that if tile=-1, it will ! set the values for the global grid. !----------------------------------------------------------------------- ! integer :: Istr, IstrB, IstrP, IstrR, IstrT, IstrM, IstrU integer :: Iend, IendB, IendP, IendR, IendT integer :: Jstr, JstrB, JstrP, JstrR, JstrT, JstrM, JstrV integer :: Jend, JendB, JendP, JendR, JendT integer :: Istrm3, Istrm2, Istrm1, IstrUm2, IstrUm1 integer :: Iendp1, Iendp2, Iendp2i, Iendp3 integer :: Jstrm3, Jstrm2, Jstrm1, JstrVm2, JstrVm1 integer :: Jendp1, Jendp2, Jendp2i, Jendp3 ! Istr =BOUNDS(ng) % Istr (tile) IstrB =BOUNDS(ng) % IstrB (tile) IstrM =BOUNDS(ng) % IstrM (tile) IstrP =BOUNDS(ng) % IstrP (tile) IstrR =BOUNDS(ng) % IstrR (tile) IstrT =BOUNDS(ng) % IstrT (tile) IstrU =BOUNDS(ng) % IstrU (tile) Iend =BOUNDS(ng) % Iend (tile) IendB =BOUNDS(ng) % IendB (tile) IendP =BOUNDS(ng) % IendP (tile) IendR =BOUNDS(ng) % IendR (tile) IendT =BOUNDS(ng) % IendT (tile) Jstr =BOUNDS(ng) % Jstr (tile) JstrB =BOUNDS(ng) % JstrB (tile) JstrM =BOUNDS(ng) % JstrM (tile) JstrP =BOUNDS(ng) % JstrP (tile) JstrR =BOUNDS(ng) % JstrR (tile) JstrT =BOUNDS(ng) % JstrT (tile) JstrV =BOUNDS(ng) % JstrV (tile) Jend =BOUNDS(ng) % Jend (tile) JendB =BOUNDS(ng) % JendB (tile) JendP =BOUNDS(ng) % JendP (tile) JendR =BOUNDS(ng) % JendR (tile) JendT =BOUNDS(ng) % JendT (tile) ! Istrm3 =BOUNDS(ng) % Istrm3 (tile) ! Istr-3 Istrm2 =BOUNDS(ng) % Istrm2 (tile) ! Istr-2 Istrm1 =BOUNDS(ng) % Istrm1 (tile) ! Istr-1 IstrUm2=BOUNDS(ng) % IstrUm2(tile) ! IstrU-2 IstrUm1=BOUNDS(ng) % IstrUm1(tile) ! IstrU-1 Iendp1 =BOUNDS(ng) % Iendp1 (tile) ! Iend+1 Iendp2 =BOUNDS(ng) % Iendp2 (tile) ! Iend+2 Iendp2i=BOUNDS(ng) % Iendp2i(tile) ! Iend+2 interior Iendp3 =BOUNDS(ng) % Iendp3 (tile) ! Iend+3 Jstrm3 =BOUNDS(ng) % Jstrm3 (tile) ! Jstr-3 Jstrm2 =BOUNDS(ng) % Jstrm2 (tile) ! Jstr-2 Jstrm1 =BOUNDS(ng) % Jstrm1 (tile) ! Jstr-1 JstrVm2=BOUNDS(ng) % JstrVm2(tile) ! JstrV-2 JstrVm1=BOUNDS(ng) % JstrVm1(tile) ! JstrV-1 Jendp1 =BOUNDS(ng) % Jendp1 (tile) ! Jend+1 Jendp2 =BOUNDS(ng) % Jendp2 (tile) ! Jend+2 Jendp2i=BOUNDS(ng) % Jendp2i(tile) ! Jend+2 interior Jendp3 =BOUNDS(ng) % Jendp3 (tile) ! Jend+3 ! ! Set array initialization range. ! Imin=BOUNDS(ng)%LBi(tile) Imax=BOUNDS(ng)%UBi(tile) Jmin=BOUNDS(ng)%LBj(tile) Jmax=BOUNDS(ng)%UBj(tile) ! !----------------------------------------------------------------------- ! Initialize module variables. !----------------------------------------------------------------------- ! ! Nonlinear model state. ! IF ((model.eq.0).or.(model.eq.iNLM)) THEN DO j=Jmin,Jmax DO i=Imin,Imax FORCES(ng) % sustr(i,j) = IniVal FORCES(ng) % svstr(i,j) = IniVal FORCES(ng) % sustrG(i,j,1) = IniVal FORCES(ng) % sustrG(i,j,2) = IniVal FORCES(ng) % svstrG(i,j,1) = IniVal FORCES(ng) % svstrG(i,j,2) = IniVal FORCES(ng) % bustr(i,j) = IniVal FORCES(ng) % bvstr(i,j) = IniVal FORCES(ng) % Pair(i,j) = IniVal FORCES(ng) % PairG(i,j,1) = IniVal FORCES(ng) % PairG(i,j,2) = IniVal FORCES(ng) % srflx(i,j) = IniVal FORCES(ng) % srflxG(i,j,1) = IniVal FORCES(ng) % srflxG(i,j,2) = IniVal FORCES(ng) % lhflx(i,j) = IniVal FORCES(ng) % lrflx(i,j) = IniVal FORCES(ng) % shflx(i,j) = IniVal FORCES(ng) % Hair(i,j) = IniVal FORCES(ng) % Tair(i,j) = IniVal FORCES(ng) % Uwind(i,j) = IniVal FORCES(ng) % Vwind(i,j) = IniVal FORCES(ng) % rain(i,j) = IniVal FORCES(ng) % evap(i,j) = IniVal FORCES(ng) % lrflxG(i,j,1) = IniVal FORCES(ng) % lrflxG(i,j,2) = IniVal FORCES(ng) % HairG(i,j,1) = IniVal FORCES(ng) % HairG(i,j,2) = IniVal FORCES(ng) % TairG(i,j,1) = IniVal FORCES(ng) % TairG(i,j,2) = IniVal FORCES(ng) % UwindG(i,j,1) = IniVal FORCES(ng) % UwindG(i,j,2) = IniVal FORCES(ng) % VwindG(i,j,1) = IniVal FORCES(ng) % VwindG(i,j,2) = IniVal FORCES(ng) % rainG(i,j,1) = IniVal FORCES(ng) % rainG(i,j,2) = IniVal DO itrc=1,NT(ng) FORCES(ng) % stflux(i,j,itrc) = IniVal FORCES(ng) % stfluxG(i,j,1,itrc) = IniVal FORCES(ng) % stfluxG(i,j,2,itrc) = IniVal FORCES(ng) % stflx(i,j,itrc) = IniVal FORCES(ng) % btflux(i,j,itrc) = IniVal FORCES(ng) % btfluxG(i,j,1,itrc) = IniVal FORCES(ng) % btfluxG(i,j,2,itrc) = IniVal FORCES(ng) % btflx(i,j,itrc) = IniVal END DO END DO END DO END IF ! ! Tangent linear model state. ! IF ((model.eq.0).or.(model.eq.iTLM).or.(model.eq.iRPM)) THEN DO j=Jmin,Jmax DO i=Imin,Imax FORCES(ng) % tl_sustr(i,j) = IniVal FORCES(ng) % tl_svstr(i,j) = IniVal FORCES(ng) % tl_bustr(i,j) = IniVal FORCES(ng) % tl_bvstr(i,j) = IniVal END DO DO i=Imin,Imax FORCES(ng) % tl_srflx(i,j) = IniVal END DO DO i=Imin,Imax FORCES(ng) % tl_lhflx(i,j) = IniVal FORCES(ng) % tl_lrflx(i,j) = IniVal FORCES(ng) % tl_shflx(i,j) = IniVal FORCES(ng) % tl_evap(i,j) = IniVal END DO DO itrc=1,NT(ng) DO i=Imin,Imax FORCES(ng) % tl_stflux(i,j,itrc) = IniVal FORCES(ng) % tl_stflx (i,j,itrc) = IniVal FORCES(ng) % tl_btflux(i,j,itrc) = IniVal FORCES(ng) % tl_btflx (i,j,itrc) = IniVal END DO END DO END DO END IF ! ! Adjoint model state. ! IF ((model.eq.0).or.(model.eq.iADM)) THEN DO j=Jmin,Jmax DO i=Imin,Imax FORCES(ng) % ad_sustr(i,j) = IniVal FORCES(ng) % ad_svstr(i,j) = IniVal FORCES(ng) % ad_bustr(i,j) = IniVal FORCES(ng) % ad_bvstr(i,j) = IniVal FORCES(ng) % ad_bustr_sol(i,j) = IniVal FORCES(ng) % ad_bvstr_sol(i,j) = IniVal END DO DO i=Imin,Imax FORCES(ng) % ad_srflx(i,j) = IniVal END DO DO i=Imin,Imax FORCES(ng) % ad_lhflx(i,j) = IniVal FORCES(ng) % ad_lrflx(i,j) = IniVal FORCES(ng) % ad_shflx(i,j) = IniVal FORCES(ng) % ad_evap(i,j) = IniVal END DO DO itrc=1,NT(ng) DO i=Imin,Imax FORCES(ng) % ad_stflx(i,j,itrc) = IniVal FORCES(ng) % ad_btflx(i,j,itrc) = IniVal END DO END DO END DO END IF ! RETURN END SUBROUTINE initialize_forces END MODULE mod_forces