Back to home page

MITgcm

 
 

    


File indexing completed on 2018-03-02 18:37:38 UTC

view on githubraw file Latest commit 672ec4d5 on 2016-12-01 17:05:48 UTC
672ec4d5cf Jean*0001 <HTML>
                0002 <TITLE>module monin_obukhov</TITLE>
                0003 <BODY BGCOLOR="#AABBCC" TEXT="#332211" >
                0004 
                0005 <DIV ALIGN="CENTER"> <FONT SIZE=1>
                0006 <A HREF="#INTERFACE">PUBLIC INTERFACE</A> / 
                0007 <A HREF="#ROUTINES">ROUTINES</A> / 
                0008 <A HREF="#NAMELIST">NAMELIST</A> / 
                0009 <A HREF="#CHANGES">CHANGES</A> / 
                0010 <A HREF="#ERRORS">ERRORS</A> / 
                0011 <A HREF="#REFERENCES">REFERENCES</A> / 
                0012 <A HREF="#NOTES">NOTES</A> 
                0013 </FONT>
                0014 <BR><BR></DIV><HR>
                0015 
                0016 
                0017 <H2>module monin_obukhov</H2>
                0018 <A NAME="HEADER">
                0019 <PRE>
                0020      <B>Contact:</B>   Isaac Held
                0021      <B>Reviewers:</B>
                0022 
                0023      <B><A HREF=".doc.log#monin_obukhov.f90">Tags/Status</A></B>
                0024 </PRE>
                0025 </A><!-- END HEADER -->
                0026 <!--------------------------------------------------------------------->
                0027 <A NAME="OVERVIEW">
                0028 <HR>
                0029 <H4>OVERVIEW</H4>
                0030 <!-- BEGIN OVERVIEW -->
                0031 <PRE>
                0032 
                0033       Routines for computing surface drag coefficients from
                0034       data at the lowest model level using Monin-Obukhov scaling,
                0035       for estimating the wind, temperature, and buoyancy profiles
                0036       between the lowest model level and surface, and for computing
                0037       diffusivities consistent these profiles
                0038 
                0039 </PRE>
                0040 </A><!-- END OVERVIEW -->
                0041 <!--------------------------------------------------------------------->
                0042 <A NAME="DESCRIPTION">
                0043 <!-- BEGIN DESCRIPTION -->
                0044 <PRE>
                0045 
                0046   Monin-Obukhov similarity (MOS) theory is the standard method 
                0047   for computing  surface fluxes from the lowest level winds, 
                0048   temperatures, and tracer mixing ratios in GCMs.  The lowest level 
                0049   is assumed to lie within the "surface layer" in which turbulent 
                0050   fluxes have negligible vertical variation, and in which MOS assumes 
                0051   that the wind and buoyancy profiles are a function only of the 
                0052   surface stress, the surface buoyancy flux, and the height z. 
                0053   A good reference is Garratt, J. R. "The Atmospheric Boundary Layer", 
                0054   Cambridge University Press, 1992.  For a  discussion of the detailed 
                0055   form of the similarity theory used in this module, see  
                0056   <A HREF="monin_obukhov.tech.ps">monin_obukhov.tech.ps</A>
                0057 
                0058   The particular form of the monin-obukhov similarity theory utilized
                0059   is defined by the similarity functions phi_m and phi_t chosen, where
                0060           du/dzeta  = phi_m * u_star/( vonkarm * zeta )
                0061       db/dzeta( = phi_t * b_star/( vonkarm * zeta )
                0062         (where b = buoyancy or specific humidity, 
                0063                vonkarm = Von Karman's constant
                0064            u_star = friction velocity (stress = density*u_star**2)
                0065            b_star = buoyancy scale (flux = density*u_star*b_star)
                0066            zeta = z/L
                0067            z = height of lowest model level
                0068            L = monin_obukhov length )
                0069   We use:
                0070      on the unstable side 
                0071              phi_m = (1 - 16.0*zeta)**(-0.25)
                0072              phi_t = (1 - 16.0*zeta)**(-0.5)
                0073 
                0074      on the stable side -- phi_t = phi_m
                0075        there are two options
                0076        
                0077        option_1: phi =  1.0 + zeta*(5.0 + b_stab*zeta)/(1.0 + zeta)
                0078               where b_stab = 1/rich_crit 
                0079               (rich_crit is a namelist parameter)
                0080           
                0081        option_2: phi = 1 + 5.0*zeta (zeta < zeta_trans)
                0082                      = 1.0 + (5.0 - b_stab)*zeta_trans + b_stab*zeta
                0083               (zeta_trans also a namelist parameter)
                0084           
                0085   In order to change these similarity functions, one would need to 
                0086      1) modify the defintion of phi_m and phi_t in the internal
                0087         subroutines mo_derivative_m and mo_derivative_t
                0088      2) provide analytical versions of the integral similarity
                0089         functions  (PHI_m = int(zeta_0 to zeta) of phi/zeta)
                0090     (see notes) in the subroutines mo_integral_m and mo_integral_tq
                0091      3) modify stable_mix to be consistent with the new similarity functions
                0092         (see notes) (mod_diff will take care of itself)
                0093 
                0094 
                0095 </PRE>
                0096 </A><!-- END DESCRIPTION -->
                0097 <!--------------------------------------------------------------------->
                0098 <A NAME="MODULES_USED">
                0099 <HR>
                0100 <H4>OTHER MODULES USED</H4>
                0101 <!-- BEGIN MODULES_USED -->
                0102 <PRE>
                0103 
                0104     This modules uses 
                0105        constants_mod, only : grav, vonkarm
                0106        utilities_mod, only : error_mesg, FATAL, file_exist,   &
                0107                              check_nml_error, open_file,      &
                0108                              get_my_pe, close_file
                0109 
                0110 
                0111 </PRE>
                0112 </A><!-- END MODULES_USED -->
                0113 <!--------------------------------------------------------------------->
                0114 <A NAME="INTERFACE">
                0115 <HR>
                0116 <H4>PUBLIC INTERFACE</H4>
                0117 <!-- BEGIN INTERFACE -->
                0118 <PRE>
                0119 
                0120     use monin_obukhov_mod [,only: mo_drag    , &
                0121                                   mo_profile , &
                0122                                   mo_diff    , &
                0123                   stable_mix   ]
                0124 
                0125     mo_drag   : computes the drag coefficients for momentum, heat, and moisture;
                0126                 also returns u_star (the friction velocity) and
                0127                 b_star (the buoyancy scale)
                0128 
                0129     mo_profile :  provides the profile of momentum, temperature, and 
                0130                   specific humidity in the constant flux layer 
                0131           -- useful when one needs
                0132                   model output at a standard level below the lowest
                0133                   model level.  For example, the term "surface wind" 
                0134           in meteorology often refers to 
                0135                   the wind 10 meters above the surface.
                0136 
                0137     mo_diff : computes the difusivity that, acting in isolation,
                0138               would produce the correct profiles in ths surface layer
                0139          
                0140     stable_mix : returns f(Ri) under stable conditions where Ri is the 
                0141              local Richardson's number,  from which on can compute the 
                0142          diffusivity that  produces the correct profile in the constant 
                0143          flux layer from D = l**2 |dv\dz| F(Ri) where l = vonkarm*z .
                0144          Useful for matching a local Ri-dependent diffusivity in the free 
                0145          atmosphere under stable conditions o that implied by the 
                0146          surface layer profile
                0147 
                0148 
                0149 
                0150   Notes:
                0151 
                0152        all of the interfaces can be called with either 2d, 1d or 0d in/output.  
                0153        The base routines for mo_drag and mo_profile are 1d, while 
                0154        those for mo_diff and stable_mix are 2d, as these are the 
                0155        versions called within our GCMs.  (mo_drag and mo_profile are
                0156        called from the exchange_grid, which is 1d).  Other versions call these
                0157        base routines -- they are included for convenience -- the presumption is
                0158        that efficiency is not of particular importance for these alternative
                0159        versions)  
                0160 
                0161 
                0162 
                0163 
                0164 </PRE>
                0165 </A><!-- END INTERFACE -->
                0166 <!--------------------------------------------------------------------->
                0167 <A NAME="ROUTINES">
                0168 <HR>
                0169 <H4>PUBLIC ROUTINES</H4>
                0170 <!-- BEGIN ROUTINES -->
                0171 <PRE>
                0172 
                0173 ==========================================================================
                0174 
                0175  call  mo_drag (pt, pt0, z, z0, zt, zq, speed, drag_m, drag_t, drag_q, &
                0176                 u_star, b_star, [mask])
                0177 
                0178    input: 
                0179 
                0180       pt,  real, dimension(:) 
                0181          virtual potential temperature at lowest model level
                0182          degrees Kelvin
                0183 
                0184       pt0, real, dimension(:)
                0185          virtual potential temperature at surface
                0186          degrees Kelvin
                0187 
                0188       z, real, dimension(:) 
                0189          height above surface of lowest model layer 
                0190          meters
                0191 
                0192       z0, real, dimension(:)
                0193           surface roughness for momentum 
                0194           meters
                0195 
                0196       zt, real, dimension(:)
                0197           surface roughness for temperature
                0198           meters
                0199 
                0200       zq, real, dimension(:) 
                0201           surface roughness for moisture
                0202           meters
                0203 
                0204       speed, real, dimension(:) 
                0205           wind speed at lowest model level with respect to surface 
                0206              (any "gustiness" factor should be included in speed)
                0207           meters/sec
                0208 
                0209    output:
                0210 
                0211         drag_m, real, dimension(:) 
                0212               non-dimensional drag coefficient for momentum
                0213            
                0214         drag_t, real, dimension(:) 
                0215               non-dimensional drag coefficient for temperature
                0216            
                0217         drag_q, real, dimension(:) 
                0218               non-dimensional drag coefficient for specific humidity
                0219 
                0220         u_star, real, dimension(:) 
                0221            friction velocity 
                0222            meters/sec
                0223 
                0224         b_star, real, dimension(:) 
                0225            buoyancy scale 
                0226            (meters/sec)**2
                0227 
                0228 
                0229             The magnitude of the wind stress is 
                0230                  density*(ustar**2)
                0231             The drag coefficient for momentum is 
                0232                  (u_star/speed)**2
                0233             The buoyancy flux is
                0234                  density*ustar*bstar
                0235             The drag coefficient for heat etc is
                0236                  (u_star/speed)*(b_star/delta_b)
                0237                  where delta_b is the buoyancy difference between
                0238                   the surface and the lowest model level
                0239             The buoyancy == grav*(p0 - p)/p0
                0240 
                0241       
                0242     optional:
                0243        mask    : logical, dimension(:) 
                0244                    computation performed only where mask = .true.
                0245            
                0246           
                0247   Also:
                0248   
                0249      can be called with all 1d arrays replaced by 2d arrays or by 0d scalars
                0250      NOTE:  in the 0d and 2d cases, there is no avail optional argument
                0251 ==========================================================================
                0252 
                0253  subroutine mo_profile(zref_m, zref_t, z, z0, zt, zq, &
                0254       u_star, b_star, q_star, del_m, del_t, del_q, [mask])
                0255 
                0256   input:
                0257 
                0258      zref_m, real
                0259              height above surface to which interpolation is requested
                0260          for horizontal components of momentum
                0261              meters
                0262        
                0263      zref_t, real
                0264              height above surface to which interpolation is requested
                0265          for temperature and specific humidity
                0266              meters
                0267 
                0268      z, real, dimension(:) 
                0269         height of lowest model layer above surface
                0270         meters
                0271 
                0272      z0, real, dimension(:)
                0273          surface roughness for momentum 
                0274          meters
                0275 
                0276      zt, real, dimension(:) 
                0277          surface roughness for temperature
                0278          meters
                0279 
                0280      zq, real, dimension(:) 
                0281          surface roughness for moisture
                0282          meters
                0283 
                0284      u_star, real, dimension(:) 
                0285              friction velocity 
                0286              meters/sec
                0287 
                0288      b_star, real, dimension(:) 
                0289              buoyancy scale
                0290              (meters/sec)**2
                0291 
                0292      q_star, real, dimension(:)
                0293              moisture scale
                0294              kg/kg
                0295 
                0296           (Note:  u_star and b_star are output from mo_drag,
                0297                   q_star = flux_q/(u_star * density)
                0298 
                0299     optional input:
                0300 
                0301        mask, logical, dimension(:) 
                0302                    computation performed only where mask = .true.
                0303 
                0304 
                0305     output:
                0306 
                0307        del_m, real, dimension(:)
                0308               dimensionless ratio, as defined below, for momentum
                0309 
                0310        del_t, real, dimension(:)
                0311               dimensionless ratio, as defined below, for temperature
                0312 
                0313        del_q, real, dimension(:) o
                0314               dimensionless ratio, as defined below, for specific humidity
                0315 
                0316           Ratios are  (f(zref) - f_surf)/(f(z) - f_surf)
                0317       
                0318   Also:
                0319   
                0320      can be called with all 1d arrays replaced by 2d arrays or by 0d scalars
                0321      NOTE:  in the 0d and 2d cases, there is no avail optional argument
                0322      
                0323      in addition, this routine can be called with the interface 
                0324      
                0325      subroutine mo_profile(zref, z, z0, zt, zq, &
                0326            u_star, b_star, q_star, del_m, del_h, del_q)
                0327        
                0328       in which zref is a 1d array of the heights at which the profiles
                0329         are requested, and all other input is 0d, 1d, or 2d as above,
                0330     while the output del_m, del_h, del_q is dimensioned
                0331     in the 0d case:  dimension (size(zref))
                0332     in the 1d case:  dimension (size(input), size(zref))
                0333     in the 2d case:  dimension (size(input,1), size(input,2), size(zref))
                0334     (note that the two inputs zref_m and zref_h and here replaced by
                0335     the input vector zref)
                0336 
                0337 ==========================================================================
                0338                
                0339   subroutine mo_diff(z, u_star, b_star, k_m, k_h)
                0340 
                0341     input: 
                0342 
                0343         u_star, real, dimension(:,:) -- 2d horizontal array
                0344                 surface friction velocity
                0345                 (meters/sec)
                0346 
                0347         b_star, real, dimension(:,:) -- 2d horizontal array
                0348                 buoyancy scale 
                0349                 (meters/sec**2)
                0350 
                0351           (Note:  u_star and b_star are output from mo_drag)
                0352       
                0353     z, real, dimension(:,:,:)  
                0354             (first two dimensions conform to u_star, b_star)
                0355                 height above surface at which diffusivities 
                0356                 are desired -- i.e., diffusivities returned at 
                0357         z(i,j,:) for each point (i,j) 
                0358                 (meters)
                0359 
                0360     output:
                0361 
                0362         k_m   : real, dimension(:,:,:) conforms to z
                0363                 kinematic diffusivity for momentum 
                0364                 (meters**2/sec)
                0365 
                0366         k_h   : real, dimension(:,:,:)
                0367                 kinematic diffusivity for temperature
                0368                 (meters**2/sec)
                0369 
                0370 Also 0d and 1d versions available
                0371   0d:  u_star, b_star = scalars; 
                0372        z, k_m, k_h = dimension(:)
                0373   1d   u_star, b_star = dimension(:), 
                0374        z, k_m, k_h = dimension(size(u_star), :)
                0375        
                0376 Also, can be called so as to return diffusivities at one level 
                0377        (the height of this level can still vary from column to column)
                0378   0d:  u_star, b_star, z, k_m, k_h  = scalars
                0379   1d:  u_star, b_star, z, k_m, k_h  = dimension(:)
                0380   2d   u_star, b_star, z, k_m, k_h  = dimension(:,:)
                0381   
                0382 ==========================================================================
                0383                
                0384   subroutine stable_mix(rich, mix)
                0385 
                0386     input: 
                0387 
                0388         rich, real, dimension(:,:,:) 
                0389           local Richardson's number
                0390                 (none)
                0391 
                0392     output:
                0393 
                0394         mix,  real, dimension(:,:,:)
                0395           dimensional factor that produces the diffusivity
                0396           consistent with the similarity profile
                0397           under stable conditions when multiplied by 
                0398           (L**2) * |dv/dz| where v is the vector horizontal wind
                0399           and L = vonkarm*z
                0400 
                0401 Also 0d, 1d, and 2d versions vailable
                0402 
                0403 </PRE>
                0404 </A><!-- END ROUTINES -->
                0405 <!--------------------------------------------------------------------->
                0406 <A NAME="NAMELIST">
                0407 <HR>
                0408 <H4>NAMELIST</H4>
                0409 <!-- BEGIN NAMELIST -->
                0410 <PRE>
                0411 
                0412    (default values provided)
                0413    
                0414    logical :: neutral    = .false.
                0415 
                0416       If set to .true., all stability dependence is suppressed and 
                0417        neutral logarithmic profiles are used for all calculations;
                0418        all other namelist parameters ignored
                0419        
                0420    real  :: drag_min   = 1.e-05   (non-dimensional)
                0421 
                0422        The drag coefficients (for both momentum and temperature) 
                0423        are not allowed to fall below drag_min
                0424 
                0425    
                0426    integer :: stable_option = 1
                0427    
                0428        must equal either 1 or 2
                0429        two version for the similarity functions on the stable side
                0430 
                0431    real :: rich_crit  = 2.0      (non-dimensional)
                0432 
                0433        The first step in applying the similarity theory is computing
                0434        the bulk Richardson's  number Ri between the lowest model level
                0435        and the surface. If Ri is greater than rich_crit, the drag 
                0436        coefficients are set to drag_min. 
                0437        
                0438        This value also affects the drag coefficients under stable conditions
                0439        for either choice of stable_option, as the drag is arranged to
                0440        becomes very small as rich_crit is approached
                0441        
                0442    real :: zeta_trans = 0.5 (non-dimensional)
                0443     
                0444        A parameter in the stable-option = 2 piecewise linear similarity
                0445        function for the stable side.  Increasing zeta_trans reduces
                0446        the drag and the implied diffusivities
                0447          
                0448 
                0449 </PRE>
                0450 </A><!-- END NAMELIST -->
                0451 <!--------------------------------------------------------------------->
                0452 <A NAME="CHANGES">
                0453 <HR>
                0454 <H4>CHANGE HISTORY</H4>
                0455 <!-- BEGIN CHANGES -->
                0456 <PRE>
                0457 <B><A HREF=".doc.log#monin_obukhov.f90">Revision history</A></B>
                0458 
                0459 <b>changes</b> (10/4/1999)
                0460 
                0461      MPP version created. Minor changes for open_file, error_mesg,
                0462      and Fortran write statements. Answers should reproduce the
                0463      previous version.
                0464 
                0465 <b>more changes</b>  (10/01)
                0466 
                0467     -- time smoothing option removed
                0468     -- second optional similarity function on stable side added
                0469     -- stable_mix interface added
                0470     -- code rearranged substantially for clarity, removing initial
                0471        aggregation of points into unstable and stable sets
                0472     -- answers reproduce for stable_option = 1 except in the very 
                0473        rare instance of points that are very close to neutral.  
                0474        These need to be treated separately to avoid a divide by sero.  
                0475        In this version, the drag coefficients at these points are simply 
                0476        set to the neutral value
                0477 
                0478 </PRE>
                0479 </A><!-- END CHANGES -->
                0480 <!--------------------------------------------------------------------->
                0481 <A NAME="ERRORS">
                0482 <HR>
                0483 <H4>ERROR MESSAGES</H4>
                0484 <!-- BEGIN ERRORS -->
                0485 <PRE>
                0486 
                0487 FATAL ERRORS in MONIN_OBUKHOV_INIT in MONIN_OBUKHOV_MOD
                0488 
                0489     -- rich_crit in monin_obukhov_mod must be > 0.25
                0490 
                0491     -- drag_min in monin_obukhov_mod must be >= 0.0
                0492     
                0493     -- the only allowable values of stable_option are 1 and 2
                0494     
                0495     -- zeta_trans must be positive
                0496     
                0497 
                0498 FATAL ERROR in solve_zeta in monin_obukhov_mod
                0499 
                0500     -- surface drag iteration did not converge
                0501 
                0502 </PRE>
                0503 </A><!-- END ERRORS -->
                0504 <!--------------------------------------------------------------------->
                0505 <A NAME="REFERENCES">
                0506 <HR>
                0507 <H4>REFERENCES</H4>
                0508 <!-- BEGIN REFERENCES -->
                0509 <PRE>
                0510 
                0511 </PRE>
                0512 </A><!-- END REFERENCES -->
                0513 <!--------------------------------------------------------------------->
                0514 <A NAME="BUGS">
                0515 <HR>
                0516 <H4>KNOWN BUGS</H4>
                0517 <!-- BEGIN BUGS -->
                0518 <PRE>
                0519 None known
                0520 </PRE>
                0521 </A><!-- END BUGS -->
                0522 <!--------------------------------------------------------------------->
                0523 <A NAME="NOTES">
                0524 <HR>
                0525 <H4>NOTES</H4>
                0526 <!-- BEGIN NOTES -->
                0527 <PRE>
                0528 
                0529 If iteration convergence is ever a problem, one can consider
                0530 precomputation and table interpolation
                0531 
                0532 It would be convenient if changing the derivative similarity functions
                0533 automatically modified the integral functions and stable_mix
                0534 (this would presumably require numerical integration and table look-ups 
                0535 as opposed to analytical expressions)
                0536 
                0537 
                0538 </PRE>
                0539 </A><!-- END NOTES -->
                0540 <!--------------------------------------------------------------------->
                0541 <A NAME="PLANS">
                0542 <HR>
                0543 <H4>FUTURE PLANS</H4>
                0544 <!-- BEGIN PLANS -->
                0545 <PRE>
                0546 
                0547 
                0548 </PRE>
                0549 </A><!-- END PLANS -->
                0550 <!--------------------------------------------------------------------->
                0551 
                0552 <HR>
                0553 </BODY>
                0554 </HTML>