FV3 Bundle
monin_obukhov.F90
Go to the documentation of this file.
1 !***********************************************************************
2 !* GNU Lesser General Public License
3 !*
4 !* This file is part of the GFDL Flexible Modeling System (FMS).
5 !*
6 !* FMS is free software: you can redistribute it and/or modify it under
7 !* the terms of the GNU Lesser General Public License as published by
8 !* the Free Software Foundation, either version 3 of the License, or (at
9 !* your option) any later version.
10 !*
11 !* FMS is distributed in the hope that it will be useful, but WITHOUT
12 !* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13 !* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
14 !* for more details.
15 !*
16 !* You should have received a copy of the GNU Lesser General Public
17 !* License along with FMS. If not, see <http://www.gnu.org/licenses/>.
18 !***********************************************************************
19 
21 
22 
23 !=======================================================================
24 !
25 ! MONIN-OBUKHOV MODULE
26 !
27 ! Routines for computing surface drag coefficients
28 ! from data at the lowest model level
29 ! and for computing the profile of fields
30 ! between the lowest model level and the ground
31 ! using Monin-Obukhov scaling
32 !
33 !=======================================================================
34 
35 
36 use constants_mod, only: grav, vonkarm
37 use mpp_mod, only: input_nml_file
38 use fms_mod, only: error_mesg, fatal, file_exist, &
39  check_nml_error, open_namelist_file, &
40  mpp_pe, mpp_root_pe, close_file, stdlog, &
41  write_version_number
42 
43 implicit none
44 private
45 
46 !=======================================================================
47  public monin_obukhov_init
48  public monin_obukhov_end
49  public mo_drag
50  public mo_profile
51  public mo_diff
52  public stable_mix
53 !=======================================================================
54 
55 interface mo_drag
56  module procedure mo_drag_0d, mo_drag_1d, mo_drag_2d
57 end interface
58 
59 
60 interface mo_profile
61  module procedure mo_profile_0d, mo_profile_1d, mo_profile_2d, &
63 end interface
64 
65 interface mo_diff
66  module procedure mo_diff_0d_n, mo_diff_0d_1, &
69 end interface
70 
71 interface stable_mix
72  module procedure stable_mix_0d, stable_mix_1d, &
74 end interface
75 
76 
77 !--------------------- version number ---------------------------------
78 
79 character(len=128) :: version = '$Id$'
80 character(len=128) :: tagname = '$Name$'
81 
82 !=======================================================================
83 
84 ! DEFAULT VALUES OF NAMELIST PARAMETERS:
85 
86 real :: rich_crit = 2.0
87 real :: drag_min_heat = 1.e-05
88 real :: drag_min_moist = 1.e-05
89 real :: drag_min_mom = 1.e-05
90 logical :: neutral = .false.
91 integer :: stable_option = 1
92 real :: zeta_trans = 0.5
93 logical :: new_mo_option = .false.
94 
95 
96 namelist /monin_obukhov_nml/ rich_crit, neutral, drag_min_heat, &
99 
100 !=======================================================================
101 
102 ! MODULE VARIABLES
103 
104 real, parameter :: small = 1.e-04
107 logical :: module_is_initialized = .false.
108 
109 
110 contains
111 
112 !=======================================================================
113 
114 subroutine monin_obukhov_init
116 integer :: unit, ierr, io, logunit
117 
118 !------------------- read namelist input -------------------------------
119 
120 #ifdef INTERNAL_FILE_NML
121  read (input_nml_file, nml=monin_obukhov_nml, iostat=io)
122  ierr = check_nml_error(io,"monin_obukhov_nml")
123 #else
124  if (file_exist('input.nml')) then
125  unit = open_namelist_file()
126  ierr=1; do while (ierr /= 0)
127  read (unit, nml=monin_obukhov_nml, iostat=io, end=10)
128  ierr = check_nml_error(io,'monin_obukhov_nml')
129  enddo
130  10 call close_file (unit)
131  endif
132 #endif
133 
134 !---------- output namelist to log-------------------------------------
135 
136  if ( mpp_pe() == mpp_root_pe() ) then
137  call write_version_number(version, tagname)
138  logunit = stdlog()
139  write (logunit, nml=monin_obukhov_nml)
140  endif
141 
142 !----------------------------------------------------------------------
143 
144 if(rich_crit.le.0.25) call error_mesg( &
145  'MONIN_OBUKHOV_INIT in MONIN_OBUKHOV_MOD', &
146  'rich_crit in monin_obukhov_mod must be > 0.25', fatal)
147 
148 if(drag_min_heat.le.0.0) call error_mesg( &
149  'MONIN_OBUKHOV_INIT in MONIN_OBUKHOV_MOD', &
150  'drag_min_heat in monin_obukhov_mod must be >= 0.0', fatal)
151 
152 if(drag_min_moist.le.0.0) call error_mesg( &
153  'MONIN_OBUKHOV_INIT in MONIN_OBUKHOV_MOD', &
154  'drag_min_moist in monin_obukhov_mod must be >= 0.0', fatal)
155 
156 if(drag_min_mom.le.0.0) call error_mesg( &
157  'MONIN_OBUKHOV_INIT in MONIN_OBUKHOV_MOD', &
158  'drag_min_mom in monin_obukhov_mod must be >= 0.0', fatal)
159 
160 if(stable_option < 1 .or. stable_option > 2) call error_mesg( &
161  'MONIN_OBUKHOV_INIT in MONIN_OBUKHOV_MOD', &
162  'the only allowable values of stable_option are 1 and 2', fatal)
163 
164 if(stable_option == 2 .and. zeta_trans < 0) call error_mesg( &
165  'MONIN_OBUKHOV_INIT in MONIN_OBUKHOV_MOD', &
166  'zeta_trans must be positive', fatal)
167 
168 b_stab = 1.0/rich_crit
169 r_crit = 0.95*rich_crit ! convergence can get slow if one is
170  ! close to rich_crit
171 
172 sqrt_drag_min_heat = 0.0
174 
177 
178 sqrt_drag_min_mom = 0.0
180 
181 lambda = 1.0 + (5.0 - b_stab)*zeta_trans ! used only if stable_option = 2
182 rich_trans = zeta_trans/(1.0 + 5.0*zeta_trans) ! used only if stable_option = 2
183 
184 module_is_initialized = .true.
185 
186 return
187 end subroutine monin_obukhov_init
188 
189 !=======================================================================
190 
191 subroutine monin_obukhov_end
193 module_is_initialized = .false.
194 
195 end subroutine monin_obukhov_end
196 
197 !=======================================================================
198 
199 subroutine mo_drag_1d &
200  (pt, pt0, z, z0, zt, zq, speed, drag_m, drag_t, drag_q, &
201  u_star, b_star, avail)
203 real, intent(in) , dimension(:) :: pt, pt0, z, z0, zt, zq, speed
204 real, intent(inout), dimension(:) :: drag_m, drag_t, drag_q, u_star, b_star
205 logical, intent(in), optional, dimension(:) :: avail
206 
207 logical :: lavail, avail_dummy(1)
208 integer :: n, ier
209 
210 integer, parameter :: max_iter = 20
211 real , parameter :: error=1.e-04, zeta_min=1.e-06, small=1.e-04
212 
213 ! #include "monin_obukhov_interfaces.h"
214 
215 if(.not.module_is_initialized) call error_mesg('mo_drag_1d in monin_obukhov_mod', &
216  'monin_obukhov_init has not been called', fatal)
217 
218 n = size(pt)
219 lavail = .false.
220 if(present(avail)) lavail = .true.
221 
222 
223 if(lavail) then
224  if (count(avail) .eq. 0) return
225  call monin_obukhov_drag_1d(grav, vonkarm, &
226  & error, zeta_min, max_iter, small, &
229  & n, pt, pt0, z, z0, zt, zq, speed, drag_m, drag_t, &
230  & drag_q, u_star, b_star, lavail, avail, ier)
231 else
232  call monin_obukhov_drag_1d(grav, vonkarm, &
233  & error, zeta_min, max_iter, small, &
236  & n, pt, pt0, z, z0, zt, zq, speed, drag_m, drag_t, &
237  & drag_q, u_star, b_star, lavail, avail_dummy, ier)
238 endif
239 
240 end subroutine mo_drag_1d
241 
242 
243 !=======================================================================
244 
245 subroutine mo_profile_1d(zref, zref_t, z, z0, zt, zq, u_star, b_star, q_star, &
246  del_m, del_t, del_q, avail)
248 real, intent(in) :: zref, zref_t
249 real, intent(in) , dimension(:) :: z, z0, zt, zq, u_star, b_star, q_star
250 real, intent(out), dimension(:) :: del_m, del_t, del_q
251 logical, intent(in) , optional, dimension(:) :: avail
252 
253 logical :: dummy_avail(1)
254 integer :: n, ier
255 
256 ! #include "monin_obukhov_interfaces.h"
257 
258 if(.not.module_is_initialized) call error_mesg('mo_profile_1d in monin_obukhov_mod', &
259  'monin_obukhov_init has not been called', fatal)
260 
261 n = size(z)
262 if(present(avail)) then
263 
264  if (count(avail) .eq. 0) return
265 
266  call monin_obukhov_profile_1d(vonkarm, &
268  & n, zref, zref_t, z, z0, zt, zq, u_star, b_star, q_star, &
269  & del_m, del_t, del_q, .true., avail, ier)
270 
271 else
272 
273  call monin_obukhov_profile_1d(vonkarm, &
275  & n, zref, zref_t, z, z0, zt, zq, u_star, b_star, q_star, &
276  & del_m, del_t, del_q, .false., dummy_avail, ier)
277 
278 endif
279 
280 end subroutine mo_profile_1d
281 
282 !=======================================================================
283 
284 subroutine stable_mix_3d(rich, mix)
286 real, intent(in) , dimension(:,:,:) :: rich
287 real, intent(out), dimension(:,:,:) :: mix
288 
289 integer :: n, ier
290 
291 if(.not.module_is_initialized) call error_mesg('stable_mix_3d in monin_obukhov_mod', &
292  'monin_obukhov_init has not been called', fatal)
293 
294 n = size(rich,1)*size(rich,2)*size(rich,3)
295 call monin_obukhov_stable_mix(stable_option, rich_crit, zeta_trans, &
296  & n, rich, mix, ier)
297 
298 
299 end subroutine stable_mix_3d
300 
301 !=======================================================================
302 
303 subroutine mo_diff_2d_n(z, u_star, b_star, k_m, k_h)
305 real, intent(in), dimension(:,:,:) :: z
306 real, intent(in), dimension(:,:) :: u_star, b_star
307 real, intent(out), dimension(:,:,:) :: k_m, k_h
308 
309 integer :: ni, nj, nk, ier
310 real, parameter :: ustar_min = 1.e-10
311 
312 if(.not.module_is_initialized) call error_mesg('mo_diff_2d_n in monin_obukhov_mod', &
313  'monin_obukhov_init has not been called', fatal)
314 
315 ni = size(z, 1); nj = size(z, 2); nk = size(z, 3)
317  & ustar_min, &
319  & ni, nj, nk, z, u_star, b_star, k_m, k_h, ier)
320 
321 end subroutine mo_diff_2d_n
322 
323 !=======================================================================
324 ! The following routines are used by the public interfaces above
325 !=======================================================================
326 
327 subroutine solve_zeta(rich, z, z0, zt, zq, f_m, f_t, f_q, mask)
329 real , intent(in) , dimension(:) :: rich, z, z0, zt, zq
330 logical, intent(in) , dimension(:) :: mask
331 real , intent(out), dimension(:) :: f_m, f_t, f_q
332 
333 
334 real, parameter :: error = 1.e-04
335 real, parameter :: zeta_min = 1.e-06
336 integer, parameter :: max_iter = 20
337 
338 real :: max_cor
339 integer :: iter
340 
341 real, dimension(size(rich(:))) :: &
342  d_rich, rich_1, correction, corr, z_z0, z_zt, z_zq, &
343  ln_z_z0, ln_z_zt, ln_z_zq, zeta, &
344  phi_m, phi_m_0, phi_t, phi_t_0, rzeta, &
345  zeta_0, zeta_t, zeta_q, df_m, df_t
346 
347 logical, dimension(size(rich(:))) :: mask_1
348 
349 
350 z_z0 = z/z0
351 z_zt = z/zt
352 z_zq = z/zq
353 ln_z_z0 = log(z_z0)
354 ln_z_zt = log(z_zt)
355 ln_z_zq = log(z_zq)
356 
357 corr = 0.0
358 mask_1 = mask
359 
360 ! initial guess
361 
362 where(mask_1)
363  zeta = rich*ln_z_z0*ln_z_z0/ln_z_zt
364 elsewhere
365  zeta = 0.0
366 end where
367 
368 where (mask_1 .and. rich >= 0.0)
369  zeta = zeta/(1.0 - rich/rich_crit)
370 end where
371 
372 iter_loop: do iter = 1, max_iter
373 
374  where (mask_1 .and. abs(zeta).lt.zeta_min)
375  zeta = 0.0
376  f_m = ln_z_z0
377  f_t = ln_z_zt
378  f_q = ln_z_zq
379  mask_1 = .false. ! don't do any more calculations at these pts
380  end where
381 
382  where (mask_1)
383  rzeta = 1.0/zeta
384  zeta_0 = zeta/z_z0
385  zeta_t = zeta/z_zt
386  zeta_q = zeta/z_zq
387  elsewhere
388  zeta_0 = 0.0
389  zeta_t = 0.0
390  zeta_q = 0.0
391  end where
392 
393  call mo_derivative_m(phi_m , zeta , mask_1)
394  call mo_derivative_m(phi_m_0, zeta_0, mask_1)
395  call mo_derivative_t(phi_t , zeta , mask_1)
396  call mo_derivative_t(phi_t_0, zeta_t, mask_1)
397 
398  call mo_integral_m(f_m, zeta, zeta_0, ln_z_z0, mask_1)
399  call mo_integral_tq(f_t, f_q, zeta, zeta_t, zeta_q, ln_z_zt, ln_z_zq, mask_1)
400 
401  where (mask_1)
402  df_m = (phi_m - phi_m_0)*rzeta
403  df_t = (phi_t - phi_t_0)*rzeta
404  rich_1 = zeta*f_t/(f_m*f_m)
405  d_rich = rich_1*( rzeta + df_t/f_t - 2.0 *df_m/f_m)
406  correction = (rich - rich_1)/d_rich
407  corr = min(abs(correction),abs(correction/zeta))
408  ! the criterion corr < error seems to work ok, but is a bit arbitrary
409  ! when zeta is small the tolerance is reduced
410  end where
411 
412  max_cor= maxval(corr)
413 
414  if(max_cor > error) then
415  mask_1 = mask_1 .and. (corr > error)
416  ! change the mask so computation proceeds only on non-converged points
417  where(mask_1)
418  zeta = zeta + correction
419  end where
420  cycle iter_loop
421  else
422  return
423  end if
424 
425 end do iter_loop
426 
427 call error_mesg ('solve_zeta in monin_obukhov_mod', &
428  'surface drag iteration did not converge', fatal)
429 
430 end subroutine solve_zeta
431 
432 !=======================================================================
433 
434 subroutine mo_derivative_m(phi_m, zeta, mask)
436 ! the differential similarity function for momentum
437 
438 real , intent(out), dimension(:) :: phi_m
439 real , intent(in), dimension(:) :: zeta
440 logical , intent(in), dimension(:) :: mask
441 
442 logical, dimension(size(zeta(:))) :: stable, unstable
443 real , dimension(size(zeta(:))) :: x
444 
445 stable = mask .and. zeta >= 0.0
446 unstable = mask .and. zeta < 0.0
447 
448 where (unstable)
449  x = (1 - 16.0*zeta )**(-0.5)
450  phi_m = sqrt(x) ! phi_m = (1 - 16.0*zeta)**(-0.25)
451 end where
452 
453 if(stable_option == 1) then
454 
455  where (stable)
456  phi_m = 1.0 + zeta *(5.0 + b_stab*zeta)/(1.0 + zeta)
457  end where
458 
459 else if(stable_option == 2) then
460 
461  where (stable .and. zeta < zeta_trans)
462  phi_m = 1 + 5.0*zeta
463  end where
464  where (stable .and. zeta >= zeta_trans)
465  phi_m = lambda + b_stab*zeta
466  end where
467 
468 endif
469 
470 return
471 end subroutine mo_derivative_m
472 
473 !=======================================================================
474 
475 subroutine mo_derivative_t(phi_t, zeta, mask)
477 ! the differential similarity function for buoyancy and tracers
478 
479 real , intent(out), dimension(:) :: phi_t
480 real , intent(in), dimension(:) :: zeta
481 logical , intent(in), dimension(:) :: mask
482 
483 logical, dimension(size(zeta(:))) :: stable, unstable
484 
485 stable = mask .and. zeta >= 0.0
486 unstable = mask .and. zeta < 0.0
487 
488 where (unstable)
489  phi_t = (1 - 16.0*zeta)**(-0.5)
490 end where
491 
492 if(stable_option == 1) then
493 
494  where (stable)
495  phi_t = 1.0 + zeta*(5.0 + b_stab*zeta)/(1.0 + zeta)
496  end where
497 
498 else if(stable_option == 2) then
499 
500  where (stable .and. zeta < zeta_trans)
501  phi_t = 1 + 5.0*zeta
502  end where
503  where (stable .and. zeta >= zeta_trans)
504  phi_t = lambda + b_stab*zeta
505  end where
506 
507 endif
508 
509 return
510 end subroutine mo_derivative_t
511 
512 !=======================================================================
513 
514 subroutine mo_integral_tq (psi_t, psi_q, zeta, zeta_t, zeta_q, &
515  ln_z_zt, ln_z_zq, mask)
517 ! the integral similarity function for moisture and tracers
518 
519 real , intent(out), dimension(:) :: psi_t, psi_q
520 real , intent(in), dimension(:) :: zeta, zeta_t, zeta_q, ln_z_zt, ln_z_zq
521 logical , intent(in), dimension(:) :: mask
522 
523 real, dimension(size(zeta(:))) :: x, x_t, x_q
524 
525 logical, dimension(size(zeta(:))) :: stable, unstable, &
526  weakly_stable, strongly_stable
527 
528 stable = mask .and. zeta >= 0.0
529 unstable = mask .and. zeta < 0.0
530 
531 where(unstable)
532 
533  x = sqrt(1 - 16.0*zeta)
534  x_t = sqrt(1 - 16.0*zeta_t)
535  x_q = sqrt(1 - 16.0*zeta_q)
536 
537  psi_t = ln_z_zt - 2.0*log( (1.0 + x)/(1.0 + x_t) )
538  psi_q = ln_z_zq - 2.0*log( (1.0 + x)/(1.0 + x_q) )
539 
540 end where
541 
542 if( stable_option == 1) then
543 
544  where (stable)
545 
546  psi_t = ln_z_zt + (5.0 - b_stab)*log((1.0 + zeta)/(1.0 + zeta_t)) &
547  + b_stab*(zeta - zeta_t)
548  psi_q = ln_z_zq + (5.0 - b_stab)*log((1.0 + zeta)/(1.0 + zeta_q)) &
549  + b_stab*(zeta - zeta_q)
550 
551  end where
552 
553 else if (stable_option == 2) then
554 
555  weakly_stable = stable .and. zeta <= zeta_trans
556  strongly_stable = stable .and. zeta > zeta_trans
557 
558  where (weakly_stable)
559  psi_t = ln_z_zt + 5.0*(zeta - zeta_t)
560  psi_q = ln_z_zq + 5.0*(zeta - zeta_q)
561  end where
562 
563  where(strongly_stable)
564  x = (lambda - 1.0)*log(zeta/zeta_trans) + b_stab*(zeta - zeta_trans)
565  endwhere
566 
567  where (strongly_stable .and. zeta_t <= zeta_trans)
568  psi_t = ln_z_zt + x + 5.0*(zeta_trans - zeta_t)
569  end where
570  where (strongly_stable .and. zeta_t > zeta_trans)
571  psi_t = lambda*ln_z_zt + b_stab*(zeta - zeta_t)
572  endwhere
573 
574  where (strongly_stable .and. zeta_q <= zeta_trans)
575  psi_q = ln_z_zq + x + 5.0*(zeta_trans - zeta_q)
576  end where
577  where (strongly_stable .and. zeta_q > zeta_trans)
578  psi_q = lambda*ln_z_zq + b_stab*(zeta - zeta_q)
579  endwhere
580 
581 end if
582 
583 return
584 end subroutine mo_integral_tq
585 
586 !=======================================================================
587 
588 subroutine mo_integral_m (psi_m, zeta, zeta_0, ln_z_z0, mask)
590 ! the integral similarity function for momentum
591 
592 real , intent(out), dimension(:) :: psi_m
593 real , intent(in), dimension(:) :: zeta, zeta_0, ln_z_z0
594 logical , intent(in), dimension(:) :: mask
595 
596 real, dimension(size(zeta(:))) :: x, x_0, x1, x1_0, num, denom, y
597 
598 logical, dimension(size(zeta(:))) :: stable, unstable, &
599  weakly_stable, strongly_stable
600 
601 stable = mask .and. zeta >= 0.0
602 unstable = mask .and. zeta < 0.0
603 
604 where(unstable)
605 
606  x = sqrt(1 - 16.0*zeta)
607  x_0 = sqrt(1 - 16.0*zeta_0)
608 
609  x = sqrt(x)
610  x_0 = sqrt(x_0)
611 
612  x1 = 1.0 + x
613  x1_0 = 1.0 + x_0
614 
615  num = x1*x1*(1.0 + x*x)
616  denom = x1_0*x1_0*(1.0 + x_0*x_0)
617  y = atan(x) - atan(x_0)
618  psi_m = ln_z_z0 - log(num/denom) + 2*y
619 
620 end where
621 
622 if( stable_option == 1) then
623 
624  where (stable)
625  psi_m = ln_z_z0 + (5.0 - b_stab)*log((1.0 + zeta)/(1.0 + zeta_0)) &
626  + b_stab*(zeta - zeta_0)
627  end where
628 
629 else if (stable_option == 2) then
630 
631  weakly_stable = stable .and. zeta <= zeta_trans
632  strongly_stable = stable .and. zeta > zeta_trans
633 
634  where (weakly_stable)
635  psi_m = ln_z_z0 + 5.0*(zeta - zeta_0)
636  end where
637 
638  where(strongly_stable)
639  x = (lambda - 1.0)*log(zeta/zeta_trans) + b_stab*(zeta - zeta_trans)
640  endwhere
641 
642  where (strongly_stable .and. zeta_0 <= zeta_trans)
643  psi_m = ln_z_z0 + x + 5.0*(zeta_trans - zeta_0)
644  end where
645  where (strongly_stable .and. zeta_0 > zeta_trans)
646  psi_m = lambda*ln_z_z0 + b_stab*(zeta - zeta_0)
647  endwhere
648 
649 end if
650 
651 return
652 end subroutine mo_integral_m
653 
654 
655 !=======================================================================
656 ! The following routines allow the public interfaces to be used
657 ! with different dimensions of the input and output
658 !
659 !=======================================================================
660 
661 
662 subroutine mo_drag_2d &
663  (pt, pt0, z, z0, zt, zq, speed, drag_m, drag_t, drag_q, u_star, b_star)
665 real, intent(in) , dimension(:,:) :: z, speed, pt, pt0, z0, zt, zq
666 real, intent(out) , dimension(:,:) :: drag_m, drag_t, drag_q
667 real, intent(inout), dimension(:,:) :: u_star, b_star
668 
669 integer :: j
670 
671 do j = 1, size(pt,2)
672  call mo_drag_1d (pt(:,j), pt0(:,j), z(:,j), z0(:,j), zt(:,j), zq(:,j), &
673  speed(:,j), drag_m(:,j), drag_t(:,j), drag_q(:,j), &
674  u_star(:,j), b_star(:,j))
675 end do
676 
677 
678 return
679 end subroutine mo_drag_2d
680 
681 !=======================================================================
682 subroutine mo_drag_0d &
683  (pt, pt0, z, z0, zt, zq, speed, drag_m, drag_t, drag_q, u_star, b_star)
685 real, intent(in) :: z, speed, pt, pt0, z0, zt, zq
686 real, intent(out) :: drag_m, drag_t, drag_q, u_star, b_star
687 
688 real, dimension(1) :: pt_1, pt0_1, z_1, z0_1, zt_1, zq_1, speed_1, &
689  drag_m_1, drag_t_1, drag_q_1, u_star_1, b_star_1
690 
691 pt_1(1) = pt
692 pt0_1(1) = pt0
693 z_1(1) = z
694 z0_1(1) = z0
695 zt_1(1) = zt
696 zq_1(1) = zq
697 speed_1(1) = speed
698 
699 call mo_drag_1d (pt_1, pt0_1, z_1, z0_1, zt_1, zq_1, speed_1, &
700  drag_m_1, drag_t_1, drag_q_1, u_star_1, b_star_1)
701 
702 drag_m = drag_m_1(1)
703 drag_t = drag_t_1(1)
704 drag_q = drag_q_1(1)
705 u_star = u_star_1(1)
706 b_star = b_star_1(1)
707 
708 return
709 end subroutine mo_drag_0d
710 !=======================================================================
711 
712 subroutine mo_profile_2d(zref, zref_t, z, z0, zt, zq, u_star, b_star, q_star, &
713  del_m, del_h, del_q)
715 real, intent(in) :: zref, zref_t
716 real, intent(in) , dimension(:,:) :: z, z0, zt, zq, u_star, b_star, q_star
717 real, intent(out), dimension(:,:) :: del_m, del_h, del_q
718 
719 integer :: j
720 
721 do j = 1, size(z,2)
722  call mo_profile_1d (zref, zref_t, z(:,j), z0(:,j), zt(:,j), &
723  zq(:,j), u_star(:,j), b_star(:,j), q_star(:,j), &
724  del_m(:,j), del_h(:,j), del_q(:,j))
725 enddo
726 
727 return
728 end subroutine mo_profile_2d
729 
730 !=======================================================================
731 
732 subroutine mo_profile_0d(zref, zref_t, z, z0, zt, zq, u_star, b_star, q_star, &
733  del_m, del_h, del_q)
735 real, intent(in) :: zref, zref_t
736 real, intent(in) :: z, z0, zt, zq, u_star, b_star, q_star
737 real, intent(out) :: del_m, del_h, del_q
738 
739 real, dimension(1) :: z_1, z0_1, zt_1, zq_1, u_star_1, b_star_1, q_star_1, &
740  del_m_1, del_h_1, del_q_1
741 
742 z_1(1) = z
743 z0_1(1) = z0
744 zt_1(1) = zt
745 zq_1(1) = zq
746 u_star_1(1) = u_star
747 b_star_1(1) = b_star
748 q_star_1(1) = q_star
749 
750 call mo_profile_1d (zref, zref_t, z_1, z0_1, zt_1, zq_1, &
751  u_star_1, b_star_1, q_star_1, &
752  del_m_1, del_h_1, del_q_1)
753 
754 del_m = del_m_1(1)
755 del_h = del_h_1(1)
756 del_q = del_q_1(1)
757 
758 
759 return
760 end subroutine mo_profile_0d
761 
762 !=======================================================================
763 
764 subroutine mo_profile_1d_n(zref, z, z0, zt, zq, u_star, b_star, q_star, &
765  del_m, del_t, del_q, avail)
767 real, intent(in), dimension(:) :: zref
768 real, intent(in) , dimension(:) :: z, z0, zt, zq, u_star, b_star, q_star
769 real, intent(out), dimension(:,:) :: del_m, del_t, del_q
770 logical, intent(in) , optional, dimension(:) :: avail
771 
772 integer :: k
773 
774 do k = 1, size(zref(:))
775  if(present(avail)) then
776  call mo_profile_1d (zref(k), zref(k), z, z0, zt, zq, &
777  u_star, b_star, q_star, del_m(:,k), del_t(:,k), del_q(:,k), avail)
778  else
779  call mo_profile_1d (zref(k), zref(k), z, z0, zt, zq, &
780  u_star, b_star, q_star, del_m(:,k), del_t(:,k), del_q(:,k))
781  endif
782 enddo
783 
784 return
785 end subroutine mo_profile_1d_n
786 
787 !=======================================================================
788 
789 subroutine mo_profile_0d_n(zref, z, z0, zt, zq, u_star, b_star, q_star, &
790  del_m, del_t, del_q)
792 real, intent(in), dimension(:) :: zref
793 real, intent(in) :: z, z0, zt, zq, u_star, b_star, q_star
794 real, intent(out), dimension(:) :: del_m, del_t, del_q
795 
796 integer :: k
797 
798 do k = 1, size(zref(:))
799  call mo_profile_0d (zref(k), zref(k), z, z0, zt, zq, &
800  u_star, b_star, q_star, del_m(k), del_t(k), del_q(k))
801 enddo
802 
803 return
804 end subroutine mo_profile_0d_n
805 
806 !=======================================================================
807 
808 subroutine mo_profile_2d_n(zref, z, z0, zt, zq, u_star, b_star, q_star, &
809  del_m, del_t, del_q)
811 real, intent(in), dimension(:) :: zref
812 real, intent(in), dimension(:,:) :: z, z0, zt, zq, u_star, b_star, q_star
813 real, intent(out), dimension(:,:,:) :: del_m, del_t, del_q
814 
815 integer :: k
816 
817 do k = 1, size(zref(:))
818  call mo_profile_2d (zref(k), zref(k), z, z0, zt, zq, &
819  u_star, b_star, q_star, del_m(:,:,k), del_t(:,:,k), del_q(:,:,k))
820 enddo
821 
822 return
823 end subroutine mo_profile_2d_n
824 
825 !=======================================================================
826 
827 subroutine mo_diff_2d_1(z, u_star, b_star, k_m, k_h)
829 real, intent(in), dimension(:,:) :: z, u_star, b_star
830 real, intent(out), dimension(:,:) :: k_m, k_h
831 
832 real , dimension(size(z,1),size(z,2),1) :: z_n, k_m_n, k_h_n
833 
834 z_n(:,:,1) = z
835 
836 call mo_diff_2d_n(z_n, u_star, b_star, k_m_n, k_h_n)
837 
838 k_m = k_m_n(:,:,1)
839 k_h = k_h_n(:,:,1)
840 
841 return
842 end subroutine mo_diff_2d_1
843 
844 
845 !=======================================================================
846 
847 subroutine mo_diff_1d_1(z, u_star, b_star, k_m, k_h)
849 real, intent(in), dimension(:) :: z, u_star, b_star
850 real, intent(out), dimension(:) :: k_m, k_h
851 
852 real, dimension(size(z),1,1) :: z_n, k_m_n, k_h_n
853 real, dimension(size(z),1) :: u_star_n, b_star_n
854 
855 z_n(:,1,1) = z
856 u_star_n(:,1) = u_star
857 b_star_n(:,1) = b_star
858 
859 call mo_diff_2d_n(z_n, u_star_n, b_star_n, k_m_n, k_h_n)
860 
861 k_m = k_m_n(:,1,1)
862 k_h = k_h_n(:,1,1)
863 
864 return
865 end subroutine mo_diff_1d_1
866 
867 !=======================================================================
868 
869 subroutine mo_diff_1d_n(z, u_star, b_star, k_m, k_h)
871 real, intent(in), dimension(:,:) :: z
872 real, intent(in), dimension(:) :: u_star, b_star
873 real, intent(out), dimension(:,:) :: k_m, k_h
874 
875 real, dimension(size(z,1),1) :: u_star2, b_star2
876 real, dimension(size(z,1),1, size(z,2)) :: z2, k_m2, k_h2
877 
878 integer :: n
879 
880 do n = 1, size(z,2)
881  z2(:,1,n) = z(:,n)
882 enddo
883 u_star2(:,1) = u_star
884 b_star2(:,1) = b_star
885 
886 call mo_diff_2d_n(z2, u_star2, b_star2, k_m2, k_h2)
887 
888 do n = 1, size(z,2)
889  k_m(:,n) = k_m2(:,1,n)
890  k_h(:,n) = k_h2(:,1,n)
891 enddo
892 
893 return
894 end subroutine mo_diff_1d_n
895 
896 !=======================================================================
897 
898 subroutine mo_diff_0d_1(z, u_star, b_star, k_m, k_h)
900 real, intent(in) :: z, u_star, b_star
901 real, intent(out) :: k_m, k_h
902 
903 integer :: ni, nj, nk, ier
904 real, parameter :: ustar_min = 1.e-10
905 
906 if(.not.module_is_initialized) call error_mesg('mo_diff_0d_1 in monin_obukhov_mod', &
907  'monin_obukhov_init has not been called', fatal)
908 
909 ni = 1; nj = 1; nk = 1
911  & ustar_min, &
913  & ni, nj, nk, z, u_star, b_star, k_m, k_h, ier)
914 
915 end subroutine mo_diff_0d_1
916 
917 !=======================================================================
918 
919 subroutine mo_diff_0d_n(z, u_star, b_star, k_m, k_h)
921 real, intent(in), dimension(:) :: z
922 real, intent(in) :: u_star, b_star
923 real, intent(out), dimension(:) :: k_m, k_h
924 
925 integer :: ni, nj, nk, ier
926 real, parameter :: ustar_min = 1.e-10
927 
928 if(.not.module_is_initialized) call error_mesg('mo_diff_0d_n in monin_obukhov_mod', &
929  'monin_obukhov_init has not been called', fatal)
930 
931 ni = 1; nj = 1; nk = size(z(:))
933  & ustar_min, &
935  & ni, nj, nk, z, u_star, b_star, k_m, k_h, ier)
936 
937 end subroutine mo_diff_0d_n
938 
939 !=======================================================================
940 
941 subroutine stable_mix_2d(rich, mix)
943 real, intent(in) , dimension(:,:) :: rich
944 real, intent(out), dimension(:,:) :: mix
945 
946 real, dimension(size(rich,1),size(rich,2),1) :: rich_3d, mix_3d
947 
948 rich_3d(:,:,1) = rich
949 
950 call stable_mix_3d(rich_3d, mix_3d)
951 
952 mix = mix_3d(:,:,1)
953 
954 return
955 end subroutine stable_mix_2d
956 
957 
958 !=======================================================================
959 
960 subroutine stable_mix_1d(rich, mix)
962 real, intent(in) , dimension(:) :: rich
963 real, intent(out), dimension(:) :: mix
964 
965 real, dimension(size(rich),1,1) :: rich_3d, mix_3d
966 
967 rich_3d(:,1,1) = rich
968 
969 call stable_mix_3d(rich_3d, mix_3d)
970 
971 mix = mix_3d(:,1,1)
972 
973 return
974 end subroutine stable_mix_1d
975 
976 !=======================================================================
977 
978 subroutine stable_mix_0d(rich, mix)
980 real, intent(in) :: rich
981 real, intent(out) :: mix
982 
983 real, dimension(1,1,1) :: rich_3d, mix_3d
984 
985 rich_3d(1,1,1) = rich
986 
987 call stable_mix_3d(rich_3d, mix_3d)
988 
989 mix = mix_3d(1,1,1)
990 
991 return
992 end subroutine stable_mix_0d
993 !=======================================================================
994 
995 end module monin_obukhov_mod
996 
Definition: fms.F90:20
real, parameter small
subroutine mo_drag_0d(pt, pt0, z, z0, zt, zq, speed, drag_m, drag_t, drag_q, u_star, b_star)
subroutine mo_profile_0d_n(zref, z, z0, zt, zq, u_star, b_star, q_star, del_m, del_t, del_q)
real, parameter, public vonkarm
Von Karman constant [dimensionless].
Definition: constants.F90:131
logical module_is_initialized
subroutine mo_integral_tq(psi_t, psi_q, zeta, zeta_t, zeta_q, ln_z_zt, ln_z_zq, mask)
subroutine mo_drag_1d(pt, pt0, z, z0, zt, zq, speed, drag_m, drag_t, drag_q, u_star, b_star, avail)
subroutine stable_mix_1d(rich, mix)
Definition: mpp.F90:39
integer function, public check_nml_error(IOSTAT, NML_NAME)
Definition: fms.F90:658
character(len=128) tagname
subroutine mo_profile_2d(zref, zref_t, z, z0, zt, zq, u_star, b_star, q_star, del_m, del_h, del_q)
subroutine stable_mix_2d(rich, mix)
subroutine, public monin_obukhov_end
*f90 *! $Id interface _PURE subroutine monin_obukhov_diff(vonkarm, &&ustar_min, &&neutral, stable_option, new_mo_option, rich_crit, zeta_trans, &!miz &ni, nj, nk, z, u_star, b_star, k_m, k_h, ier) real
subroutine mo_profile_1d_n(zref, z, z0, zt, zq, u_star, b_star, q_star, del_m, del_t, del_q, avail)
character(len=input_str_length), dimension(:), allocatable, target, public input_nml_file
Definition: mpp.F90:1378
subroutine stable_mix_3d(rich, mix)
subroutine mo_profile_1d(zref, zref_t, z, z0, zt, zq, u_star, b_star, q_star, del_m, del_t, del_q, avail)
subroutine mo_drag_2d(pt, pt0, z, z0, zt, zq, speed, drag_m, drag_t, drag_q, u_star, b_star)
subroutine mo_derivative_t(phi_t, zeta, mask)
subroutine mo_diff_2d_1(z, u_star, b_star, k_m, k_h)
subroutine stable_mix_0d(rich, mix)
subroutine solve_zeta(rich, z, z0, zt, zq, f_m, f_t, f_q, mask)
subroutine mo_profile_2d_n(zref, z, z0, zt, zq, u_star, b_star, q_star, del_m, del_t, del_q)
subroutine mo_diff_1d_n(z, u_star, b_star, k_m, k_h)
real, parameter, public grav
Acceleration due to gravity [m/s^2].
Definition: constants.F90:76
character(len=128) version
subroutine mo_diff_1d_1(z, u_star, b_star, k_m, k_h)
subroutine, public monin_obukhov_init
subroutine mo_integral_m(psi_m, zeta, zeta_0, ln_z_z0, mask)
subroutine mo_profile_0d(zref, zref_t, z, z0, zt, zq, u_star, b_star, q_star, del_m, del_h, del_q)
subroutine mo_diff_0d_n(z, u_star, b_star, k_m, k_h)
subroutine mo_derivative_m(phi_m, zeta, mask)
#define min(a, b)
Definition: mosaic_util.h:32
subroutine mo_diff_0d_1(z, u_star, b_star, k_m, k_h)
subroutine, public error_mesg(routine, message, level)
Definition: fms.F90:529
subroutine mo_diff_2d_n(z, u_star, b_star, k_m, k_h)