FMS 2025.01-dev
Flexible Modeling System
Loading...
Searching...
No Matches
monin_obukhov.inc
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
20
21subroutine mo_drag_1d_ &
22 (pt, pt0, z, z0, zt, zq, speed, drag_m, drag_t, drag_q, &
23 u_star, b_star, avail)
24
25real(kind=fms_mo_kind_), intent(in) , dimension(:) :: pt, pt0, z, z0, zt, zq, speed
26real(kind=fms_mo_kind_), intent(inout), dimension(:) :: drag_m, drag_t, drag_q, u_star, b_star
27logical, intent(in), optional, dimension(:) :: avail
28
29logical :: lavail, avail_dummy(1)
30integer :: n, ier
31
32integer, parameter :: max_iter = 20
33integer, parameter :: lkind = fms_mo_kind_
34real(kind=fms_mo_kind_), parameter :: error = 1.0e-04_lkind, &
35 zeta_min = 1.0e-06_lkind, &
36 small = 1.0e-04_lkind
37
38! #include "monin_obukhov_interfaces.h"
39
40if(.not.module_is_initialized) call error_mesg('mo_drag_1d in monin_obukhov_mod', &
41 'monin_obukhov_init has not been called', fatal)
42
43n = size(pt)
44lavail = .false.
45if(present(avail)) lavail = .true.
46
47
48if(lavail) then
49 if (count(avail) .eq. 0) return
50 call monin_obukhov_drag_1d(real(grav, fms_mo_kind_), real(vonkarm, fms_mo_kind_), &
51 & error, zeta_min, max_iter, real(small, fms_mo_kind_), neutral, stable_option, &
52 & new_mo_option, real(rich_crit, fms_mo_kind_), real(zeta_trans, fms_mo_kind_), &!miz
53 & real(drag_min_heat, fms_mo_kind_), real(drag_min_moist, fms_mo_kind_), &
54 & real(drag_min_mom, fms_mo_kind_), n, pt, pt0, z, z0, zt, &
55 & zq, speed, drag_m, drag_t, drag_q, u_star, b_star, lavail, avail, ier)
56else
57call monin_obukhov_drag_1d(real(grav, fms_mo_kind_), real(vonkarm, fms_mo_kind_), &
58 & error, zeta_min, max_iter, real(small, fms_mo_kind_), neutral, stable_option, &
59 & new_mo_option, real(rich_crit, fms_mo_kind_), real(zeta_trans, fms_mo_kind_), &!miz
60 & real(drag_min_heat, fms_mo_kind_), real(drag_min_moist, fms_mo_kind_), &
61 & real(drag_min_mom, fms_mo_kind_), n, pt, pt0, z, z0, zt, zq, speed, drag_m, drag_t, &
62 & drag_q, u_star, b_star, lavail, avail_dummy, ier)
63endif
64
65end subroutine mo_drag_1d_
66
67
68!=======================================================================
69
70subroutine mo_profile_1d_(zref, zref_t, z, z0, zt, zq, u_star, b_star, q_star, &
71 del_m, del_t, del_q, avail)
72
73real(kind=fms_mo_kind_), intent(in) :: zref, zref_t
74real(kind=fms_mo_kind_), intent(in) , dimension(:) :: z, z0, zt, zq, u_star, b_star, q_star
75real(kind=fms_mo_kind_), intent(out), dimension(:) :: del_m, del_t, del_q
76logical, intent(in) , optional, dimension(:) :: avail
77
78logical :: dummy_avail(1)
79integer :: n, ier
80
81! #include "monin_obukhov_interfaces.h"
82
83if(.not.module_is_initialized) call error_mesg('mo_profile_1d in monin_obukhov_mod', &
84 'monin_obukhov_init has not been called', fatal)
85
86n = size(z)
87if(present(avail)) then
88
89 if (count(avail) .eq. 0) return
90
91 call monin_obukhov_profile_1d(real(vonkarm, fms_mo_kind_), &
92 & neutral, stable_option, new_mo_option, real(rich_crit, fms_mo_kind_), &
93 & real(zeta_trans, fms_mo_kind_), n, zref, zref_t, z, z0, zt, zq, u_star, &
94 & b_star, q_star, del_m, del_t, del_q, .true., avail, ier)
95
96else
97
98 call monin_obukhov_profile_1d(real(vonkarm, fms_mo_kind_), &
99 & neutral, stable_option, new_mo_option, real(rich_crit, fms_mo_kind_), &
100 & real(zeta_trans, fms_mo_kind_), n, zref, zref_t, z, z0, zt, zq, u_star, &
101 & b_star, q_star, del_m, del_t, del_q, .false., dummy_avail, ier)
102
103endif
104
105end subroutine mo_profile_1d_
106
107!=======================================================================
108
109subroutine stable_mix_3d_(rich, mix)
110
111real(kind=fms_mo_kind_), intent(in) , dimension(:,:,:) :: rich
112real(kind=fms_mo_kind_), intent(out), dimension(:,:,:) :: mix
113integer :: n2 !< Size of dimension 2 of mix and rich
114integer :: n3 !< Size of dimension 3 of mix and rich
115integer :: i, j !< Loop indices
116
117n2 = size(mix, 2)
118n3 = size(mix, 3)
119
120do j=1, n3
121 do i=1, n2
122 call stable_mix(rich(:, i, j), mix(:, i, j))
123 enddo
124enddo
125
126end subroutine stable_mix_3d_
127
128!=======================================================================
129
130subroutine mo_diff_2d_n_(z, u_star, b_star, k_m, k_h)
131
132real(kind=fms_mo_kind_), intent(in), dimension(:,:,:) :: z
133real(kind=fms_mo_kind_), intent(in), dimension(:,:) :: u_star, b_star
134real(kind=fms_mo_kind_), intent(out), dimension(:,:,:) :: k_m, k_h
135
136integer :: ni, nj, nk, ier
137integer, parameter :: lkind = fms_mo_kind_
138real(kind=fms_mo_kind_), parameter :: ustar_min = 1.0e-10_lkind
139
140if(.not.module_is_initialized) call error_mesg('mo_diff_2d_n in monin_obukhov_mod', &
141 'monin_obukhov_init has not been called', fatal)
142
143ni = size(z, 1); nj = size(z, 2); nk = size(z, 3)
144call monin_obukhov_diff(real(vonkarm, fms_mo_kind_), ustar_min, neutral, &
145 & stable_option, new_mo_option, real(rich_crit, fms_mo_kind_), &
146 & real(zeta_trans, fms_mo_kind_), ni, nj, nk, z, u_star, b_star, &
147 & k_m, k_h, ier)
148
149end subroutine mo_diff_2d_n_
150
151!=======================================================================
152! The following routines are used by the public interfaces above
153!=======================================================================
154
155subroutine solve_zeta_(rich, z, z0, zt, zq, f_m, f_t, f_q, mask)
156
157real(kind=fms_mo_kind_), intent(in) , dimension(:) :: rich, z, z0, zt, zq
158logical, intent(in) , dimension(:) :: mask
159real(kind=fms_mo_kind_), intent(out), dimension(:) :: f_m, f_t, f_q
160
161integer, parameter :: lkind = fms_mo_kind_
162real(kind=fms_mo_kind_), parameter :: error = 1.0e-04_lkind
163real(kind=fms_mo_kind_), parameter :: zeta_min = 1.0e-06_lkind
164integer, parameter :: max_iter = 20
165
166real(kind=fms_mo_kind_) :: max_cor
167integer :: iter
168
169real(kind=fms_mo_kind_), dimension(size(rich(:))) :: &
170 d_rich, rich_1, correction, corr, z_z0, z_zt, z_zq, &
171 ln_z_z0, ln_z_zt, ln_z_zq, zeta, &
172 phi_m, phi_m_0, phi_t, phi_t_0, rzeta, &
173 zeta_0, zeta_t, zeta_q, df_m, df_t
174
175logical, dimension(size(rich(:))) :: mask_1
176
177z_z0 = z/z0
178z_zt = z/zt
179z_zq = z/zq
180ln_z_z0 = log(z_z0)
181ln_z_zt = log(z_zt)
182ln_z_zq = log(z_zq)
183
184corr = 0.0_lkind
185mask_1 = mask
186
187! initial guess
188
189where(mask_1)
190 zeta = rich*ln_z_z0*ln_z_z0/ln_z_zt
191elsewhere
192 zeta = 0.0_lkind
193end where
194
195where (mask_1 .and. rich >= 0.0_lkind)
196 zeta = zeta/(1.0_lkind - rich/real(rich_crit, fms_mo_kind_))
197end where
198
199iter_loop: do iter = 1, max_iter
200
201 where (mask_1 .and. abs(zeta).lt.zeta_min)
202 zeta = 0.0_lkind
203 f_m = ln_z_z0
204 f_t = ln_z_zt
205 f_q = ln_z_zq
206 mask_1 = .false. ! don't do any more calculations at these pts
207 end where
208
209 where (mask_1)
210 rzeta = 1.0_lkind/zeta
211 zeta_0 = zeta/z_z0
212 zeta_t = zeta/z_zt
213 zeta_q = zeta/z_zq
214 elsewhere
215 zeta_0 = 0.0_lkind
216 zeta_t = 0.0_lkind
217 zeta_q = 0.0_lkind
218 end where
219
220 call mo_derivative_m(phi_m , zeta , mask_1)
221 call mo_derivative_m(phi_m_0, zeta_0, mask_1)
222 call mo_derivative_t(phi_t , zeta , mask_1)
223 call mo_derivative_t(phi_t_0, zeta_t, mask_1)
224
225 call mo_integral_m(f_m, zeta, zeta_0, ln_z_z0, mask_1)
226 call mo_integral_tq(f_t, f_q, zeta, zeta_t, zeta_q, ln_z_zt, ln_z_zq, mask_1)
227
228 where (mask_1)
229 df_m = (phi_m - phi_m_0)*rzeta
230 df_t = (phi_t - phi_t_0)*rzeta
231 rich_1 = zeta*f_t/(f_m*f_m)
232 d_rich = rich_1*( rzeta + df_t/f_t - 2.0_lkind *df_m/f_m)
233 correction = (rich - rich_1)/d_rich
234 corr = min(abs(correction),abs(correction/zeta))
235 ! the criterion corr < error seems to work ok, but is a bit arbitrary
236 ! when zeta is small the tolerance is reduced
237 end where
238
239 max_cor = maxval(corr)
240
241 if(max_cor > error) then
242 mask_1 = mask_1 .and. (corr > error)
243 ! change the mask so computation proceeds only on non-converged points
244 where(mask_1)
245 zeta = zeta + correction
246 end where
247 cycle iter_loop
248 else
249 return
250 end if
251
252end do iter_loop
253
254call error_mesg ('solve_zeta in monin_obukhov_mod', &
255 'surface drag iteration did not converge', fatal)
256
257end subroutine solve_zeta_
258
259!=======================================================================
260
261subroutine mo_derivative_m_(phi_m, zeta, mask)
262
263! the differential similarity function for momentum
264
265real(kind=fms_mo_kind_), intent(out), dimension(:) :: phi_m
266real(kind=fms_mo_kind_), intent(in), dimension(:) :: zeta
267logical, intent(in), dimension(:) :: mask
268
269logical, dimension(size(zeta(:))) :: stable, unstable
270real(kind=fms_mo_kind_), dimension(size(zeta(:))) :: x
271integer, parameter :: lkind = fms_mo_kind_
272
273stable = mask .and. zeta >= 0.0_lkind
274unstable = mask .and. zeta < 0.0_lkind
275
276where (unstable)
277 x = (1.0_lkind - 16.0_lkind*zeta )**(-0.5_lkind)
278 phi_m = sqrt(x) ! phi_m = (1 - 16.0*zeta)**(-0.25)
279end where
280
281if(stable_option == 1) then
282
283 where (stable)
284 phi_m = 1.0_lkind + zeta*(5.0_lkind + real(b_stab, fms_mo_kind_) &
285 *zeta)/(1.0_lkind + zeta)
286 end where
287
288else if(stable_option == 2) then
289
290 where (stable .and. zeta < real(zeta_trans,fms_mo_kind_))
291 phi_m = 1.0_lkind + 5.0_lkind*zeta
292 end where
293 where (stable .and. zeta >= real(zeta_trans,fms_mo_kind_))
294 phi_m = real(lambda, fms_mo_kind_) + real(b_stab, fms_mo_kind_)*zeta
295 end where
296
297endif
298
299return
300end subroutine mo_derivative_m_
301
302!=======================================================================
303
304subroutine mo_derivative_t_(phi_t, zeta, mask)
305
306! the differential similarity function for buoyancy and tracers
307
308real(kind=fms_mo_kind_), intent(out), dimension(:) :: phi_t
309real(kind=fms_mo_kind_), intent(in), dimension(:) :: zeta
310logical , intent(in), dimension(:) :: mask
311
312logical, dimension(size(zeta(:))) :: stable, unstable
313integer, parameter :: lkind = fms_mo_kind_
314
315stable = mask .and. zeta >= 0.0_lkind
316unstable = mask .and. zeta < 0.0_lkind
317
318where (unstable)
319 phi_t = (1.0_lkind - 16.0_lkind*zeta)**(-0.5_lkind)
320end where
321
322if(stable_option == 1) then
323
324 where (stable)
325 phi_t = 1.0_lkind + zeta * (5.0_lkind + real(b_stab, fms_mo_kind_)&
326 * zeta)/(1.0_lkind + zeta)
327 end where
328
329else if(stable_option == 2) then
330
331 where (stable .and. zeta < real(zeta_trans,fms_mo_kind_))
332 phi_t = 1.0_lkind + 5.0_lkind*zeta
333 end where
334 where (stable .and. zeta >= real(zeta_trans,fms_mo_kind_))
335 phi_t = real(lambda, fms_mo_kind_) + real(b_stab, fms_mo_kind_)*zeta
336 end where
337
338endif
339
340return
341end subroutine mo_derivative_t_
342
343!=======================================================================
344
345subroutine mo_integral_tq_ (psi_t, psi_q, zeta, zeta_t, zeta_q, &
346 ln_z_zt, ln_z_zq, mask)
347
348! the integral similarity function for moisture and tracers
349
350real(kind=fms_mo_kind_), intent(out), dimension(:) :: psi_t, psi_q
351real(kind=fms_mo_kind_), intent(in), dimension(:) :: zeta, zeta_t, zeta_q, ln_z_zt, ln_z_zq
352logical , intent(in), dimension(:) :: mask
353
354real(kind=fms_mo_kind_), dimension(size(zeta(:))) :: x, x_t, x_q
355
356logical, dimension(size(zeta(:))) :: stable, unstable, &
357 weakly_stable, strongly_stable
358integer, parameter :: lkind = fms_mo_kind_
359
360stable = mask .and. zeta >= 0.0_lkind
361unstable = mask .and. zeta < 0.0_lkind
362
363where(unstable)
364
365 x = sqrt(1.0_lkind - 16.0_lkind*zeta)
366 x_t = sqrt(1.0_lkind - 16.0_lkind*zeta_t)
367 x_q = sqrt(1.0_lkind - 16.0_lkind*zeta_q)
368
369 psi_t = ln_z_zt - 2.0_lkind*log( (1.0_lkind + x)/(1.0_lkind + x_t) )
370 psi_q = ln_z_zq - 2.0_lkind*log( (1.0_lkind + x)/(1.0_lkind + x_q) )
371
372end where
373
374if( stable_option == 1) then
375
376 where (stable)
377
378 psi_t = ln_z_zt + (5.0_lkind - real(b_stab, fms_mo_kind_)) &
379 *log((1.0_lkind + zeta)/(1.0_lkind + zeta_t)) &
380 + real(b_stab, fms_mo_kind_)*(zeta - zeta_t)
381 psi_q = ln_z_zq + (5.0_lkind - real(b_stab, fms_mo_kind_)) &
382 *log((1.0_lkind + zeta)/(1.0_lkind + zeta_q)) &
383 + real(b_stab, fms_mo_kind_)*(zeta - zeta_q)
384
385 end where
386
387else if (stable_option == 2) then
388
389 weakly_stable = stable .and. zeta <= real(zeta_trans,fms_mo_kind_)
390 strongly_stable = stable .and. zeta > real(zeta_trans,fms_mo_kind_)
391
392 where (weakly_stable)
393 psi_t = ln_z_zt + 5.0_lkind*(zeta - zeta_t)
394 psi_q = ln_z_zq + 5.0_lkind*(zeta - zeta_q)
395 end where
396
397 where(strongly_stable)
398 x = (real(lambda, fms_mo_kind_) - 1.0_lkind)*log(zeta/real(zeta_trans, fms_mo_kind_)) + &
399 real(b_stab, FMS_MO_KIND_)*(zeta - real(zeta_trans, FMS_MO_KIND_))
400 endwhere
401
402 where (strongly_stable .and. zeta_t <= real(zeta_trans,fms_mo_kind_))
403 psi_t = ln_z_zt + x + 5.0_lkind * (real(zeta_trans, fms_mo_kind_) - zeta_t)
404 end where
405
406 where (strongly_stable .and. zeta_t > real(zeta_trans,fms_mo_kind_))
407 psi_t = real(lambda, fms_mo_kind_)* ln_z_zt &
408 + real(b_stab, fms_mo_kind_)*(zeta - zeta_t)
409 endwhere
410
411 where (strongly_stable .and. zeta_q <= real(zeta_trans,fms_mo_kind_))
412 psi_q = ln_z_zq + x + 5.0_lkind &
413 *(real(zeta_trans, fms_mo_kind_) - zeta_q)
414 end where
415
416 where (strongly_stable .and. zeta_q > real(zeta_trans,fms_mo_kind_))
417 psi_q = real(lambda, fms_mo_kind_)*ln_z_zq + real(b_stab, fms_mo_kind_) &
418 * (zeta - zeta_q)
419 endwhere
420
421end if
422
423return
424end subroutine mo_integral_tq_
425
426!=======================================================================
427
428subroutine mo_integral_m_ (psi_m, zeta, zeta_0, ln_z_z0, mask)
429
430! the integral similarity function for momentum
431
432real(kind=fms_mo_kind_), intent(out), dimension(:) :: psi_m
433real(kind=fms_mo_kind_), intent(in), dimension(:) :: zeta, zeta_0, ln_z_z0
434logical, intent(in), dimension(:) :: mask
435
436real(kind=fms_mo_kind_), dimension(size(zeta(:))) :: x, x_0, x1, x1_0, num, denom, y
437
438logical, dimension(size(zeta(:))) :: stable, unstable, &
439 weakly_stable, strongly_stable
440integer, parameter :: lkind = fms_mo_kind_
441
442stable = mask .and. zeta >= 0.0_lkind
443unstable = mask .and. zeta < 0.0_lkind
444
445where(unstable)
446
447 x = sqrt(1.0_lkind - 16.0_lkind*zeta)
448 x_0 = sqrt(1.0_lkind - 16.0_lkind*zeta_0)
449
450 x = sqrt(x)
451 x_0 = sqrt(x_0)
452
453 x1 = 1.0_lkind + x
454 x1_0 = 1.0_lkind + x_0
455
456 num = x1*x1*(1.0_lkind + x*x)
457 denom = x1_0*x1_0*(1.0_lkind + x_0*x_0)
458 y = atan(x) - atan(x_0)
459 psi_m = ln_z_z0 - log(num/denom) + 2.0_lkind*y
460
461end where
462
463if( stable_option == 1) then
464
465 where (stable)
466 psi_m = ln_z_z0 + (5.0_lkind - real(b_stab, fms_mo_kind_)) &
467 *log((1.0_lkind + zeta)/(1.0_lkind + zeta_0)) &
468 + real(b_stab, fms_mo_kind_)*(zeta - zeta_0)
469 end where
470
471else if (stable_option == 2) then
472
473 weakly_stable = stable .and. zeta <= real(zeta_trans,fms_mo_kind_)
474 strongly_stable = stable .and. zeta > real(zeta_trans,fms_mo_kind_)
475
476 where (weakly_stable)
477 psi_m = ln_z_z0 + 5.0_lkind*(zeta - zeta_0)
478 end where
479
480 where(strongly_stable)
481 x = (real(lambda, fms_mo_kind_) - 1.0_lkind)*log(zeta/real(zeta_trans, fms_mo_kind_)) + &
482 real(b_stab, FMS_MO_KIND_)*(zeta - real(zeta_trans, FMS_MO_KIND_))
483 endwhere
484
485 where (strongly_stable .and. zeta_0 <= real(zeta_trans,fms_mo_kind_))
486 psi_m = ln_z_z0 + x + 5.0_lkind &
487 *(real(zeta_trans, fms_mo_kind_) - zeta_0)
488 end where
489 where (strongly_stable .and. zeta_0 > real(zeta_trans,fms_mo_kind_))
490 psi_m = real(lambda, fms_mo_kind_)*ln_z_z0 + real(b_stab, fms_mo_kind_) &
491 *(zeta - zeta_0)
492 endwhere
493
494end if
495
496return
497end subroutine mo_integral_m_
498
499
500!=======================================================================
501! The following routines allow the public interfaces to be used
502! with different dimensions of the input and output
503!
504!=======================================================================
505
506
507subroutine mo_drag_2d_ &
508 (pt, pt0, z, z0, zt, zq, speed, drag_m, drag_t, drag_q, u_star, b_star)
509
510real(kind=fms_mo_kind_), intent(in) , dimension(:,:) :: z, speed, pt, pt0, z0, zt, zq
511real(kind=fms_mo_kind_), intent(out) , dimension(:,:) :: drag_m, drag_t, drag_q
512real(kind=fms_mo_kind_), intent(inout), dimension(:,:) :: u_star, b_star
513
514integer :: j
515
516do j = 1, size(pt,2)
517 call mo_drag (pt(:,j), pt0(:,j), z(:,j), z0(:,j), zt(:,j), zq(:,j), &
518 speed(:,j), drag_m(:,j), drag_t(:,j), drag_q(:,j), &
519 u_star(:,j), b_star(:,j))
520end do
521
522
523return
524end subroutine mo_drag_2d_
525
526!=======================================================================
527subroutine mo_drag_0d_ &
528 (pt, pt0, z, z0, zt, zq, speed, drag_m, drag_t, drag_q, u_star, b_star)
529
530real(kind=fms_mo_kind_), intent(in) :: z, speed, pt, pt0, z0, zt, zq
531real(kind=fms_mo_kind_), intent(out) :: drag_m, drag_t, drag_q, u_star, b_star
532
533real(kind=fms_mo_kind_), dimension(1) :: pt_1, pt0_1, z_1, z0_1, zt_1, zq_1, speed_1, &
534 drag_m_1, drag_t_1, drag_q_1, u_star_1, b_star_1
535
536pt_1(1) = pt
537pt0_1(1) = pt0
538z_1(1) = z
539z0_1(1) = z0
540zt_1(1) = zt
541zq_1(1) = zq
542speed_1(1) = speed
543
544call mo_drag (pt_1, pt0_1, z_1, z0_1, zt_1, zq_1, speed_1, &
545 drag_m_1, drag_t_1, drag_q_1, u_star_1, b_star_1)
546
547drag_m = drag_m_1(1)
548drag_t = drag_t_1(1)
549drag_q = drag_q_1(1)
550u_star = u_star_1(1)
551b_star = b_star_1(1)
552
553return
554end subroutine mo_drag_0d_
555!=======================================================================
556
557subroutine mo_profile_2d_(zref, zref_t, z, z0, zt, zq, u_star, b_star, q_star, &
558 del_m, del_h, del_q)
559
560real(kind=fms_mo_kind_), intent(in) :: zref, zref_t
561real(kind=fms_mo_kind_), intent(in) , dimension(:,:) :: z, z0, zt, zq, u_star, b_star, q_star
562real(kind=fms_mo_kind_), intent(out), dimension(:,:) :: del_m, del_h, del_q
563
564integer :: j
565
566do j = 1, size(z,2)
567 call mo_profile (zref, zref_t, z(:,j), z0(:,j), zt(:,j), &
568 zq(:,j), u_star(:,j), b_star(:,j), q_star(:,j), &
569 del_m(:,j), del_h(:,j), del_q(:,j))
570enddo
571
572return
573end subroutine mo_profile_2d_
574
575!=======================================================================
576
577subroutine mo_profile_0d_(zref, zref_t, z, z0, zt, zq, u_star, b_star, q_star, &
578 del_m, del_h, del_q)
579
580real(kind=fms_mo_kind_), intent(in) :: zref, zref_t
581real(kind=fms_mo_kind_), intent(in) :: z, z0, zt, zq, u_star, b_star, q_star
582real(kind=fms_mo_kind_), intent(out) :: del_m, del_h, del_q
583
584real(kind=fms_mo_kind_), dimension(1) :: z_1, z0_1, zt_1, zq_1, u_star_1, b_star_1, q_star_1, &
585 del_m_1, del_h_1, del_q_1
586
587z_1(1) = z
588z0_1(1) = z0
589zt_1(1) = zt
590zq_1(1) = zq
591u_star_1(1) = u_star
592b_star_1(1) = b_star
593q_star_1(1) = q_star
594
595call mo_profile (zref, zref_t, z_1, z0_1, zt_1, zq_1, &
596 u_star_1, b_star_1, q_star_1, &
597 del_m_1, del_h_1, del_q_1)
598
599del_m = del_m_1(1)
600del_h = del_h_1(1)
601del_q = del_q_1(1)
602
603
604return
605end subroutine mo_profile_0d_
606
607!=======================================================================
608
609subroutine mo_profile_1d_n_(zref, z, z0, zt, zq, u_star, b_star, q_star, &
610 del_m, del_t, del_q, avail)
611
612real(kind=fms_mo_kind_), intent(in), dimension(:) :: zref
613real(kind=fms_mo_kind_), intent(in) , dimension(:) :: z, z0, zt, zq, u_star, b_star, q_star
614real(kind=fms_mo_kind_), intent(out), dimension(:,:) :: del_m, del_t, del_q
615logical, intent(in) , optional, dimension(:) :: avail
616
617integer :: k
618
619do k = 1, size(zref(:))
620 if(present(avail)) then
621 call mo_profile (zref(k), zref(k), z, z0, zt, zq, &
622 u_star, b_star, q_star, del_m(:,k), del_t(:,k), del_q(:,k), avail)
623 else
624 call mo_profile (zref(k), zref(k), z, z0, zt, zq, &
625 u_star, b_star, q_star, del_m(:,k), del_t(:,k), del_q(:,k))
626 endif
627enddo
628
629return
630end subroutine mo_profile_1d_n_
631
632!=======================================================================
633
634subroutine mo_profile_0d_n_(zref, z, z0, zt, zq, u_star, b_star, q_star, &
635 del_m, del_t, del_q)
636
637real(kind=fms_mo_kind_), intent(in), dimension(:) :: zref
638real(kind=fms_mo_kind_), intent(in) :: z, z0, zt, zq, u_star, b_star, q_star
639real(kind=fms_mo_kind_), intent(out), dimension(:) :: del_m, del_t, del_q
640
641integer :: k
642
643do k = 1, size(zref(:))
644 call mo_profile (zref(k), zref(k), z, z0, zt, zq, &
645 u_star, b_star, q_star, del_m(k), del_t(k), del_q(k))
646enddo
647
648return
649end subroutine mo_profile_0d_n_
650
651!=======================================================================
652
653subroutine mo_profile_2d_n_(zref, z, z0, zt, zq, u_star, b_star, q_star, &
654 del_m, del_t, del_q)
655
656real(kind=fms_mo_kind_), intent(in), dimension(:) :: zref
657real(kind=fms_mo_kind_), intent(in), dimension(:,:) :: z, z0, zt, zq, u_star, b_star, q_star
658real(kind=fms_mo_kind_), intent(out), dimension(:,:,:) :: del_m, del_t, del_q
659
660integer :: k
661
662do k = 1, size(zref(:))
663 call mo_profile (zref(k), zref(k), z, z0, zt, zq, &
664 u_star, b_star, q_star, del_m(:,:,k), del_t(:,:,k), del_q(:,:,k))
665enddo
666
667return
668end subroutine mo_profile_2d_n_
669
670!=======================================================================
671
672subroutine mo_diff_2d_1_(z, u_star, b_star, k_m, k_h)
673
674real(kind=fms_mo_kind_), intent(in), dimension(:,:) :: z, u_star, b_star
675real(kind=fms_mo_kind_), intent(out), dimension(:,:) :: k_m, k_h
676
677real(kind=fms_mo_kind_), dimension(size(z,1),size(z,2),1) :: z_n, k_m_n, k_h_n
678
679z_n(:,:,1) = z
680
681call mo_diff(z_n, u_star, b_star, k_m_n, k_h_n)
682
683k_m = k_m_n(:,:,1)
684k_h = k_h_n(:,:,1)
685
686return
687end subroutine mo_diff_2d_1_
688
689
690!=======================================================================
691
692subroutine mo_diff_1d_1_(z, u_star, b_star, k_m, k_h)
693
694real(kind=fms_mo_kind_), intent(in), dimension(:) :: z, u_star, b_star
695real(kind=fms_mo_kind_), intent(out), dimension(:) :: k_m, k_h
696
697real(kind=fms_mo_kind_), dimension(size(z),1,1) :: z_n, k_m_n, k_h_n
698real(kind=fms_mo_kind_), dimension(size(z),1) :: u_star_n, b_star_n
699
700z_n(:,1,1) = z
701u_star_n(:,1) = u_star
702b_star_n(:,1) = b_star
703
704call mo_diff(z_n, u_star_n, b_star_n, k_m_n, k_h_n)
705
706k_m = k_m_n(:,1,1)
707k_h = k_h_n(:,1,1)
708
709return
710end subroutine mo_diff_1d_1_
711
712!=======================================================================
713
714subroutine mo_diff_1d_n_(z, u_star, b_star, k_m, k_h)
715
716real(kind=fms_mo_kind_), intent(in), dimension(:,:) :: z
717real(kind=fms_mo_kind_), intent(in), dimension(:) :: u_star, b_star
718real(kind=fms_mo_kind_), intent(out), dimension(:,:) :: k_m, k_h
719
720real(kind=fms_mo_kind_), dimension(size(z,1),1) :: u_star2, b_star2
721real(kind=fms_mo_kind_), dimension(size(z,1),1, size(z,2)) :: z2, k_m2, k_h2
722
723integer :: n
724
725do n = 1, size(z,2)
726 z2(:,1,n) = z(:,n)
727enddo
728u_star2(:,1) = u_star
729b_star2(:,1) = b_star
730
731call mo_diff(z2, u_star2, b_star2, k_m2, k_h2)
732
733do n = 1, size(z,2)
734 k_m(:,n) = k_m2(:,1,n)
735 k_h(:,n) = k_h2(:,1,n)
736enddo
737
738return
739end subroutine mo_diff_1d_n_
740
741!=======================================================================
742
743subroutine mo_diff_0d_1_(z, u_star, b_star, k_m, k_h)
744
745real(kind=fms_mo_kind_), intent(in) :: z, u_star, b_star
746real(kind=fms_mo_kind_), intent(out) :: k_m, k_h
747
748integer :: ni, nj, nk, ier
749integer, parameter :: lkind = fms_mo_kind_
750real(kind=fms_mo_kind_), parameter :: ustar_min = 1.0e-10_lkind
751real(kind=fms_mo_kind_), dimension(1,1,1) :: z_a, k_m_a, k_h_a
752real(kind=fms_mo_kind_), dimension(1,1) :: u_star_a, b_star_a
753
754if(.not.module_is_initialized) call error_mesg('mo_diff_0d_1 in monin_obukhov_mod', &
755 'monin_obukhov_init has not been called', fatal)
756
757ni = 1; nj = 1; nk = 1
758z_a(1,1,1) = z
759u_star_a(1,1) = u_star
760b_star_a(1,1) = b_star
761call monin_obukhov_diff(real(vonkarm, fms_mo_kind_), ustar_min, neutral, &
762 & stable_option, new_mo_option, real(rich_crit, fms_mo_kind_), &
763 & real(zeta_trans, fms_mo_kind_), ni, nj, nk, z_a, u_star_a, & !miz
764 & b_star_a, k_m_a, k_h_a, ier)
765k_m = k_m_a(1,1,1)
766k_h = k_h_a(1,1,1)
767
768end subroutine mo_diff_0d_1_
769
770!=======================================================================
771
772subroutine mo_diff_0d_n_(z, u_star, b_star, k_m, k_h)
773
774real(kind=fms_mo_kind_), intent(in), dimension(:) :: z
775real(kind=fms_mo_kind_), intent(in) :: u_star, b_star
776real(kind=fms_mo_kind_), intent(out), dimension(:) :: k_m, k_h
777
778integer :: ni, nj, nk, ier
779integer, parameter :: lkind = fms_mo_kind_
780real(kind=fms_mo_kind_), parameter :: ustar_min = 1.0e-10_lkind
781real(kind=fms_mo_kind_), dimension(1,1,size(z)) :: z_a, k_m_a, k_h_a
782real(kind=fms_mo_kind_), dimension(1,1) :: u_star_a, b_star_a
783
784if(.not.module_is_initialized) call error_mesg('mo_diff_0d_n in monin_obukhov_mod', &
785 'monin_obukhov_init has not been called', fatal)
786
787ni = 1; nj = 1; nk = size(z(:))
788z_a(1,1,:) = z(:)
789u_star_a(1,1) = u_star
790b_star_a(1,1) = b_star
791call monin_obukhov_diff(real(vonkarm, fms_mo_kind_), ustar_min, neutral, &
792 & stable_option, new_mo_option, real(rich_crit, fms_mo_kind_), &
793 & real(zeta_trans, fms_mo_kind_), ni, nj, nk, z_a, u_star_a, &
794 & b_star_a, k_m_a, k_h_a, ier)
795k_m(:) = k_m_a(1,1,:)
796k_h(:) = k_h_a(1,1,:)
797end subroutine mo_diff_0d_n_
798
799!=======================================================================
800
801subroutine stable_mix_2d_(rich, mix)
802
803real(kind=fms_mo_kind_), intent(in) , dimension(:,:) :: rich
804real(kind=fms_mo_kind_), intent(out), dimension(:,:) :: mix
805integer :: n2 !< Size of dimension 2 of mix and rich
806integer :: i !< Loop index
807
808n2 = size(mix, 2)
809
810do i=1, n2
811 call stable_mix(rich(:, i), mix(:, i))
812enddo
813
814end subroutine stable_mix_2d_
815
816
817!=======================================================================
818
819subroutine stable_mix_1d_(rich, mix)
820
821real(kind=fms_mo_kind_), intent(in) , dimension(:) :: rich
822real(kind=fms_mo_kind_), intent(out), dimension(:) :: mix
823integer :: n !< Size of mix and rich
824integer :: ierr !< Error code set by monin_obukhov_stable_mix
825
826if (.not.module_is_initialized) call error_mesg('stable_mix in monin_obukhov_mod', &
827 'monin_obukhov_init has not been called', fatal)
828
829n = size(mix)
830
831call monin_obukhov_stable_mix(stable_option, real(rich_crit,fms_mo_kind_), &
832 & real(zeta_trans,fms_mo_kind_), n, rich, mix, ierr)
833
834end subroutine stable_mix_1d_
835
836!=======================================================================
837
838subroutine stable_mix_0d_(rich, mix)
839
840real(kind=fms_mo_kind_), intent(in) :: rich
841real(kind=fms_mo_kind_), intent(out) :: mix
842real(kind=fms_mo_kind_), dimension(1) :: mix_1d !< Representation of mix as a dimension(1) array
843
844call stable_mix([rich], mix_1d)
845
846mix = mix_1d(1)
847
848end subroutine stable_mix_0d_
849!=======================================================================
850
851!> @}
852! close documentation grouping