buoyancy.F90 5.31 KB
 dumoda01 committed Jan 19, 2011 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 !$Id: buoyancy.F90,v 1.8 2006-11-06 13:36:45 hb Exp$ #include"cppdefs.h" !----------------------------------------------------------------------- !BOP ! ! !ROUTINE: The buoyancy equation ! ! !INTERFACE: subroutine buoyancy(nlev,dt,cnpar,nub,gamb) ! ! !DESCRIPTION: ! This subroutine solves a transport equation for the mean ! potential buoyancy, ! ! \label{DefBuoyancy} ! B=-g\frac{\mean{\rho}-\rho_0}{\rho_0} ! \comma ! ! where $g$ is the accelaration of gravity, and $\mean{\rho}$ and $\rho_0$ ! are the mean potential density and the reference density, respectively. ! A simplified transport equation for $B$ can be written as ! ! \label{bEq} ! \dot{B} ! = {\cal D}_B ! \comma ! ! where $\dot{B}$ denotes the material derivative of $B$, and ! ${\cal D}_b$ is the sum of the turbulent and viscous transport ! terms modelled according to ! ! \label{Db} ! {\cal D}_B ! = \frstder{z} ! \left( (\nu^B_t+\nu^B) \partder{B}{z} - \tilde{\Gamma}_B \right) ! \point ! ! In this equation, $\nu^B_t$ and $\nu^B$ are the turbulent and molecular ! diffusivities of buoyancy, respectively, and $\tilde{\Gamma}_B$ ! denotes the non-local flux of buoyancy, see ! \sect{sec:turbulenceIntro}. In the current version of GOTM, ! we set $\nu^B_t = \nu^\Theta_t$ for simplicity. Source and sink ! terms are completely disregarded, and thus \eq{bEq} mainly serves ! as a convenient tool for some idealized test cases in ! GOTM. ! ! Diffusion is treated implicitly in space (see equations (\ref{sigmafirst})- ! (\ref{sigmalast})), and then solved by a ! simplified Gauss elimination. ! Vertical advection is included, and it must be non-conservative, ! which is ensured by setting the local variable {\tt adv\_mode=0}, ! see section \ref{sec:advectionMean} on page \pageref{sec:advectionMean}. ! ! !USES: use meanflow, only: h,w,buoy,T,avh use meanflow, only: w_grid,grid_method use observations, only: b_obs_NN,b_obs_surf,b_obs_sbf use observations, only: w_adv_discr,w_adv_method use util, only: Dirichlet,Neumann use util, only: oneSided,zeroDivergence ! IMPLICIT NONE ! ! !INPUT PARAMETERS: ! number of vertical layers integer, intent(in) :: nlev ! time step (s) REALTYPE, intent(in) :: dt ! numerical "implicitness" parameter REALTYPE, intent(in) :: cnpar ! diffusivity of buoyancy (m^2/s) REALTYPE, intent(in) :: nub(0:nlev) ! non-local buoyancy flux (m^2/s^3) REALTYPE, intent(in) :: gamb(0:nlev) ! ! ! !REVISION HISTORY: ! Original author(s): Hans Burchard & Karsten Bolding ! ! $Log: buoyancy.F90,v$ ! Revision 1.8 2006-11-06 13:36:45 hb ! Option for conservative vertical advection added to adv_center ! ! Revision 1.7 2005-11-17 09:58:20 hb ! explicit argument for positive definite variables in diff_center() ! ! Revision 1.6 2005/06/27 13:44:07 kbk ! modified + removed traling blanks ! ! Revision 1.5 2003/03/28 09:20:35 kbk ! added new copyright to files ! ! Revision 1.4 2003/03/28 08:56:56 kbk ! removed tabs ! ! Revision 1.3 2003/03/10 08:50:06 gotm ! Improved documentation and cleaned up code ! ! Revision 1.1.1.1 2001/02/12 15:55:57 gotm ! initial import into CVS ! !EOP ! ! !LOCAL VARIABLES: integer :: adv_mode=0 integer :: posconc=0 integer :: i integer :: DiffBcup,DiffBcdw integer :: AdvBcup,AdvBcdw REALTYPE :: DiffBup,DiffBdw REALTYPE :: AdvBup,AdvBdw REALTYPE :: Lsour(0:nlev) REALTYPE :: Qsour(0:nlev) REALTYPE :: BRelaxTau(0:nlev) REALTYPE :: zz logical, save :: first=.true. ! !----------------------------------------------------------------------- !BOC ! no diffusive flux at bottom and surface DiffBcup = Neumann DiffBcdw = Neumann DiffBup = _ZERO_ DiffBdw = _ZERO_ AdvBcup = zeroDivergence AdvBcdw = oneSided AdvBup = _ZERO_ AdvBdw = _ZERO_ ! construct initial linear profile from information in namelist if (first) then zz=_ZERO_ do i=nlev,1,-1 zz=zz+0.5*h(i) buoy(i) = b_obs_surf - zz*b_obs_NN zz=zz+0.5*h(i) end do first=.false. end if ! compose source term do i=1,nlev Lsour(i) = _ZERO_ Qsour(i) = - (gamB(i)-gamB(i-1))/h(i) end do ! no relaxation to observed values BRelaxTau = 1.e15 ! do advection step if (w_adv_method .ne. 0) then call adv_center(nlev,dt,h,h,w,AdvBcup,AdvBcdw, & AdvBup,AdvBdw,w_adv_discr,adv_mode,buoy) end if ! do diffusion step call diff_center(nlev,dt,cnpar,posconc,h,DiffBcup,DiffBcdw, & DiffBup,DiffBdw,avh,Lsour,Qsour, & BRelaxTau,buoy,buoy) ! T = buoy return end subroutine buoyancy !EOC !----------------------------------------------------------------------- ! Copyright by the GOTM-team under the GNU Public License - www.gnu.org !-----------------------------------------------------------------------