MODULE mod_tides ! !svn $Id: mod_tides.F 300 2009-01-22 18:26:41Z arango $ !================================================== Hernan G. Arango === ! Copyright (c) 2002-2009 The ROMS/TOMS Group ! ! Licensed under a MIT/X style license ! ! See License_ROMS.txt ! !======================================================================= ! ! ! Tidal Components: ! ! ! ! Each of the following arrays has a dimension in tidal components ! ! classified by period: ! ! ! ! semi-diurnal: M2, S2, N2, K2 (12.42, 12.00, 12.66, 11.97h) ! ! diurnal: K1, O1, P1, Q1 (23.93, 25.82, 24.07, 26.87h) ! ! ! ! and other longer periods. The order of these tidal components is ! ! irrelevant here. The number of components to USE is depends on ! ! the regional application. ! ! ! ! CosOmega Cosine tidal harmonics for current omega(t). ! ! SinOmega Sine tidal harmonics for current omega(t). ! ! SSH_Tamp Tidal elevation amplitude (m) at RHO-points. ! ! SSH_Tphase Tidal elevation phase (degrees/360) at RHO-points. ! ! Tperiod Tidal period (s). ! ! UV_Tangle Tidal current angle (radians; counterclockwise ! ! from EAST and rotated to curvilinear grid) at ! ! RHO-points. ! ! UV_Tmajor Maximum tidal current: tidal ellipse major axis ! ! (m/s) at RHO-points. ! ! UV_Tminor Minimum tidal current: tidal ellipse minor axis ! ! (m/s) at RHO-points. ! ! UV_Tphase Tidal current phase (degrees/360) at RHO-points. ! ! ! !======================================================================= ! USE mod_kinds implicit none TYPE T_TIDES real(r8), pointer :: Tperiod(:) real(r8), pointer :: SSH_Tamp(:,:,:) real(r8), pointer :: SSH_Tphase(:,:,:) real(r8), pointer :: UV_Tangle(:,:,:) real(r8), pointer :: UV_Tmajor(:,:,:) real(r8), pointer :: UV_Tminor(:,:,:) real(r8), pointer :: UV_Tphase(:,:,:) END TYPE T_TIDES TYPE (T_TIDES), allocatable :: TIDES(:) CONTAINS SUBROUTINE allocate_tides (ng, LBi, UBi, LBj, UBj) ! !======================================================================= ! ! ! This routine allocates all variables in the module for all nested ! ! grids. ! ! ! !======================================================================= ! USE mod_param USE mod_parallel USE mod_iounits USE mod_ncparam USE mod_netcdf USE mod_scalars USE mod_stepping ! ! Inported variable declarations. ! integer, intent(in) :: ng, LBi, UBi, LBj, UBj ! ! Local variable declarations. ! logical :: foundit integer :: Vid, i, ifile, mg, nvatt, nvdim ! !----------------------------------------------------------------------- ! Allocate module variables. !----------------------------------------------------------------------- ! ! Inquire about the maximum number of tidal components. ! IF (ng.eq.1) THEN MTC=0 DO mg=1,Ngrids foundit=.FALSE. QUERY : DO ifile=1,nFfiles(mg) CALL netcdf_inq_var (ng, iNLM, FRCname(ifile,mg), & & MyVarName = TRIM(Vname(1,idTper)), & & SearchVar = foundit, & & VarID = Vid, & & nVardim = nvdim, & & nVarAtt = nvatt) IF (exit_flag.ne.NoError) RETURN ! ! Set maximun numbert of tidal components. ! IF (foundit) THEN MTC=MAX(MTC,var_Dsize(1)) ! first dimension NTC(mg)=var_Dsize(1) TIDEname(ng)=TRIM(FRCname(ifile,mg)) EXIT QUERY END IF END DO QUERY END DO END IF ! ! Allocate structure. ! IF (ng.eq.1) allocate ( TIDES(Ngrids) ) ! ! Allocate tidal forcing variables. ! allocate ( TIDES(ng) % Tperiod(MTC) ) allocate ( TIDES(ng) % SSH_Tamp(LBi:UBi,LBj:UBj,MTC) ) allocate ( TIDES(ng) % SSH_Tphase(LBi:UBi,LBj:UBj,MTC) ) allocate ( TIDES(ng) % UV_Tangle(LBi:UBi,LBj:UBj,MTC) ) allocate ( TIDES(ng) % UV_Tmajor(LBi:UBi,LBj:UBj,MTC) ) allocate ( TIDES(ng) % UV_Tminor(LBi:UBi,LBj:UBj,MTC) ) allocate ( TIDES(ng) % UV_Tphase(LBi:UBi,LBj:UBj,MTC) ) RETURN END SUBROUTINE allocate_tides SUBROUTINE initialize_tides (ng, tile) ! !======================================================================= ! ! ! 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 ! ! Local variable declarations. ! integer :: Imin, Imax, Jmin, Jmax integer :: i, itide, j, jtide, 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. !----------------------------------------------------------------------- ! IF (.TRUE.) THEN DO itide=1,MTC TIDES(ng) % Tperiod(itide) = IniVal END DO END IF DO itide=1,MTC DO j=Jmin,Jmax DO i=Imin,Imax TIDES(ng) % SSH_Tamp(i,j,itide) = IniVal TIDES(ng) % SSH_Tphase(i,j,itide) = IniVal END DO END DO DO j=Jmin,Jmax DO i=Imin,Imax TIDES(ng) % UV_Tangle(i,j,itide) = IniVal TIDES(ng) % UV_Tmajor(i,j,itide) = IniVal TIDES(ng) % UV_Tminor(i,j,itide) = IniVal TIDES(ng) % UV_Tphase(i,j,itide) = IniVal END DO END DO END DO RETURN END SUBROUTINE initialize_tides END MODULE mod_tides