40 mpp_pe, mpp_root_pe, close_file, stdlog, &
116 integer :: unit, ierr, io, logunit
120 #ifdef INTERNAL_FILE_NML 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)
130 10
call close_file (unit)
136 if ( mpp_pe() == mpp_root_pe() )
then 139 write (logunit, nml=monin_obukhov_nml)
145 'MONIN_OBUKHOV_INIT in MONIN_OBUKHOV_MOD', &
146 'rich_crit in monin_obukhov_mod must be > 0.25', fatal)
149 'MONIN_OBUKHOV_INIT in MONIN_OBUKHOV_MOD', &
150 'drag_min_heat in monin_obukhov_mod must be >= 0.0', fatal)
153 'MONIN_OBUKHOV_INIT in MONIN_OBUKHOV_MOD', &
154 'drag_min_moist in monin_obukhov_mod must be >= 0.0', fatal)
157 'MONIN_OBUKHOV_INIT in MONIN_OBUKHOV_MOD', &
158 'drag_min_mom in monin_obukhov_mod must be >= 0.0', fatal)
161 'MONIN_OBUKHOV_INIT in MONIN_OBUKHOV_MOD', &
162 'the only allowable values of stable_option are 1 and 2', fatal)
165 'MONIN_OBUKHOV_INIT in MONIN_OBUKHOV_MOD', &
166 'zeta_trans must be positive', fatal)
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
207 logical :: lavail, avail_dummy(1)
210 integer,
parameter :: max_iter = 20
211 real ,
parameter :: error=1.e-04, zeta_min=1.e-06,
small=1.e-04
216 'monin_obukhov_init has not been called', fatal)
220 if(
present(avail)) lavail = .true.
224 if (count(avail) .eq. 0)
return 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)
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)
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
253 logical :: dummy_avail(1)
259 'monin_obukhov_init has not been called', fatal)
262 if(
present(avail))
then 264 if (count(avail) .eq. 0)
return 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)
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)
286 real,
intent(in) ,
dimension(:,:,:) :: rich
287 real,
intent(out),
dimension(:,:,:) :: mix
292 'monin_obukhov_init has not been called', fatal)
294 n =
size(rich,1)*
size(rich,2)*
size(rich,3)
305 real,
intent(in),
dimension(:,:,:) :: z
306 real,
intent(in),
dimension(:,:) :: u_star, b_star
307 real,
intent(out),
dimension(:,:,:) :: k_m, k_h
309 integer :: ni, nj, nk, ier
310 real,
parameter :: ustar_min = 1.e-10
313 'monin_obukhov_init has not been called', fatal)
315 ni =
size(z, 1); nj =
size(z, 2); nk =
size(z, 3)
319 & ni, nj, nk, z, u_star, b_star, k_m, k_h, ier)
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
334 real,
parameter :: error = 1.e-04
335 real,
parameter :: zeta_min = 1.e-06
336 integer,
parameter :: max_iter = 20
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
347 logical,
dimension(size(rich(:))) :: mask_1
363 zeta = rich*ln_z_z0*ln_z_z0/ln_z_zt
368 where (mask_1 .and. rich >= 0.0)
372 iter_loop:
do iter = 1, max_iter
374 where (mask_1 .and. abs(zeta).lt.zeta_min)
399 call mo_integral_tq(f_t, f_q, zeta, zeta_t, zeta_q, ln_z_zt, ln_z_zq, 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))
412 max_cor= maxval(corr)
414 if(max_cor > error)
then 415 mask_1 = mask_1 .and. (corr > error)
418 zeta = zeta + correction
427 call error_mesg (
'solve_zeta in monin_obukhov_mod', &
428 'surface drag iteration did not converge', fatal)
438 real ,
intent(out),
dimension(:) :: phi_m
439 real ,
intent(in),
dimension(:) :: zeta
440 logical ,
intent(in),
dimension(:) :: mask
442 logical,
dimension(size(zeta(:))) :: stable, unstable
443 real ,
dimension(size(zeta(:))) :: x
445 stable = mask .and. zeta >= 0.0
446 unstable = mask .and. zeta < 0.0
449 x = (1 - 16.0*zeta )**(-0.5)
456 phi_m = 1.0 + zeta *(5.0 +
b_stab*zeta)/(1.0 + zeta)
479 real ,
intent(out),
dimension(:) :: phi_t
480 real ,
intent(in),
dimension(:) :: zeta
481 logical ,
intent(in),
dimension(:) :: mask
483 logical,
dimension(size(zeta(:))) :: stable, unstable
485 stable = mask .and. zeta >= 0.0
486 unstable = mask .and. zeta < 0.0
489 phi_t = (1 - 16.0*zeta)**(-0.5)
495 phi_t = 1.0 + zeta*(5.0 +
b_stab*zeta)/(1.0 + zeta)
515 ln_z_zt, ln_z_zq, mask)
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
523 real,
dimension(size(zeta(:))) :: x, x_t, x_q
525 logical,
dimension(size(zeta(:))) :: stable, unstable, &
526 weakly_stable, strongly_stable
528 stable = mask .and. zeta >= 0.0
529 unstable = mask .and. zeta < 0.0
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)
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) )
546 psi_t = ln_z_zt + (5.0 -
b_stab)*log((1.0 + zeta)/(1.0 + zeta_t)) &
548 psi_q = ln_z_zq + (5.0 -
b_stab)*log((1.0 + zeta)/(1.0 + zeta_q)) &
555 weakly_stable = stable .and. zeta <=
zeta_trans 556 strongly_stable = stable .and. zeta >
zeta_trans 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)
563 where(strongly_stable)
567 where (strongly_stable .and. zeta_t <=
zeta_trans)
568 psi_t = ln_z_zt + x + 5.0*(
zeta_trans - zeta_t)
570 where (strongly_stable .and. zeta_t >
zeta_trans)
574 where (strongly_stable .and. zeta_q <=
zeta_trans)
575 psi_q = ln_z_zq + x + 5.0*(
zeta_trans - zeta_q)
577 where (strongly_stable .and. zeta_q >
zeta_trans)
588 subroutine mo_integral_m (psi_m, zeta, zeta_0, ln_z_z0, mask)
592 real ,
intent(out),
dimension(:) :: psi_m
593 real ,
intent(in),
dimension(:) :: zeta, zeta_0, ln_z_z0
594 logical ,
intent(in),
dimension(:) :: mask
596 real,
dimension(size(zeta(:))) :: x, x_0, x1, x1_0, num, denom, y
598 logical,
dimension(size(zeta(:))) :: stable, unstable, &
599 weakly_stable, strongly_stable
601 stable = mask .and. zeta >= 0.0
602 unstable = mask .and. zeta < 0.0
606 x = sqrt(1 - 16.0*zeta)
607 x_0 = sqrt(1 - 16.0*zeta_0)
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
625 psi_m = ln_z_z0 + (5.0 -
b_stab)*log((1.0 + zeta)/(1.0 + zeta_0)) &
631 weakly_stable = stable .and. zeta <=
zeta_trans 632 strongly_stable = stable .and. zeta >
zeta_trans 634 where (weakly_stable)
635 psi_m = ln_z_z0 + 5.0*(zeta - zeta_0)
638 where(strongly_stable)
642 where (strongly_stable .and. zeta_0 <=
zeta_trans)
643 psi_m = ln_z_z0 + x + 5.0*(
zeta_trans - zeta_0)
645 where (strongly_stable .and. zeta_0 >
zeta_trans)
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
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))
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
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
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)
712 subroutine mo_profile_2d(zref, zref_t, z, z0, zt, zq, u_star, b_star, q_star, &
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
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))
732 subroutine mo_profile_0d(zref, zref_t, z, z0, zt, zq, u_star, b_star, q_star, &
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
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
751 u_star_1, b_star_1, q_star_1, &
752 del_m_1, del_h_1, del_q_1)
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
774 do k = 1,
size(zref(:))
775 if(
present(avail))
then 777 u_star, b_star, q_star, del_m(:,k), del_t(:,k), del_q(:,k), avail)
780 u_star, b_star, q_star, del_m(:,k), del_t(:,k), del_q(:,k))
789 subroutine mo_profile_0d_n(zref, z, z0, zt, zq, u_star, b_star, q_star, &
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
798 do k = 1,
size(zref(:))
800 u_star, b_star, q_star, del_m(k), del_t(k), del_q(k))
808 subroutine mo_profile_2d_n(zref, z, z0, zt, zq, u_star, b_star, q_star, &
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
817 do k = 1,
size(zref(:))
819 u_star, b_star, q_star, del_m(:,:,k), del_t(:,:,k), del_q(:,:,k))
829 real,
intent(in),
dimension(:,:) :: z, u_star, b_star
830 real,
intent(out),
dimension(:,:) :: k_m, k_h
832 real ,
dimension(size(z,1),size(z,2),1) :: z_n, k_m_n, k_h_n
849 real,
intent(in),
dimension(:) :: z, u_star, b_star
850 real,
intent(out),
dimension(:) :: k_m, k_h
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
856 u_star_n(:,1) = u_star
857 b_star_n(:,1) = b_star
859 call mo_diff_2d_n(z_n, u_star_n, b_star_n, k_m_n, k_h_n)
871 real,
intent(in),
dimension(:,:) :: z
872 real,
intent(in),
dimension(:) :: u_star, b_star
873 real,
intent(out),
dimension(:,:) :: k_m, k_h
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
883 u_star2(:,1) = u_star
884 b_star2(:,1) = b_star
889 k_m(:,n) = k_m2(:,1,n)
890 k_h(:,n) = k_h2(:,1,n)
900 real,
intent(in) :: z, u_star, b_star
901 real,
intent(out) :: k_m, k_h
903 integer :: ni, nj, nk, ier
904 real,
parameter :: ustar_min = 1.e-10
907 'monin_obukhov_init has not been called', fatal)
909 ni = 1; nj = 1; nk = 1
913 & ni, nj, nk, z, u_star, b_star, k_m, k_h, ier)
921 real,
intent(in),
dimension(:) :: z
922 real,
intent(in) :: u_star, b_star
923 real,
intent(out),
dimension(:) :: k_m, k_h
925 integer :: ni, nj, nk, ier
926 real,
parameter :: ustar_min = 1.e-10
929 'monin_obukhov_init has not been called', fatal)
931 ni = 1; nj = 1; nk =
size(z(:))
935 & ni, nj, nk, z, u_star, b_star, k_m, k_h, ier)
943 real,
intent(in) ,
dimension(:,:) :: rich
944 real,
intent(out),
dimension(:,:) :: mix
946 real,
dimension(size(rich,1),size(rich,2),1) :: rich_3d, mix_3d
948 rich_3d(:,:,1) = rich
962 real,
intent(in) ,
dimension(:) :: rich
963 real,
intent(out),
dimension(:) :: mix
965 real,
dimension(size(rich),1,1) :: rich_3d, mix_3d
967 rich_3d(:,1,1) = rich
980 real,
intent(in) :: rich
981 real,
intent(out) :: mix
983 real,
dimension(1,1,1) :: rich_3d, mix_3d
985 rich_3d(1,1,1) = rich
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].
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)
integer function, public check_nml_error(IOSTAT, NML_NAME)
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
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].
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)
subroutine mo_diff_0d_1(z, u_star, b_star, k_m, k_h)
subroutine, public error_mesg(routine, message, level)
subroutine mo_diff_2d_n(z, u_star, b_star, k_m, k_h)