XBeach
C:/repositories/XBeach/trunk/src/xbeachlibrary/roelvink.F90
Go to the documentation of this file.
00001 module roelvink_module
00002    implicit none
00003    save
00004    private
00005    public roelvink, baldock, janssen_battjes
00006    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00007    ! Copyright (C) 2007 UNESCO-IHE, WL|Delft Hydraulics and Delft University !
00008    ! Dano Roelvink, Ap van Dongeren, Ad Reniers, Jamie Lescinski,            !
00009    ! Jaap van Thiel de Vries, Robert McCall                                  !
00010    !                                                                         !
00011    ! d.roelvink@unesco-ihe.org                                               !
00012    ! UNESCO-IHE Institute for Water Education                                !
00013    ! P.O. Box 3015                                                           !
00014    ! 2601 DA Delft                                                           !
00015    ! The Netherlands                                                         !
00016    !                                                                         !
00017    ! This library is free software; you can redistribute it and/or           !
00018    ! modify it under the terms of the GNU Lesser General Public              !
00019    ! License as published by the Free Software Foundation; either            !
00020    ! version 2.1 of the License, or (at your option) any later version.      !
00021    !                                                                         !
00022    ! This library is distributed in the hope that it will be useful,         !
00023    ! but WITHOUT ANY WARRANTY; without even the implied warranty of          !
00024    ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU        !
00025    ! Lesser General Public License for more details.                         !
00026    !                                                                         !
00027    ! You should have received a copy of the GNU Lesser General Public        !
00028    ! License along with this library; if not, write to the Free Software     !
00029    ! Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307     !
00030    ! USA                                                                     !
00031    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00032 
00033    interface roelvink
00034       module procedure roelvink_1D
00035       module procedure roelvink_2D
00036    end interface roelvink
00037 
00038    interface baldock
00039       module procedure baldock_1D
00040       module procedure baldock_2D
00041    end interface baldock
00042 
00043    interface janssen_battjes
00044       module procedure janssen_battjes_1D
00045       module procedure janssen_battjes_2D
00046    end interface janssen_battjes
00047 
00048 contains
00049 
00050    subroutine roelvink_1D(par,s,i)
00051 
00052       use params,         only: parameters
00053       use spaceparams
00054       use paramsconst
00055 
00056       implicit none
00057 
00058       type(spacepars)                                 :: s
00059       type(parameters)                                :: par
00060       integer,intent(in)                              :: i
00061 
00062       integer                                         :: j
00063 
00064       real*8, dimension(s%ny+1)                       ::kmr,hr,H,arg
00065 
00066       ! Dissipation according to Roelvink (1993)
00067 
00068       H   = s%H(i,:)
00069       hr  = s%hhw(i,:)
00070 
00071       kmr = min(max(s%k(i,:), 0.01d0), 100.d0)
00072 
00073       if (par%break /= BREAK_ROELVINK_DALY) then
00074          if (par%wci==1) then
00075             arg = -( H / (par%gamma*tanh(kmr*hr)/kmr))**par%n
00076          else
00077             arg = -( H / (par%gamma*hr              ))**par%n
00078          endif
00079 
00080          s%Qb(i,:) = min(1.d0 - exp(max(arg,-100.d0)), 1.d0)
00081       else
00082          do j=1,s%ny+1
00083             if (H(j) > par%gamma *hr(j)) s%Qb(i,j) = 1.d0
00084             if (H(j) < par%gamma2*hr(j)) s%Qb(i,j) = 0.d0
00085          enddo
00086 
00087          s%Qb(i,:) = max(s%Qb(i,:), 0.d0)
00088       endif
00089 
00090       s%D(i,:) = s%Qb(i,:) * 2.d0 * par%alpha * s%E(i,:)
00091 
00092       if (par%break == BREAK_ROELVINK1) then
00093          if (par%wci==1) then
00094             s%D(i,:) = s%D(i,:) * s%sigm(i,:)/2.d0/par%px;
00095          else
00096             s%D(i,:) = s%D(i,:) / par%Trep
00097          endif
00098       elseif (par%break == BREAK_ROELVINK2 .or. par%break == BREAK_ROELVINK_DALY) then
00099          ! Jaap: also wci switch for roelvink2
00100          if (par%wci==1) then
00101             s%D(i,:) = s%D(i,:) * s%sigm(i,:)/2.d0/par%px * H/s%hh(i,:);
00102          else
00103             s%D(i,:) = s%D(i,:) / par%Trep * H/s%hh(i,:)
00104          endif
00105       end if
00106 
00107    end subroutine roelvink_1D
00108 
00109    subroutine roelvink_2D(par,s)
00110 
00111       use params,         only: parameters
00112       use spaceparams
00113 
00114       implicit none
00115 
00116       type(spacepars)                                 :: s
00117       type(parameters)                                :: par
00118 
00119       integer                                         :: i
00120 
00121       do i = 1,s%nx+1
00122          call roelvink_1D(par,s,i)
00123       enddo
00124 
00125    end subroutine roelvink_2D
00126 
00127    subroutine baldock_1D(par,s,i)
00128 
00129       use params,        only: parameters
00130       use spaceparams
00131 
00132       implicit none
00133 
00134       type(spacepars)                 :: s
00135       type(parameters)                :: par
00136 
00137       integer,intent(in)              :: i
00138 
00139       real*8, dimension(s%ny+1)       :: kh,f,k,H,Hb,R,gamma
00140 
00141       ! Dissipation according to Baldock et al. (1998)
00142 
00143       k = s%k(i,:)
00144       if (par%wci==1) then
00145          f = s%sigm(i,:) / 2.d0 / par%px
00146       else
00147          f = 1.d0 / par%Trep
00148       endif
00149 
00150       kh  = s%k(i,:) * s%hhw(i,:)
00151 
00152       if (par%wci == 1) then
00153          gamma = 0.76d0*kh + 0.29d0 !Jaap: spatial varying gamma according to Ruessink et al., 1998
00154       else
00155          gamma = par%gamma
00156       endif
00157 
00158       H   = s%H(i,:)
00159       Hb  = tanh(gamma*kh/0.88d0)*(0.88d0/max(k,1d-10))
00160       R   = Hb/max(H,0.00001d0)
00161 
00162       s%Qb(i,:)   = exp(-R**2)
00163       s%D (i,:)   = 0.25d0 * par%alpha * f * par%rho * par%g * (Hb**2+H**2) * s%Qb(i,:)
00164 
00165    end subroutine baldock_1D
00166 
00167    subroutine baldock_2D(par,s)
00168 
00169       use params,        only: parameters
00170       use spaceparams
00171 
00172       implicit none
00173 
00174       type(spacepars)                                 :: s
00175       type(parameters)                                :: par
00176 
00177       integer                                         :: i
00178 
00179       do i = 1,s%nx+1
00180          call baldock_1D(par,s,i)
00181       enddo
00182 
00183    end subroutine baldock_2D
00184 
00185    subroutine janssen_battjes_1D(par,s,i)
00186 
00187       use params,         only: parameters
00188       use spaceparams
00189       use math_tools, only: xerf
00190 
00191       implicit none
00192 
00193       type(spacepars)                                 :: s
00194       type(parameters)                                :: par
00195 
00196       integer,intent(in)                              :: i
00197       real*8                                          :: B
00198 
00199       real*8, dimension(s%ny+1)                       :: kh,f,k,H,Hb,R
00200 
00201       ! Dissipation according to Janssen and Battjes (2007)
00202 
00203       B   = par%alpha
00204 
00205       k = s%k(i,:)
00206       if (par%wci==1) then
00207          f = s%sigm(i,:) / 2.d0 / par%px
00208       else
00209          f = 1.d0 / par%Trep
00210       endif
00211 
00212       kh  = s%k(i,:) * s%hhw(i,:)
00213 
00214       H   = s%H(i,:)
00215       Hb  = tanh(par%gamma*kh/0.88d0)*(0.88d0/k)
00216       R   = Hb/max(H,0.00001d0)
00217 
00218       s%Qb(i,:)   = 1 + 4/(3*sqrt(par%px)) * (R**3 + 1.5d0*R) * exp(-R**2) - xerf(R)
00219       ! wwvv was:
00220       !s%Qb(i,:)   = 1 + 4/(3*sqrt(pi)) * (R**3 + 3/2*R) * exp(-R**2) - xerf(R)
00221       s%D (i,:)   = 3*sqrt(par%px)/16      * B * f * par%rho * par%g * H**3/s%hh(i,:) * s%Qb(i,:)
00222 
00223    end subroutine janssen_battjes_1D
00224 
00225    subroutine janssen_battjes_2D(par,s)
00226 
00227       use params,         only: parameters
00228       use spaceparams
00229 
00230       implicit none
00231 
00232       type(spacepars)                                 :: s
00233       type(parameters)                                :: par
00234 
00235       integer                                         :: i
00236 
00237       do i = 1,s%nx+1
00238          call janssen_battjes_1D(par,s,i)
00239       enddo
00240 
00241    end subroutine janssen_battjes_2D
00242 
00243 end module roelvink_module
 All Classes Files Functions Variables Defines