#include "cppdefs.h" MODULE frc_inw_mod #if defined INWAVE_MODEL # if defined WEC_MELLOR || defined WEC_VF ! !svn $Id: celer_inw.F 732 2008-09-07 01:55:51Z jcwarner $ !================================================== John C. Warner ===== ! ! ! This routine computes the energy dissipation ! ! ! !======================================================================= ! implicit none PRIVATE PUBLIC :: frc_inw CONTAINS ! !*********************************************************************** SUBROUTINE frc_inw (ng, tile) !*********************************************************************** ! USE mod_param USE mod_inwave_params USE mod_inwave_vars USE mod_ocean USE mod_coupling USE mod_stepping USE mod_forces ! ! Imported variable declarations. ! integer, intent(in) :: ng, tile ! ! Local variable declarations. ! # include "tile.h" ! !# ifdef PROFILE ! CALL wclock_on (ng, iNLM, 35) !# endif CALL frc_inw_tile(ng, tile, & & LBi, UBi, LBj, UBj, & & IminS, ImaxS, JminS, JmaxS, & & FORCES(ng) % Hwave, & & FORCES(ng) % Dwave, & & FORCES(ng) % Lwave, & # if defined WAVES_UB & FORCES(ng) % Uwave_rms, & # endif # if defined WAVES_TOP_PERIOD & FORCES(ng) % Pwave_top, & # endif # if defined WAVES_BOT_PERIOD & FORCES(ng) % Pwave_bot, & # endif & WAVEP(ng)% AC, WAVEP(ng)% Tr, & & WAVEP(ng)% Ta, WAVEP(ng)% kwc, & & WAVEG(ng)% WD, & & WAVEP(ng) % h_tot) !# ifdef PROFILE ! CALL wclock_off (ng, iNLM, 35) !# endif RETURN END SUBROUTINE frc_inw ! !*********************************************************************** SUBROUTINE frc_inw_tile(ng, tile, & & LBi, UBi, LBj, UBj, & & IminS, ImaxS, JminS, JmaxS, & & Hwave, Dwave, Lwave, & # if defined WAVES_UB & Uwave_rms, & # endif # if defined WAVES_TOP_PERIOD & Pwave_top, & # endif # if defined WAVES_BOT_PERIOD & Pwave_bot, & # endif & AC, Tr, Ta, kwc, & & WD, h_tot) !*********************************************************************** ! USE mod_param USE mod_inwave_params USE mod_boundary USE mod_grid USE mod_scalars USE mod_stepping USE exchange_2d_mod # ifdef DISTRIBUTE USE mp_exchange_mod, ONLY : mp_exchange2d # endif USE bc_2d_mod ! ! Imported variable declarations. ! integer, intent(in) :: ng, tile integer, intent(in) :: LBi, UBi, LBj, UBj integer, intent(in) :: IminS, ImaxS, JminS, JmaxS ! # ifdef ASSUMED_SHAPE real(r8), intent(inout) :: Hwave(LBi:,LBj:) real(r8), intent(inout) :: Dwave(LBi:,LBj:) real(r8), intent(inout) :: Lwave(LBi:,LBj:) # if defined WAVES_UB real(r8), intent(inout) :: Uwave_rms(LBi:,LBj:) # endif # if defined WAVES_TOP_PERIOD real(r8), intent(inout) :: Pwave_top(LBi:,LBj:) # endif # if defined WAVES_BOT_PERIOD real(r8), intent(inout) :: Pwave_bot(LBi:,LBj:) # endif real(r8), intent(in) :: AC(LBi:,LBj:,:,:) real(r8), intent(in) :: Tr(LBi:,LBj:,:) real(r8), intent(in) :: Ta(LBi:,LBj:,:) real(r8), intent(in) :: kwc(LBi:,LBj:,:) real(r8), intent(in) :: WD(:) real(r8), intent(in) :: h_tot(LBi:,LBj:) # else real(r8), intent(inout) :: Hwave(LBi:UBi,LBj:UBj) real(r8), intent(inout) :: Dwave(LBi:UBi,LBj:UBj) real(r8), intent(inout) :: Lwave(LBi:UBi,LBj:UBj) # if defined WAVES_UB real(r8), intent(inout) :: Uwave_rms(LBi:UBi,LBj:UBj) # endif # if defined WAVES_TOP_PERIOD real(r8), intent(inout) :: Pwave_top(LBi:UBi,LBj:UBj) # endif # if defined WAVES_BOT_PERIOD real(r8), intent(inout) :: Pwave_bot(LBi:UBi,LBj:UBj) # endif real(r8), intent(in) :: AC(LBi:UBi,LBj:UBj,ND,3) real(r8), intent(in) :: Tr(LBi:UBi,LBj:UBj,ND) real(r8), intent(in) :: Ta(LBi:UBi,LBj:UBj,ND) real(r8), intent(in) :: kwc(LBi:UBi,LBj:UBj,ND) real(r8), intent(in) :: WD(ND) real(r8), intent(in) :: h_tot(LBi:UBi,LBj:UBj) # endif ! ! Local variable declarations. ! integer :: i, j, d real(r8) :: TRM, H, DIRE, DIRN, twopi real(r8) :: EN, EW, oEW, KM, DwaveE, DwaveN real(r8), parameter :: Lmin=10.0_r8 real(r8), parameter :: Lmax=1000.0_r8 real(r8), parameter :: Pmin=2.0_r8 real(r8), parameter :: Trmin=1.0_r8 # if defined WAVES_UB real(r8), parameter :: eps = 1.0E-10_r8 real(r8), parameter :: K1 = 0.6666666666_r8 ! Coefficients for real(r8), parameter :: K2 = 0.3555555555_r8 ! explicit real(r8), parameter :: K3 = 0.1608465608_r8 ! wavenumber real(r8), parameter :: K4 = 0.0632098765_r8 ! calculation real(r8), parameter :: K5 = 0.0217540484_r8 ! (Dean and real(r8), parameter :: K6 = 0.0065407983_r8 ! Dalrymple, 1991) real(r8) :: Kbh, Kbh2, Kdh real(r8) :: Fwave_bot, Ab # endif # include "set_bounds.h" twopi=2.0_r8*pi DO j=Jstr,Jend DO i=Istr,Iend EW=0.0_r8 KM=0.0_r8 TRM=0.0_r8 DIRE=0.0_r8 DIRN=0.0_r8 IF (h_tot(i,j).ge.Dcrit(ng)) THEN DO d=1,ND !======================================================================= ! Compute the energy from the action balance and the wave height !======================================================================= EN=MAX(0.0_r8,AC(i,j,d,nnew(ng))*twopi/ & & MAX(Trmin,Tr(i,j,d))) !======================================================================= ! Compute the total energy !======================================================================= EW=EW+EN !======================================================================= ! Compute the mean wave number and intrinsic periods ! What we do is give more importance to those wave ! numbers with more energy !======================================================================= KM=KM+kwc(i,j,d)*EN TRM=TRM+Tr(i,j,d)*EN ! DIRM=DIRM+WD(d)*EN DIRE=DIRE+sin(WD(d))*EN DIRN=DIRN+cos(WD(d))*EN ENDDO IF (EW.ge.0.001_r8) THEN oEW=1.0_r8/EW KM=KM*oEW TRM=TRM*oEW Lwave(i,j)=MIN(twopi/KM,Lmax) DwaveE=DIRE*oEW DwaveN=DIRN*oEW Dwave(i,j)=ATAN2(DwaveE,DwaveN) IF (Dwave(i,j).lt.0.0_r8) THEN Dwave(i,j)=Dwave(i,j)+2.0_r8*pi END IF !======================================================================= ! Compute the wave height Hwave (Significant). !======================================================================= ! Hwave(i,j)=(8.0_r8*EW/(g*rho0))**0.5_r8 Hwave(i,j)=(16.0_r8*EW/(g*rho0))**0.5_r8 # if defined WAVES_TOP_PERIOD Pwave_top(i,j)=TRM # endif # if defined WAVES_BOT_PERIOD Pwave_bot(i,j)=TRM # endif # if defined WAVES_UB Fwave_bot=twopi/MAX(Pwave_bot(i,j),0.05_r8) Kdh=h_tot(i,j)*Fwave_bot**2/g Kbh2=Kdh*Kdh+ & & Kdh/(1.0_r8+Kdh*(K1+Kdh*(K2+Kdh*(K3+Kdh*(K4+ & & Kdh*(K5+K6*Kdh)))))) Kbh=SQRT(Kbh2) Ab=0.5_r8*Hwave(i,j)/SINH(Kbh)+eps Uwave_rms(i,j)=Fwave_bot*Ab+eps # endif ELSE Hwave(i,j)=0.0_r8 Dwave(i,j)=0.0_r8 Lwave(i,j)=Lmin # if defined WAVES_UB Uwave_rms(i,j)=0.0_r8 # endif # if defined WAVES_TOP_PERIOD Pwave_top(i,j)=Pmin # endif # if defined WAVES_BOT_PERIOD Pwave_bot(i,j)=Pmin # endif ENDIF ELSE Hwave(i,j)=0.0_r8 Dwave(i,j)=0.0_r8 Lwave(i,j)=Lmin # if defined WAVES_UB Uwave_rms(i,j)=0.0_r8 # endif # if defined WAVES_TOP_PERIOD Pwave_top(i,j)=Pmin # endif # if defined WAVES_BOT_PERIOD Pwave_bot(i,j)=Pmin # endif ENDIF ENDDO ENDDO ! ! Apply boundary conditions. ! CALL bc_r2d_tile (ng, tile, & & LBi, UBi, LBj, UBj, & & Hwave) CALL bc_r2d_tile (ng, tile, & & LBi, UBi, LBj, UBj, & & Dwave) CALL bc_r2d_tile (ng, tile, & & LBi, UBi, LBj, UBj, & & Lwave) # if defined WAVES_UB CALL bc_r2d_tile (ng, tile, & & LBi, UBi, LBj, UBj, & & Uwave_rms) # endif # if defined WAVES_TOP_PERIOD CALL bc_r2d_tile (ng, tile, & & LBi, UBi, LBj, UBj, & & Pwave_top) # endif # if defined WAVES_BOT_PERIOD CALL bc_r2d_tile (ng, tile, & & LBi, UBi, LBj, UBj, & & Pwave_bot) # endif IF (EWperiodic(ng).or.NSperiodic(ng)) THEN CALL exchange_r2d_tile (ng, tile, & & LBi, UBi, LBj, UBj, & & Hwave) CALL exchange_r2d_tile (ng, tile, & & LBi, UBi, LBj, UBj, & & Dwave) CALL exchange_r2d_tile (ng, tile, & & LBi, UBi, LBj, UBj, & & Lwave) # if defined WAVES_UB CALL exchange_r2d_tile (ng, tile, & & LBi, UBi, LBj, UBj, & & Uwave_rms) # endif # if defined WAVES_TOP_PERIOD CALL exchange_r2d_tile (ng, tile, & & LBi, UBi, LBj, UBj, & & Pwave_top) # endif # if defined WAVES_BOT_PERIOD CALL exchange_r2d_tile (ng, tile, & & LBi, UBi, LBj, UBj, & & Pwave_bot) # endif END IF # ifdef DISTRIBUTE ! ! Exchange boundary data. ! CALL mp_exchange2d (ng, tile, iNLM, 3, & & LBi, UBi, LBj, UBj, & & NghostPoints, & & EWperiodic(ng), NSperiodic(ng), & & Hwave, Dwave, Lwave) # if defined WAVES_UB CALL mp_exchange2d (ng, tile, iNLM, 1, & & LBi, UBi, LBj, UBj, & & NghostPoints, & & EWperiodic(ng), NSperiodic(ng), & & Uwave_rms) # endif # if defined WAVES_TOP_PERIOD CALL mp_exchange2d (ng, tile, iNLM, 1, & & LBi, UBi, LBj, UBj, & & NghostPoints, & & EWperiodic(ng), NSperiodic(ng), & & Pwave_top) # endif # if defined WAVES_BOT_PERIOD CALL mp_exchange2d (ng, tile, iNLM, 1, & & LBi, UBi, LBj, UBj, & & NghostPoints, & & EWperiodic(ng), NSperiodic(ng), & & Pwave_bot) # endif # endif RETURN END SUBROUTINE frc_inw_tile # endif #endif END MODULE frc_inw_mod