MODULE mod_ocean ! !svn $Id: mod_ocean.F 407 2009-11-02 21:27:07Z arango $ !================================================== Hernan G. Arango === ! Copyright (c) 2002-2009 The ROMS/TOMS Group ! ! Licensed under a MIT/X style license ! ! See License_ROMS.txt ! !======================================================================= ! ! ! 2D Primitive Variables. ! ! ! ! rubar Right-hand-side of 2D U-momentum equation (m4/s2). ! ! rvbar Right-hand-side of 2D V-momentum equation (m4/s2). ! ! rzeta Right-hand-side of free surface equation (m3/s). ! ! ubar Vertically integrated U-momentum component (m/s). ! ! vbar Vertically integrated V-momentum component (m/s). ! ! zeta Free surface (m). ! ! ! ! 3D Primitive Variables. ! ! ! ! pden Potential Density anomaly (kg/m3). ! ! rho Density anomaly (kg/m3). ! ! ru Right-hand-side of 3D U-momentum equation (m4/s2). ! ! rv Right hand side of 3D V-momentum equation (m4/s2). ! ! t Tracer type variables (active and passive). ! ! u 3D U-momentum component (m/s). ! ! v 3D V-momentum component (m/s). ! ! W S-coordinate (omega*Hz/mn) vertical velocity (m3/s). ! ! ! ! Biology Variables. ! ! ! ! pH Surface concentration of hydrogen ions. ! ! ! ! Sediment Variables. ! ! ! ! bed Sediment properties in each bed layer: ! ! bed(:,:,:,ithck) => layer thickness ! ! bed(:,:,:,iaged) => layer age ! ! bed(:,:,:,iporo) => layer porosity ! ! bed(:,:,:,idiff) => layer bio-diffusivity ! ! bedldu Bed load u-transport (kg/m/s) ! ! bedldv Bed load v-transport (kg/m/s) ! ! bed_frac Sediment fraction of each size class in each bed layer ! ! (nondimensional: 0-1.0). Sum of bed_frac = 1.0 ! ! bed_mass Sediment mass of each size class in each bed layer ! ! (kg/m2). ! ! bottom Exposed sediment layer properties: ! ! bottom(:,:,isd50) => mean grain diameter ! ! bottom(:,:,idens) => mean grain density ! ! bottom(:,:,iwsed) => mean settling velocity ! ! bottom(:,:,itauc) => mean critical erosion stress ! ! bottom(:,:,irlen) => ripple length ! ! bottom(:,:,irhgt) => ripple height ! ! bottom(:,:,ibwav) => bed wave excursion amplitude ! ! bottom(:,:,izNik) => Nikuradse bottom roughness ! ! bottom(:,:,izbio) => biological bottom roughness ! ! bottom(:,:,izbfm) => bed form bottom roughness ! ! bottom(:,:,izbld) => bed load bottom roughness ! ! bottom(:,:,izapp) => apparent bottom roughness ! ! bottom(:,:,izwbl) => wave bottom roughness ! ! bottom(:,:,izdef) => default bottom roughness ! ! bottom(:,:,iactv) => active layer thickness ! ! bottom(:,:,ishgt) => saltation height ! ! ero_flux Flux from erosion. ! ! settling_flux Flux from settling. ! ! ! !======================================================================= ! USE mod_kinds implicit none TYPE T_OCEAN ! ! Nonlinear model state. ! real(r8), pointer :: rubar(:,:,:) real(r8), pointer :: rvbar(:,:,:) real(r8), pointer :: rzeta(:,:,:) real(r8), pointer :: ubar(:,:,:) real(r8), pointer :: vbar(:,:,:) real(r8), pointer :: zeta(:,:,:) real(r8), pointer :: pden(:,:,:) real(r8), pointer :: rho(:,:,:) real(r8), pointer :: ru(:,:,:,:) real(r8), pointer :: rv(:,:,:,:) real(r8), pointer :: t(:,:,:,:,:) real(r8), pointer :: u(:,:,:,:) real(r8), pointer :: v(:,:,:,:) real(r8), pointer :: W(:,:,:) real(r8), pointer :: wvel(:,:,:) real(r8), pointer :: bed(:,:,:,:) real(r8), pointer :: bed_frac(:,:,:,:) real(r8), pointer :: bed_mass(:,:,:,:,:) real(r8), pointer :: ero_flux(:,:,:) real(r8), pointer :: settling_flux(:,:,:) real(r8), pointer :: bottom(:,:,:) END TYPE T_OCEAN TYPE (T_OCEAN), allocatable :: OCEAN(:) CONTAINS SUBROUTINE allocate_ocean (ng, LBi, UBi, LBj, UBj) ! !======================================================================= ! ! ! This routine allocates all variables in the module for all nested ! ! grids. ! ! ! !======================================================================= ! USE mod_param USE mod_sediment ! ! Imported variable declarations. ! integer, intent(in) :: ng, LBi, UBi, LBj, UBj ! !----------------------------------------------------------------------- ! Allocate and initialize module variables. !----------------------------------------------------------------------- ! IF (ng.eq.1) allocate ( OCEAN(Ngrids) ) ! ! Nonlinear model state. ! allocate ( OCEAN(ng) % rubar(LBi:UBi,LBj:UBj,2) ) allocate ( OCEAN(ng) % rvbar(LBi:UBi,LBj:UBj,2) ) allocate ( OCEAN(ng) % rzeta(LBi:UBi,LBj:UBj,2) ) allocate ( OCEAN(ng) % ubar(LBi:UBi,LBj:UBj,3) ) allocate ( OCEAN(ng) % vbar(LBi:UBi,LBj:UBj,3) ) allocate ( OCEAN(ng) % zeta(LBi:UBi,LBj:UBj,3) ) allocate ( OCEAN(ng) % pden(LBi:UBi,LBj:UBj,N(ng)) ) allocate ( OCEAN(ng) % rho(LBi:UBi,LBj:UBj,N(ng)) ) allocate ( OCEAN(ng) % ru(LBi:UBi,LBj:UBj,0:N(ng),2) ) allocate ( OCEAN(ng) % rv(LBi:UBi,LBj:UBj,0:N(ng),2) ) allocate ( OCEAN(ng) % t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng)) ) allocate ( OCEAN(ng) % u(LBi:UBi,LBj:UBj,N(ng),2) ) allocate ( OCEAN(ng) % v(LBi:UBi,LBj:UBj,N(ng),2) ) allocate ( OCEAN(ng) % W(LBi:UBi,LBj:UBj,0:N(ng)) ) allocate ( OCEAN(ng) % wvel(LBi:UBi,LBj:UBj,0:N(ng)) ) allocate ( OCEAN(ng) % bed(LBi:UBi,LBj:UBj,Nbed,MBEDP) ) allocate ( OCEAN(ng) % bed_frac(LBi:UBi,LBj:UBj,Nbed,NST) ) allocate ( OCEAN(ng) % bed_mass(LBi:UBi,LBj:UBj,Nbed,2,NST) ) allocate ( OCEAN(ng) % ero_flux(LBi:UBi,LBj:UBj,NST) ) allocate ( OCEAN(ng) % settling_flux(LBi:UBi,LBj:UBj,NST) ) allocate ( OCEAN(ng) % bottom(LBi:UBi,LBj:UBj,MBOTP) ) RETURN END SUBROUTINE allocate_ocean SUBROUTINE initialize_ocean (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 USE mod_sediment ! ! Imported variable declarations. ! integer, intent(in) :: ng, tile, model ! ! Local variable declarations. ! integer :: Imin, Imax, Jmin, Jmax integer :: i, j, rec integer :: itrc, k 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, IstrR, IstrT, IstrU, Iend, IendR, IendT integer :: Jstr, JstrR, JstrT, JstrV, Jend, JendR, JendT ! Istr =BOUNDS(ng)%Istr (tile) IstrR=BOUNDS(ng)%IstrR(tile) IstrT=BOUNDS(ng)%IstrT(tile) IstrU=BOUNDS(ng)%IstrU(tile) Iend =BOUNDS(ng)%Iend (tile) IendR=BOUNDS(ng)%IendR(tile) IendT=BOUNDS(ng)%IendT(tile) Jstr =BOUNDS(ng)%Jstr (tile) JstrR=BOUNDS(ng)%JstrR(tile) JstrT=BOUNDS(ng)%JstrT(tile) JstrV=BOUNDS(ng)%JstrV(tile) Jend =BOUNDS(ng)%Jend (tile) JendR=BOUNDS(ng)%JendR(tile) JendT=BOUNDS(ng)%JendT(tile) ! ! 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 OCEAN(ng) % rubar(i,j,1) = IniVal OCEAN(ng) % rubar(i,j,2) = IniVal OCEAN(ng) % rvbar(i,j,1) = IniVal OCEAN(ng) % rvbar(i,j,2) = IniVal OCEAN(ng) % rzeta(i,j,1) = IniVal OCEAN(ng) % rzeta(i,j,2) = IniVal OCEAN(ng) % ubar(i,j,1) = IniVal OCEAN(ng) % ubar(i,j,2) = IniVal OCEAN(ng) % ubar(i,j,3) = IniVal OCEAN(ng) % vbar(i,j,1) = IniVal OCEAN(ng) % vbar(i,j,2) = IniVal OCEAN(ng) % vbar(i,j,3) = IniVal OCEAN(ng) % zeta(i,j,1) = IniVal OCEAN(ng) % zeta(i,j,2) = IniVal OCEAN(ng) % zeta(i,j,3) = IniVal END DO DO k=1,N(ng) DO i=Imin,Imax OCEAN(ng) % pden(i,j,k) = IniVal OCEAN(ng) % rho(i,j,k) = IniVal OCEAN(ng) % u(i,j,k,1) = IniVal OCEAN(ng) % u(i,j,k,2) = IniVal OCEAN(ng) % v(i,j,k,1) = IniVal OCEAN(ng) % v(i,j,k,2) = IniVal END DO END DO DO k=0,N(ng) DO i=Imin,Imax OCEAN(ng) % ru(i,j,k,1) = IniVal OCEAN(ng) % ru(i,j,k,2) = IniVal OCEAN(ng) % rv(i,j,k,1) = IniVal OCEAN(ng) % rv(i,j,k,2) = IniVal OCEAN(ng) % W(i,j,k) = IniVal OCEAN(ng) % wvel(i,j,k) = IniVal END DO END DO DO itrc=1,NT(ng) DO k=1,N(ng) DO i=Imin,Imax OCEAN(ng) % t(i,j,k,1,itrc) = IniVal OCEAN(ng) % t(i,j,k,2,itrc) = IniVal OCEAN(ng) % t(i,j,k,3,itrc) = IniVal END DO END DO END DO DO itrc=1,NST DO k=1,Nbed DO i=Imin,Imax OCEAN(ng) % bed_frac(i,j,k,itrc) = IniVal OCEAN(ng) % bed_mass(i,j,k,1,itrc) = IniVal OCEAN(ng) % bed_mass(i,j,k,2,itrc) = IniVal END DO END DO END DO DO itrc=1,MBEDP DO k=1,Nbed DO i=Imin,Imax OCEAN(ng) % bed(i,j,k,itrc) = IniVal END DO END DO END DO DO itrc=1,NST DO i=Imin,Imax OCEAN(ng) % ero_flux(i,j,itrc) = IniVal OCEAN(ng) % settling_flux(i,j,itrc) = IniVal END DO END DO DO itrc=1,MBOTP DO i=Imin,Imax OCEAN(ng) % bottom(i,j,itrc) = IniVal END DO END DO END DO END IF RETURN END SUBROUTINE initialize_ocean END MODULE mod_ocean