module roelvink_module contains subroutine roelvink1(E,hh,Trep,alpha,gamma,n,rho,g,delta,D,ntot,break) !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! Copyright (C) 2007 UNESCO-IHE, WL|Delft Hydraulics and Delft University ! ! Dano Roelvink, Ap van Dongeren, Ad Reniers, Jamie Lescinski, ! ! Jaap van Thiel de Vries, Robert McCall ! ! ! ! d.roelvink@unesco-ihe.org ! ! UNESCO-IHE Institute for Water Education ! ! P.O. Box 3015 ! ! 2601 DA Delft ! ! The Netherlands ! ! ! ! This library is free software; you can redistribute it and/or ! ! modify it under the terms of the GNU Lesser General Public ! ! License as published by the Free Software Foundation; either ! ! version 2.1 of the License, or (at your option) any later version. ! ! ! ! This library is distributed in the hope that it will be useful, ! ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU ! ! Lesser General Public License for more details. ! ! ! ! You should have received a copy of the GNU Lesser General Public ! ! License along with this library; if not, write to the Free Software ! ! Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 ! ! USA ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! IMPLICIT NONE integer :: ntot character(24) :: break real*8,dimension(ntot) :: E,hh,D real*8 :: fac,rho,g,delta,alpha,gamma,n,Trep real*8,dimension(:),allocatable,save :: H,Qb,hroelvink if (.not. allocated(H)) then allocate(H (ntot)) allocate(Qb (ntot)) allocate(hroelvink(ntot)) endif ! Dissipation acc. to Roelvink (1993) fac=8.0d0/rho/g H=sqrt(fac*E) hroelvink=hh+delta*H !add breakerdelay here, to correct for waterdepth, and cancel subroutine breakerdelay !danojaaprobertap 10/2, also in baldock, roelvink Qb=min(1-exp(-(H/gamma/hroelvink)**n),1.0d0) !gebruik hier Hb = (0.88d0/k)*tanh(gamma*kh/0.88d0) ! D=Qb*2.*alpha/Trep*E ! cjaap : two options: if (trim(break)=='roelvink1') then D=Qb*2*alpha/Trep*E elseif (trim(break)=='roelvink2') then ! D=Qb*2.*alpha/Trep*E*H/hroelvink; !breaker delay depth for dissipation in bore and fraction of breaking waves D=Qb*2*alpha/Trep*E*H/hh !breaker delay depth for fraction breaking waves and actual water depth for disipation in bore end if end subroutine roelvink1 subroutine roelvink(par,s,km) use params use spaceparams IMPLICIT NONE type(spacepars) :: s type(parameters) :: par real*8 :: fac integer :: i,j real*8,dimension(s%nx+1,s%ny+1) :: km real*8,dimension(:,:),allocatable,save :: H,hroelvink,arg,kmr if (.not. allocated(H)) then allocate(H (s%nx+1,s%ny+1)) allocate(arg (s%nx+1,s%ny+1)) allocate(kmr (s%nx+1,s%ny+1)) allocate(hroelvink(s%nx+1,s%ny+1)) endif ! Dissipation acc. to Roelvink (1993) fac=8.0d0/par%rho/par%g hroelvink=s%hh+par%delta*s%H kmr=km where(km<0.01d0) kmr=0.01d0 elsewhere(km>100.d0) kmr=100.d0 endwhere H=sqrt(fac*s%E) if (trim(par%break)/='roelvink_daly') then if (par%wci==1) then arg = -(H/(par%gamma*tanh(min(max(km,0.01d0),100.d0)*hroelvink)/min(max(km,0.01d0),100.d0)))**par%n s%Qb = min(1.d0-exp(max(arg,-100.d0)),1.d0) else s%Qb=min(1-exp(-(H/par%gamma/hroelvink)**par%n),1.0d0) endif else do j=1,s%ny+1 do i=1,s%nx+1 if (H(i,j)>par%gamma *hroelvink(i,j)) s%Qb(i,j)=1.d0 if (H(i,j)