123 type(model_state_type),
target,
intent(inout) :: current_state
127 if (.not. is_component_enabled(current_state%options_database,
"smagorinsky"))
then
128 call log_master_log(log_error,
"Subgrid model diags requested but subgrid model not enabled - check config")
131 vertical_grid=current_state%global_grid%configuration%vertical
133 total_points=(current_state%local_grid%size(y_index) * current_state%local_grid%size(x_index))
135 allocate(
uwsg_tot(current_state%local_grid%size(z_index)) &
136 ,
vwsg_tot(current_state%local_grid%size(z_index)) &
137 ,
uusg_tot(current_state%local_grid%size(z_index)) &
138 ,
vvsg_tot(current_state%local_grid%size(z_index)) &
139 ,
wwsg_tot(current_state%local_grid%size(z_index)) &
140 ,
tkesg_tot(current_state%local_grid%size(z_index)) &
141 ,
wkesg_tot(current_state%local_grid%size(z_index)) &
142 ,
vis_coef_tot(current_state%local_grid%size(z_index)) &
147 ,
ssub_tot(current_state%local_grid%size(z_index)) &
148 ,
sed_tot(current_state%local_grid%size(z_index)) &
149 ,
buoysg_tot(current_state%local_grid%size(z_index)) )
152 allocate(
ssq(current_state%local_grid%size(z_index)))
153 allocate(
elamr_sq(current_state%local_grid%size(z_index)))
155 allocate(
subgrid_tke(current_state%local_grid%size(z_index)))
158 allocate(
subke_2d(current_state%local_grid%size(y_index),current_state%local_grid%size(x_index)))
159 if (current_state%th%active)
then
160 allocate(
epsth(current_state%local_grid%size(z_index)) &
162 ,
wtsg_tot(current_state%local_grid%size(z_index)) &
163 ,
th2sg_tot(current_state%local_grid%size(z_index)))
166 if (.not. current_state%passive_q .and. current_state%number_q_fields .gt. 0)
then
167 allocate(
wqsg_tot(current_state%local_grid%size(z_index)))
168 iqv=get_q_index(standard_q_names%VAPOUR,
'subgrid_profile_diags')
169 iql=get_q_index(standard_q_names%CLOUD_LIQUID_MASS,
'subgrid_profile_diags')
170 qlcrit=options_get_real(current_state%options_database,
"qlcrit")
171 allocate(
wqv_sg_tot(current_state%local_grid%size(z_index)), &
172 wql_sg_tot(current_state%local_grid%size(z_index)))
174 if (current_state%rain_water_mixing_ratio_index > 0)
then
175 iqr = current_state%rain_water_mixing_ratio_index
176 allocate(
wqr_sg_tot(current_state%local_grid%size(z_index)))
178 if (current_state%ice_water_mixing_ratio_index > 0)
then
179 iqi = current_state%ice_water_mixing_ratio_index
180 allocate(
wqi_sg_tot(current_state%local_grid%size(z_index)))
182 if (current_state%snow_water_mixing_ratio_index > 0)
then
183 iqs = current_state%snow_water_mixing_ratio_index
184 allocate(
wqs_sg_tot(current_state%local_grid%size(z_index)))
186 if (current_state%graupel_water_mixing_ratio_index > 0)
then
187 iqg = current_state%graupel_water_mixing_ratio_index
188 allocate(
wqg_sg_tot(current_state%local_grid%size(z_index)))
197 type(model_state_type),
target,
intent(inout) :: current_state
200 integer :: jcol, icol, target_x_index, target_y_index
202 real(kind=default_precision),
dimension(current_state%local_grid%size(Z_INDEX)) :: &
203 s11, s22, s33, s12, s23, s13, s33_on_p
204 real(kind=default_precision),
dimension(current_state%local_grid%size(Z_INDEX)) :: &
205 tau11,tau22, tau33, tau12, tau23, tau13, tau33_on_p
206 real(kind=default_precision),
dimension(current_state%local_grid%size(Z_INDEX)) :: u_i_prime_tau_i
207 real(kind=default_precision) :: vistmp, vis12, vis12a, visonp2, visonp2a, vis13, vis23
208 real(kind=default_precision),
dimension(current_state%local_grid%size(Z_INDEX)) :: umean, wu_umean, vmean, wv_vmean, &
209 w_pprime_at_w, rke1, w_qvprime_at_w, w_qclprime_at_w, w_thprime_at_w, wq, rho, rec_rho, buoy_cof, &
210 uw_tot, vw_tot,qvprime_at_w,qclprime_at_w,thprime_at_w, sed_eq13, sed_eq23, sed33,uwsg,vwsg,sed_swap
211 real(kind=default_precision) :: c_virtual
212 real(kind=default_precision) :: upr_at_w,vpr_at_w, vprime_w_local, uprime_w_local, u_1, u_1_bar
213 real(kind=default_precision) :: buoy_prod, sg_shear_prod, dissipation_rate
215 logical :: use_Ri_for_buoyant_prod=.true.
217 c_virtual = (ratio_mol_wts-1.0_default_precision)
218 jcol=current_state%column_local_y
219 icol=current_state%column_local_x
220 target_y_index=jcol-current_state%local_grid%halo_size(y_index)
221 target_x_index=icol-current_state%local_grid%halo_size(x_index)
224 if (current_state%first_timestep_column)
then
225 sed_tot(:) = 0.0_default_precision
228 sed_eq23(:) = 0.0_default_precision
229 sed_eq13(:) = 0.0_default_precision
230 sed_swap(:) = 0.0_default_precision
243 subke_2d(:,:) = 0.0_default_precision
244 if (current_state%th%active)
then
249 if (.not. current_state%passive_q .and. &
250 current_state%number_q_fields .gt. 0)
then
261 if (.not. current_state%halo_column)
then
263 do k=2,current_state%local_grid%size(z_index)-1
270 s11(k)=current_state%global_grid%configuration%horizontal%cx*(&
271 (current_state%u%data(k+1, jcol, icol)-&
272 current_state%u%data(k+1, jcol, icol-1))+&
273 (current_state%u%data(k, jcol, icol)-&
274 current_state%u%data(k, jcol, icol-1)))
276 s11(k)=0.0_default_precision
280 s22(k)=current_state%global_grid%configuration%horizontal%cy*(&
281 (current_state%v%data(k+1, jcol, icol)-&
282 current_state%v%data(k+1, jcol-1, icol))+&
283 (current_state%v%data(k, jcol, icol)-&
284 current_state%v%data(k, jcol-1, icol)))
286 s22(k)=0.0_default_precision
290 s33(k)=((current_state%w%data(k, jcol, icol)-&
291 current_state%w%data(k-1, jcol, icol))*&
292 current_state%global_grid%configuration%vertical%rdz(k)) +&
293 ((current_state%w%data(k+1, jcol, icol)-&
294 current_state%w%data(k, jcol, icol))*&
295 current_state%global_grid%configuration%vertical%rdz(k+1))
298 s33_on_p(k) = 2.0_default_precision * &
299 ((current_state%w%data(k, jcol, icol)-&
300 current_state%w%data(k-1, jcol, icol))*&
301 current_state%global_grid%configuration%vertical%rdz(k))
303 s33(k)=0.0_default_precision
304 s33_on_p(k) = 0.0_default_precision
306 #if defined(U_ACTIVE) && defined(W_ACTIVE)
308 s13(k)=(((current_state%u%data(k+1, jcol, icol)-&
309 current_state%u%data(k, jcol, icol))*&
310 current_state%global_grid%configuration%vertical%rdzn(k+1)+&
311 (current_state%w%data(k, jcol, icol+1)-&
312 current_state%w%data(k, jcol, icol))*&
313 current_state%global_grid%configuration%horizontal%cx)+&
314 ((current_state%u%data(k+1, jcol, icol-1)-&
315 current_state%u%data(k, jcol, icol-1))*&
316 current_state%global_grid%configuration%vertical%rdzn(k+1)+&
317 (current_state%w%data(k, jcol, icol)-&
318 current_state%w%data(k, jcol, icol-1))*&
319 current_state%global_grid%configuration%horizontal%cx))*0.5_default_precision
320 #elif defined(U_ACTIVE) && !defined(W_ACTIVE)
321 s13(k)=(((current_state%u%data(k+1, jcol, icol)-&
322 current_state%u%data(k, jcol, icol))*&
323 current_state%global_grid%configuration%vertical%rdzn(k+1))+&
324 ((current_state%u%data(k+1, jcol, icol-1)-&
325 current_state%u%data(k, jcol, icol-1))*&
326 current_state%global_grid%configuration%vertical%rdzn(k+1)))*0.5_default_precision
327 #elif !defined(U_ACTIVE) && defined(W_ACTIVE)
328 s13(k)=(((current_state%w%data(k, jcol, icol+1)-&
329 current_state%w%data(k, jcol, icol))*&
330 current_state%global_grid%configuration%horizontal%cx)+&
331 ((current_state%w%data(k, jcol, icol)-&
332 current_state%w%data(k, jcol, icol-1))*&
333 current_state%global_grid%configuration%horizontal%cx))*0.5_default_precision
335 s13(k)=0.0_default_precision
337 #if defined(W_ACTIVE) && defined(V_ACTIVE)
339 s23(k)=(((current_state%w%data(k, jcol, icol) - &
340 current_state%w%data(k, jcol-1, icol)) * &
341 current_state%global_grid%configuration%horizontal%cy + &
342 (current_state%v%data(k+1, jcol-1, icol) - &
343 current_state%v%data(k, jcol-1, icol)) * &
344 current_state%global_grid%configuration%vertical%rdzn(k+1)) + &
345 ((current_state%w%data(k, jcol+1, icol) - &
346 current_state%w%data(k, jcol, icol)) * &
347 current_state%global_grid%configuration%horizontal%cy + &
348 (current_state%v%data(k+1, jcol, icol)-&
349 current_state%v%data(k, jcol, icol))*&
350 current_state%global_grid%configuration%vertical%rdzn(k+1)))*0.5_default_precision
351 #elif defined(W_ACTIVE) && !defined(V_ACTIVE)
352 s23(k)=(((current_state%w%data(k, jcol, icol)-&
353 current_state%w%data(k, jcol-1, icol))*&
354 current_state%global_grid%configuration%horizontal%cy)+&
355 ((current_state%w%data(k, jcol+1, icol)-&
356 current_state%w%data(k, jcol, icol))*&
357 current_state%global_grid%configuration%horizontal%cy))*0.5_default_precision
358 #elif !defined(W_ACTIVE) && defined(V_ACTIVE)
359 s23(k)=(((current_state%v%data(k+1, jcol-1, icol)-&
360 current_state%v%data(k, jcol-1, icol))*&
361 current_state%global_grid%configuration%vertical%rdzn(k+1))+&
362 ((current_state%v%data(k+1, jcol, icol)-&
363 current_state%v%data(k, jcol, icol))*&
364 current_state%global_grid%configuration%vertical%rdzn(k+1)))*0.5_default_precision
366 s23(k)=0.0_default_precision
368 #if defined(U_ACTIVE) && defined(V_ACTIVE)
370 s12(k)=(((((current_state%u%data(k, jcol, icol-1)-&
371 current_state%u%data(k, jcol-1, icol-1))*&
372 current_state%global_grid%configuration%horizontal%cy+&
373 (current_state%v%data(k, jcol-1, icol)-&
374 current_state%v%data(k, jcol-1, icol-1))*&
375 current_state%global_grid%configuration%horizontal%cx) +&
376 ((current_state%u%data(k, jcol+1, icol-1)-&
377 current_state%u%data(k, jcol, icol-1))*&
378 current_state%global_grid%configuration%horizontal%cy+&
379 (current_state%v%data(k, jcol, icol)-&
380 current_state%v%data(k, jcol, icol-1))*&
381 current_state%global_grid%configuration%horizontal%cx)) +(&
382 ((current_state%u%data(k, jcol, icol)-&
383 current_state%u%data(k, jcol-1, icol))*&
384 current_state%global_grid%configuration%horizontal%cy+&
385 (current_state%v%data(k, jcol-1, icol+1)-&
386 current_state%v%data(k, jcol-1, icol))*&
387 current_state%global_grid%configuration%horizontal%cx) +&
388 ((current_state%u%data(k, jcol+1, icol)-&
389 current_state%u%data(k, jcol, icol))*&
390 current_state%global_grid%configuration%horizontal%cy+&
391 (current_state%v%data(k, jcol, icol+1)-&
392 current_state%v%data(k, jcol, icol))*&
393 current_state%global_grid%configuration%horizontal%cx)))+((&
394 ((current_state%u%data(k+1, jcol, icol-1)-&
395 current_state%u%data(k+1, jcol-1, icol-1))*&
396 current_state%global_grid%configuration%horizontal%cy+&
397 (current_state%v%data(k+1, jcol-1, icol)-&
398 current_state%v%data(k+1, jcol-1, icol-1))*&
399 current_state%global_grid%configuration%horizontal%cx)+&
400 ((current_state%u%data(k+1, jcol+1, icol-1)-&
401 current_state%u%data(k+1, jcol, icol-1))*&
402 current_state%global_grid%configuration%horizontal%cy+&
403 (current_state%v%data(k+1, jcol, icol)-&
404 current_state%v%data(k+1, jcol, icol-1))*&
405 current_state%global_grid%configuration%horizontal%cx))+(&
406 ((current_state%u%data(k+1, jcol, icol)-&
407 current_state%u%data(k+1, jcol-1, icol))*&
408 current_state%global_grid%configuration%horizontal%cy+&
409 (current_state%v%data(k+1, jcol-1, icol+1)-&
410 current_state%v%data(k+1, jcol-1, icol))*&
411 current_state%global_grid%configuration%horizontal%cx)+&
412 ((current_state%u%data(k+1, jcol+1, icol)-&
413 current_state%u%data(k+1, jcol, icol))*&
414 current_state%global_grid%configuration%horizontal%cy+&
415 (current_state%v%data(k+1, jcol, icol+1)-&
416 current_state%v%data(k+1, jcol, icol))*&
417 current_state%global_grid%configuration%horizontal%cx))))*0.125_default_precision
419 #elif defined(U_ACTIVE) && !defined(V_ACTIVE)
421 s12(k)=(((((current_state%u%data(k, jcol, icol-1)-&
422 current_state%u%data(k, jcol-1, icol-1))*&
423 current_state%global_grid%configuration%horizontal%cy) +&
424 ((current_state%u%data(k, jcol+1, icol-1)-&
425 current_state%u%data(k, jcol, icol-1))*&
426 current_state%global_grid%configuration%horizontal%cy)) +(&
427 ((current_state%u%data(k, jcol, icol)-&
428 current_state%u%data(k, jcol-1, icol))*&
429 current_state%global_grid%configuration%horizontal%cy) +&
430 ((current_state%u%data(k, jcol+1, icol)-&
431 current_state%u%data(k, jcol, icol))*&
432 current_state%global_grid%configuration%horizontal%cy)))+((&
433 ((current_state%u%data(k+1, jcol, icol-1)-&
434 current_state%u%data(k+1, jcol-1, icol-1))*&
435 current_state%global_grid%configuration%horizontal%cy)+&
436 ((current_state%u%data(k+1, jcol+1, icol-1)-&
437 current_state%u%data(k+1, jcol, icol-1))*&
438 current_state%global_grid%configuration%horizontal%cy))+(&
439 ((current_state%u%data(k+1, jcol, icol)-&
440 current_state%u%data(k+1, jcol-1, icol))*&
441 current_state%global_grid%configuration%horizontal%cy)+&
442 ((current_state%u%data(k+1, jcol+1, icol)-&
443 current_state%u%data(k+1, jcol, icol))*&
444 current_state%global_grid%configuration%horizontal%cy))))*0.125_default_precision
446 #elif !defined(U_ACTIVE) && defined(V_ACTIVE)
448 s12(k)=(((((current_state%v%data(k, jcol-1, icol)-&
449 current_state%v%data(k, jcol-1, icol-1))*&
450 current_state%global_grid%configuration%horizontal%cx) +&
451 ((current_state%v%data(k, jcol, icol)-&
452 current_state%v%data(k, jcol, icol-1))*&
453 current_state%global_grid%configuration%horizontal%cx)) +&
454 (((current_state%v%data(k, jcol-1, icol+1)-&
455 current_state%v%data(k, jcol-1, icol))*&
456 current_state%global_grid%configuration%horizontal%cx) +&
457 ((current_state%v%data(k, jcol, icol+1)-&
458 current_state%v%data(k, jcol, icol))*&
459 current_state%global_grid%configuration%horizontal%cx)))+((&
460 ((current_state%v%data(k+1, jcol-1, icol)-&
461 current_state%v%data(k+1, jcol-1, icol-1))*&
462 current_state%global_grid%configuration%horizontal%cx)+&
463 ((current_state%v%data(k+1, jcol, icol)-&
464 current_state%v%data(k+1, jcol, icol-1))*&
465 current_state%global_grid%configuration%horizontal%cx))+(&
466 ((current_state%v%data(k+1, jcol-1, icol+1)-&
467 current_state%v%data(k+1, jcol-1, icol))*&
468 current_state%global_grid%configuration%horizontal%cx)+&
469 ((current_state%v%data(k+1, jcol, icol+1)-&
470 current_state%v%data(k+1, jcol, icol))*&
471 current_state%global_grid%configuration%horizontal%cx))))*0.125_default_precision
473 s12(k)=0.0_default_precision
477 #if defined(U_ACTIVE)
479 s13(1)=(current_state%u%data(2, jcol, icol) + &
480 current_state%u%data(2, jcol, icol-1)) * &
481 0.5_default_precision / current_state%global_grid%configuration%vertical%zn(2)
483 s13(1)=0.0_default_precision
485 #if defined(V_ACTIVE)
487 s23(1)=(current_state%v%data(2, jcol, icol) + &
488 current_state%v%data(2, jcol-1, icol)) * &
489 0.5_default_precision / current_state%global_grid%configuration%vertical%zn(2)
491 s23(1)=0.0_default_precision
493 s12(1)=0.0_default_precision
495 do k=1,current_state%local_grid%size(z_index)-1
496 tau11(k) = current_state%global_grid%configuration%vertical%rho(k) * &
497 current_state%vis_coefficient%data(k,jcol,icol) * &
499 tau22(k) = current_state%global_grid%configuration%vertical%rho(k) * &
500 current_state%vis_coefficient%data(k,jcol,icol) * &
502 tau33(k) = current_state%global_grid%configuration%vertical%rho(k) * &
503 current_state%vis_coefficient%data(k,jcol,icol) * &
505 tau12(k) = current_state%global_grid%configuration%vertical%rho(k) * &
506 current_state%vis_coefficient%data(k,jcol,icol) * &
508 tau13(k) = current_state%global_grid%configuration%vertical%rho(k) * &
509 current_state%vis_coefficient%data(k,jcol,icol) * &
511 tau23(k) = current_state%global_grid%configuration%vertical%rho(k) * &
512 current_state%vis_coefficient%data(k,jcol,icol) * &
516 do k=2,current_state%local_grid%size(z_index)-1
517 tau33_on_p(k) = current_state%global_grid%configuration%vertical%rhon(k) * 0.5 *&
518 (current_state%vis_coefficient%data(k-1,jcol,icol) + &
519 current_state%vis_coefficient%data(k, jcol,icol)) * &
524 do k=1, current_state%local_grid%size(z_index)-1
529 current_state%vis_coefficient%data(k,jcol,icol-1) * &
530 ( current_state%u%data(k+1,jcol,icol-1) - &
531 current_state%u%data(k,jcol,icol-1) ) * &
534 current_state%vis_coefficient%data(k,jcol,icol) * &
535 ( current_state%u%data(k+1,jcol,icol) - &
536 current_state%u%data(k,jcol,icol) ) * &
542 current_state%vis_coefficient%data(k,jcol-1,icol) * &
543 ( current_state%v%data(k+1,jcol-1,icol) - &
544 current_state%v%data(k,jcol-1,icol) ) * &
547 current_state%vis_coefficient%data(k,jcol,icol) * &
548 ( current_state%v%data(k+1,jcol,icol) - &
549 current_state%v%data(k,jcol,icol) ) * &
554 current_state%vis_coefficient%data(k,jcol,icol)
556 current_state%diff_coefficient%data(k,jcol,icol)
560 if (current_state%th%active)
then
561 do k=1, current_state%local_grid%size(z_index)-1
563 current_state%diff_coefficient%data(k,jcol,icol) * &
564 ( current_state%th%data(k+1,jcol,icol) + &
566 current_state%th%data(k,jcol,icol) - &
573 if (current_state%th%active .and. .not. current_state%passive_q .and. current_state%number_q_fields .gt. 0)
then
575 do k=1, current_state%local_grid%size(z_index)-1
577 current_state%diff_coefficient%data(k, jcol, icol)* &
579 (current_state%q(
iqv)%data(k,jcol,icol) - current_state%q(
iqv)%data(k+1,jcol,icol))
581 current_state%diff_coefficient%data(k, jcol, icol)* &
583 (current_state%q(
iql)%data(k,jcol,icol) - current_state%q(
iql)%data(k+1,jcol,icol))
586 do k=1, current_state%local_grid%size(z_index)-1
588 current_state%diff_coefficient%data(k, jcol, icol)* &
590 (current_state%q(
iqr)%data(k,jcol,icol) - current_state%q(
iqr)%data(k+1,jcol,icol))
594 do k=1, current_state%local_grid%size(z_index)-1
596 current_state%diff_coefficient%data(k, jcol, icol)* &
598 (current_state%q(
iqi)%data(k,jcol,icol) - current_state%q(
iqi)%data(k+1,jcol,icol))
602 do k=1, current_state%local_grid%size(z_index)-1
604 current_state%diff_coefficient%data(k, jcol, icol)* &
606 (current_state%q(
iqs)%data(k,jcol,icol) - current_state%q(
iqs)%data(k+1,jcol,icol))
610 do k=1, current_state%local_grid%size(z_index)-1
612 current_state%diff_coefficient%data(k, jcol, icol)* &
614 (current_state%q(
iqg)%data(k,jcol,icol) - current_state%q(
iqg)%data(k+1,jcol,icol))
621 u_1=sqrt(current_state%u%data(2,jcol,icol) * &
622 current_state%u%data(2,jcol,icol) + &
623 current_state%v%data(2,jcol,icol) * &
624 current_state%v%data(2,jcol,icol))
625 u_1_bar=sqrt(current_state%global_grid%configuration%vertical%olubar(2)*&
626 current_state%global_grid%configuration%vertical%olubar(2)+&
627 current_state%global_grid%configuration%vertical%olvbar(2)*&
628 current_state%global_grid%configuration%vertical%olvbar(2))
631 do k=1, current_state%local_grid%size(z_index)-1
632 rho(k) = current_state%global_grid%configuration%vertical%rho(k)
637 do k=1, current_state%local_grid%size(z_index)
638 buoy_cof(k)=g/current_state%global_grid%configuration%vertical%thref(k)
646 ri_crit = 0.25_default_precision
648 ssq=calculate_half_squared_strain_rate(current_state, current_state%u, current_state%v, current_state%w)
649 richardson_number=calculate_richardson_number(current_state,
ssq, current_state%th, current_state%q)
651 do k=2, current_state%local_grid%size(z_index)-1
656 current_state%vis_coefficient%data(k, jcol, icol) &
659 current_state%vis_coefficient%data(k, jcol, icol) / &
661 current_state%diff_coefficient%data(k, jcol, icol) / &
662 current_state%vis_coefficient%data(k, jcol, icol) )
665 dissipation_rate =
ssq(k) * ( &
666 current_state%vis_coefficient%data(k, jcol, icol) - &
668 current_state%diff_coefficient%data(k, jcol, icol) )
670 if (use_ri_for_buoyant_prod)
then
673 current_state%diff_coefficient%data(k, jcol, icol)
676 if (current_state%th%active .and. &
677 .not. current_state%passive_q .and. current_state%number_q_fields .gt. 0)
then
679 -current_state%diff_coefficient%data(k,jcol,icol) * &
681 current_state%global_grid%configuration%vertical%rdzn(k+1) * &
685 (current_state%th%data(k+1,jcol,icol) + &
686 current_state%global_grid%configuration%vertical%thref(k+1) + &
688 current_state%global_grid%configuration%vertical%thref(k+1) * &
689 (c_virtual * current_state%q(1)%data(k+1,jcol,icol) - &
690 current_state%q(2)%data(k+1,jcol,icol))) - &
693 (current_state%th%data(k,jcol,icol) + &
694 current_state%global_grid%configuration%vertical%thref(k) + &
695 current_state%global_grid%configuration%vertical%thref(k) * &
696 (c_virtual * current_state%q(1)%data(k,jcol,icol) - &
697 current_state%q(2)%data(k,jcol,icol))) )
699 call log_master_log(log_error, &
700 "Subgrid diags - buoy_prod calc needs theta and q active, STOP")
705 dissipation_rate =
ssq(k) * &
706 current_state%vis_coefficient%data(k, jcol, icol) + &
728 current_state%w%data(k,jcol,icol)
734 subke_2d(target_y_index, target_x_index) =
subke_2d(target_y_index, target_x_index) + &
745 if (current_state%th%active .and. &
746 .not. current_state%passive_q .and. current_state%number_q_fields .gt. 0)
then
749 -current_state%diff_coefficient%data(k,jcol,icol) * &
751 current_state%global_grid%configuration%vertical%rdzn(k+1) * &
755 (current_state%th%data(k+1,jcol,icol) + &
756 current_state%global_grid%configuration%vertical%thref(k+1) + &
758 current_state%global_grid%configuration%vertical%thref(k+1) * &
759 (c_virtual * current_state%q(1)%data(k+1,jcol,icol) &
760 - current_state%q(2)%data(k+1,jcol,icol))) - &
763 (current_state%th%data(k,jcol,icol) + &
764 current_state%global_grid%configuration%vertical%thref(k) + &
765 current_state%global_grid%configuration%vertical%thref(k) * &
766 (c_virtual * current_state%q(1)%data(k,jcol,icol) &
767 - current_state%q(2)%data(k,jcol,icol))) )
773 current_state%vis_coefficient%data(1, jcol, icol) &
775 (current_state%global_grid%configuration%vertical%zn(2) * &
776 current_state%global_grid%configuration%vertical%zn(2)) + buoy_prod
778 call log_master_log(log_error, &
779 "Subgrid diags - buoy_prod calc needs theta and q active, STOP")
786 subke_2d(target_y_index, target_x_index) =
subke_2d(target_y_index, target_x_index)/ &
792 if (current_state%th%active)
then
793 epsth=calculate_thermal_dissipation_rate(current_state, current_state%th)
796 do k=2, current_state%local_grid%size(z_index)-1
809 do k=2,current_state%local_grid%size(z_index)-1
813 umean(k)=(current_state%global_grid%configuration%vertical%olubar(k+1) - &
814 current_state%global_grid%configuration%vertical%olubar(k)) * &
815 current_state%global_grid%configuration%vertical%rdzn(k+1)
816 vmean(k)=(current_state%global_grid%configuration%vertical%olvbar(k+1) - &
817 current_state%global_grid%configuration%vertical%olvbar(k)) * &
818 current_state%global_grid%configuration%vertical%rdzn(k+1)
820 sg_shear_prod = current_state%vis_coefficient%data(k,jcol,icol)* &
821 (s13(k)*umean(k)+s23(k)*vmean(k))
826 uwsg(k)=0.5*(tau13(k)+tau13(k+1))/rho(k)
827 vwsg(k)=0.5*(tau23(k)+tau23(k+1))/rho(k)
833 sg_shear_prod = current_state%vis_coefficient%data(1, jcol, icol) * &
835 (current_state%u%data(2,jcol,icol-1) + &
836 current_state%u%data(2,jcol,icol) ) * &
837 current_state%global_grid%configuration%vertical%olubar(2) + &
839 (current_state%v%data(2,jcol-1,icol) + &
840 current_state%v%data(2,jcol,icol) ) * &
841 current_state%global_grid%configuration%vertical%olvbar(2) ) / &
842 (current_state%global_grid%configuration%vertical%zn(2) * &
843 current_state%global_grid%configuration%vertical%zn(2))
851 do k=2, current_state%local_grid%size(z_index)
853 u_i_prime_tau_i(k) = ( 0.5_default_precision * &
854 (current_state%u%data(k,jcol,icol-1) + &
855 current_state%u%data(k,jcol,icol) ) - &
856 current_state%global_grid%configuration%vertical%olubar(k)) * &
857 0.5_default_precision * ( tau13(k-1)+tau13(k))
859 u_i_prime_tau_i(k) = u_i_prime_tau_i(k) + ( 0.5_default_precision * &
860 (current_state%v%data(k,jcol-1,icol) + &
861 current_state%v%data(k,jcol,icol) ) - &
862 current_state%global_grid%configuration%vertical%olvbar(k)) * &
863 0.5_default_precision * ( tau23(k-1)+tau23(k))
865 u_i_prime_tau_i(k) = u_i_prime_tau_i(k) + ( 0.5_default_precision * &
866 (current_state%w%data(k-1,jcol,icol)+&
867 current_state%w%data(k,jcol,icol))) *&
872 u_i_prime_tau_i(1) = &
873 ( 0.5_default_precision * &
874 (current_state%u%data(2,jcol,icol-1) + &
875 current_state%u%data(2,jcol,icol) ) - &
876 current_state%global_grid%configuration%vertical%olubar(2) ) * &
877 ( 0.5_default_precision * &
878 (current_state%u%data(2,jcol,icol-1) + &
879 current_state%u%data(2,jcol,icol) ) )
881 u_i_prime_tau_i(1) = u_i_prime_tau_i(1) + &
882 ( 0.5_default_precision * &
883 (current_state%v%data(2,jcol-1,icol) + &
884 current_state%v%data(2,jcol,icol) ) - &
885 current_state%global_grid%configuration%vertical%olvbar(2) ) * &
886 ( 0.5_default_precision * &
887 (current_state%v%data(2,jcol-1,icol) + &
888 current_state%v%data(2,jcol,icol) ) )
892 u_i_prime_tau_i(1) = u_i_prime_tau_i(1) * &
893 current_state%vis_coefficient%data(1,jcol,icol) / &
894 current_state%global_grid%configuration%vertical%zn(2)
897 do k=2, current_state%local_grid%size(z_index)-1
898 sed_tot(k)=
sed_tot(k) + (u_i_prime_tau_i(k+1)-u_i_prime_tau_i(k)) * &
899 current_state%global_grid%configuration%vertical%rdzn(k+1) / &
900 current_state%global_grid%configuration%vertical%rho(k)
904 current_state%global_grid%configuration%vertical%zn(2)
916 type(model_state_type),
target,
intent(inout) :: current_state
917 character(len=*),
intent(in) :: name
918 type(component_field_information_type),
intent(out) :: field_information
920 field_information%field_type=component_array_field_type
921 field_information%data_type=component_double_data_type
922 if (name .eq.
'subke')
then
923 field_information%number_dimensions=2
924 field_information%dimension_sizes(1)=current_state%local_grid%size(y_index)
925 field_information%dimension_sizes(2)=current_state%local_grid%size(x_index)
927 field_information%number_dimensions=1
928 field_information%dimension_sizes(1)=current_state%local_grid%size(z_index)
931 if (name .eq.
"th2sg_total_local" .or. name .eq.
"wtsg_total_local" &
932 .or. name .eq.
"theta_dis_total_local" )
then
933 field_information%enabled=current_state%th%active
934 else if (name .eq.
"wqv_sg_total_local" .or. name .eq.
"wql_sg_total_local")
then
935 field_information%enabled= (.not.current_state%passive_q) .and. &
936 (current_state%number_q_fields .gt. 0)
937 else if (name .eq.
"wqr_sg_total_local")
then
938 field_information%enabled= current_state%rain_water_mixing_ratio_index .gt. 0
939 else if (name .eq.
"wqi_sg_total_local")
then
940 field_information%enabled= current_state%ice_water_mixing_ratio_index .gt. 0
941 else if (name .eq.
"wqs_sg_total_local")
then
942 field_information%enabled= current_state%snow_water_mixing_ratio_index .gt. 0
943 else if (name .eq.
"wqg_sg_total_local")
then
944 field_information%enabled= current_state%graupel_water_mixing_ratio_index .gt. 0
948 else if (name .eq.
"i_th2sg_total_local")
then
949 field_information%enabled=current_state%th%active
950 else if (name .eq.
"i_wqsg_total_local")
then
951 field_information%enabled= (.not.current_state%passive_q) .and. &
952 (current_state%liquid_water_mixing_ratio_index > 0)
955 field_information%enabled=.true.
964 type(model_state_type),
target,
intent(inout) :: current_state
965 character(len=*),
intent(in) :: name
966 type(component_field_value_type),
intent(out) :: field_value
970 if (name .eq.
"uwsg_total_local")
then
971 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
972 do k = 1, current_state%local_grid%size(z_index)
973 field_value%real_1d_array(k)=
uwsg_tot(k)
975 else if (name .eq.
"sed_total_local")
then
976 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
977 do k = 1, current_state%local_grid%size(z_index)
978 field_value%real_1d_array(k)=
sed_tot(k)
980 else if (name .eq.
"ssub_total_local")
then
981 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
982 do k = 1, current_state%local_grid%size(z_index)
983 field_value%real_1d_array(k)=
ssub_tot(k)
985 else if (name .eq.
"buoysg_total_local")
then
986 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
987 do k = 1, current_state%local_grid%size(z_index)
990 else if (name .eq.
"vwsg_total_local")
then
991 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
992 do k = 1, current_state%local_grid%size(z_index)
993 field_value%real_1d_array(k)=
vwsg_tot(k)
995 else if (name .eq.
"uusg_total_local")
then
996 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
997 do k = 1, current_state%local_grid%size(z_index)
998 field_value%real_1d_array(k)=
uusg_tot(k)
1000 else if (name .eq.
"vvsg_total_local")
then
1001 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1002 do k = 1, current_state%local_grid%size(z_index)
1003 field_value%real_1d_array(k)=
vvsg_tot(k)
1005 else if (name .eq.
"wwsg_total_local")
then
1006 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1007 do k = 1, current_state%local_grid%size(z_index)
1008 field_value%real_1d_array(k)=
wwsg_tot(k)
1010 else if (name .eq.
"tkesg_total_local")
then
1011 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1012 do k = 1, current_state%local_grid%size(z_index)
1013 field_value%real_1d_array(k)=
tkesg_tot(k)
1015 else if (name .eq.
"wkesg_total_local")
then
1016 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1017 do k = 1, current_state%local_grid%size(z_index)
1018 field_value%real_1d_array(k)=
wkesg_tot(k)
1020 else if (name .eq.
"viscosity_coef_total_local")
then
1021 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1022 do k = 1, current_state%local_grid%size(z_index)
1025 else if (name .eq.
"diffusion_coef_total_local")
then
1026 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1027 do k = 1, current_state%local_grid%size(z_index)
1030 else if (name .eq.
"richardson_number_total_local")
then
1031 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1032 do k = 1, current_state%local_grid%size(z_index)
1035 else if (name .eq.
"richardson_squared_total_local")
then
1036 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1037 do k = 1, current_state%local_grid%size(z_index)
1040 else if (name .eq.
"dissipation_total_local")
then
1041 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1042 do k = 1, current_state%local_grid%size(z_index)
1045 else if (name .eq.
"wtsg_total_local")
then
1046 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1047 do k = 1, current_state%local_grid%size(z_index)
1048 field_value%real_1d_array(k)=
wtsg_tot(k)
1050 else if (name .eq.
"th2sg_total_local")
then
1051 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1052 do k = 1, current_state%local_grid%size(z_index)
1053 field_value%real_1d_array(k)=
th2sg_tot(k)
1055 else if (name .eq.
"theta_dis_total_local")
then
1056 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1057 do k = 1, current_state%local_grid%size(z_index)
1060 else if (name .eq.
"wqv_sg_total_local")
then
1061 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1062 do k = 1, current_state%local_grid%size(z_index)
1065 else if (name .eq.
"wql_sg_total_local")
then
1066 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1067 do k = 1, current_state%local_grid%size(z_index)
1070 else if (name .eq.
"wqr_sg_total_local")
then
1071 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1072 do k = 1, current_state%local_grid%size(z_index)
1075 else if (name .eq.
"wqi_sg_total_local")
then
1076 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1077 do k = 1, current_state%local_grid%size(z_index)
1080 else if (name .eq.
"wqs_sg_total_local")
then
1081 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1082 do k = 1, current_state%local_grid%size(z_index)
1085 else if (name .eq.
"wqg_sg_total_local")
then
1086 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1087 do k = 1, current_state%local_grid%size(z_index)
1092 else if (name .eq.
"i_uwsg_total_local")
then
1093 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1094 do k = 1, current_state%local_grid%size(z_index)
1095 field_value%real_1d_array(k)=
uwsg_tot(k)
1097 else if (name .eq.
"i_vwsg_total_local")
then
1098 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1099 do k = 1, current_state%local_grid%size(z_index)
1100 field_value%real_1d_array(k)=
vwsg_tot(k)
1102 else if (name .eq.
"i_uusg_total_local")
then
1103 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1104 do k = 1, current_state%local_grid%size(z_index)
1105 field_value%real_1d_array(k)=
uusg_tot(k)
1107 else if (name .eq.
"i_vvsg_total_local")
then
1108 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1109 do k = 1, current_state%local_grid%size(z_index)
1110 field_value%real_1d_array(k)=
vvsg_tot(k)
1112 else if (name .eq.
"i_wwsg_total_local")
then
1113 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1114 do k = 1, current_state%local_grid%size(z_index)
1115 field_value%real_1d_array(k)=
wwsg_tot(k)
1117 else if (name .eq.
"i_tkesg_total_local")
then
1118 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1119 do k = 1, current_state%local_grid%size(z_index)
1120 field_value%real_1d_array(k)=
tkesg_tot(k)
1122 else if (name .eq.
"i_wtsg_total_local")
then
1123 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1124 do k = 1, current_state%local_grid%size(z_index)
1125 field_value%real_1d_array(k)=
wtsg_tot(k)
1127 else if (name .eq.
"i_th2sg_total_local")
then
1128 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1129 do k = 1, current_state%local_grid%size(z_index)
1130 field_value%real_1d_array(k)=
th2sg_tot(k)
1132 else if (name .eq.
"i_wqsg_total_local")
then
1133 allocate(field_value%real_1d_array(current_state%local_grid%size(z_index)))
1134 do k = 1, current_state%local_grid%size(z_index)
1135 field_value%real_1d_array(k)=
wqsg_tot(k)
1138 else if (name .eq.
'subke')
then
1139 allocate(field_value%real_2d_array(current_state%local_grid%size(y_index), &
1140 current_state%local_grid%size(x_index)))
1141 field_value%real_2d_array(:,:)=
subke_2d(:,:)