MONC
Functions/Subroutines | Variables
smagorinsky_mod Module Reference

Calculates the Smagorinsky eddy viscosity and diffusivity at l_w-points. More...

Functions/Subroutines

type(component_descriptor_type) function, public smagorinsky_get_descriptor ()
 Provides the descriptor back to the caller and is used in component registration. More...
 
subroutine initialisation_callback (current_state)
 Initialisation call back which will read in the coriolis configuration and set up the geostrophic winds. More...
 
subroutine timestep_callback (current_state)
 For each none halo cell this will calculate the subgrid terms for viscosity and diffusion. More...
 
subroutine finalisation_callback (current_state)
 Called when the model is finishing up, will finalise the halo communications represented by the state. More...
 
subroutine update_viscous_number (current_state)
 Update viscous number based upon the viscosity and diffusivity coefficients. More...
 
subroutine setfri (current_state, richardson_number, ssq)
 Calculates the eddy viscosity (VIS) and diffusivity (DIFF) depending on the Richardson Number (RI) and half squared strain rate. More...
 
real(kind=default_precision) function, dimension(current_state%local_grid%size(z_index)), public calculate_richardson_number (current_state, ssq, th, q)
 Calculates the richardson number depending upon the setup of the model and the method selected. More...
 
real(kind=default_precision) function, dimension(current_state%local_grid%size(z_index)) moist_ri_2 (current_state, ssq, th, q)
 Calculates another "moist version" of the Richardson number based on "change in subgrid buoyancy flux when a fraction EPS is exchanged". More...
 
real(kind=default_precision) function, dimension(current_state%local_grid%size(z_index)) moist_ri_1 (current_state, ssq, th, q)
 Calculates numerator of Richarson number as the difference in buoyancy between a parcel at K+1 and a parcel lifted from K to K+1. Condensation is calculated by assuming that TL=T-(Lv/CP)*QL+gz/CP and QT = QV+QL are conserved on lifting. Microphysics conversions are assumed to happen slowly compared to turbulent processes so are, therefore, neglected in the calculation. More...
 
real(kind=default_precision) function, dimension(current_state%local_grid%size(z_index)), public calculate_half_squared_strain_rate (current_state, u, v, w)
 Calculates the half squared strain rate on l_w-points which is used in determining the Richardson number. CX=1./DX, CY=1./DY, RDZ(K)=1./DZ(K), RDZN(K) =1./DZN(K) _SSQ= 0.5*^^DU_I/DX_J+DU_J/DX_I^^**2 _SSQIJ= (DU_I/DX_J+DU_J/DX_I)**2 _Hence SSQ= SUM(I,J) {0.5*(SSQIJ)}. More...
 
real(kind=default_precision) function, dimension(current_state%local_grid%size(z_index)), public calculate_thermal_dissipation_rate (current_state, th)
 
subroutine copy_vis_to_halo_buffer (current_state, neighbour_description, dim, source_index, pid_location, current_page, source_data)
 Copies the viscosity field data to halo buffers for a specific process in a dimension and halo cell. More...
 
subroutine copy_vis_corners_to_halo_buffer (current_state, neighbour_description, corner_loc, x_source_index, y_source_index, pid_location, current_page, source_data)
 Copies the viscosity field corner data to halo buffers for a specific process. More...
 
subroutine copy_diff_to_halo_buffer (current_state, neighbour_description, dim, source_index, pid_location, current_page, source_data)
 Copies the diffusion field data to halo buffers for a specific process in a dimension and halo cell. More...
 
subroutine copy_diff_corners_to_halo_buffer (current_state, neighbour_description, corner_loc, x_source_index, y_source_index, pid_location, current_page, source_data)
 Copies the diffusion field corner data to halo buffers for a specific process. More...
 

Variables

integer, parameter richardson_number_calculation =2
 
real(kind=default_precision) eps
 
real(kind=default_precision) repsh
 
real(kind=default_precision) thcona
 
real(kind=default_precision) thconb
 
real(kind=default_precision) thconap1
 
real(kind=default_precision) suba
 
real(kind=default_precision) subb
 
real(kind=default_precision) subc
 
real(kind=default_precision) subg
 
real(kind=default_precision) subh
 
real(kind=default_precision) subr
 
real(kind=default_precision) pr_n
 
real(kind=default_precision) ric
 
real(kind=default_precision) ricinv
 
integer iqv
 
integer iql
 

Detailed Description

Calculates the Smagorinsky eddy viscosity and diffusivity at l_w-points.

Function/Subroutine Documentation

◆ calculate_half_squared_strain_rate()

real(kind=default_precision) function, dimension(current_state%local_grid%size(z_index)), public smagorinsky_mod::calculate_half_squared_strain_rate ( type(model_state_type), intent(inout), target  current_state,
type(prognostic_field_type), intent(inout)  u,
type(prognostic_field_type), intent(inout)  v,
type(prognostic_field_type), intent(inout)  w 
)

Calculates the half squared strain rate on l_w-points which is used in determining the Richardson number. CX=1./DX, CY=1./DY, RDZ(K)=1./DZ(K), RDZN(K) =1./DZN(K) _SSQ= 0.5*^^DU_I/DX_J+DU_J/DX_I^^**2 _SSQIJ= (DU_I/DX_J+DU_J/DX_I)**2 _Hence SSQ= SUM(I,J) {0.5*(SSQIJ)}.

Parameters
current_stateThe current model state
Returns
The hard squared strain rate for a specific column

Definition at line 409 of file smagorinsky.F90.

409  type(model_state_type), target, intent(inout) :: current_state
410  type(prognostic_field_type), intent(inout) :: u, v, w
411  real(kind=DEFAULT_PRECISION) :: calculate_half_squared_strain_rate(current_state%local_grid%size(z_index))
412 
413  integer :: k
414  real(kind=DEFAULT_PRECISION) :: ssq11, ssq22, ssq33, ssq13, ssq23, ssq12
415 
416  do k=2,current_state%local_grid%size(z_index)-1
417 #ifdef U_ACTIVE
418  ssq11=current_state%global_grid%configuration%horizontal%cx2*(&
419  (u%data(k+1, current_state%column_local_y, current_state%column_local_x)-&
420  u%data(k+1, current_state%column_local_y, current_state%column_local_x-1))**2+&
421  (u%data(k, current_state%column_local_y, current_state%column_local_x)-&
422  u%data(k, current_state%column_local_y, current_state%column_local_x-1))**2)
423 #else
424  ssq11=0.0_default_precision
425 #endif
426 #ifdef V_ACTIVE
427  ssq22=current_state%global_grid%configuration%horizontal%cy2*(&
428  (v%data(k+1, current_state%column_local_y, current_state%column_local_x)-&
429  v%data(k+1, current_state%column_local_y-1, current_state%column_local_x))**2+&
430  (v%data(k, current_state%column_local_y, current_state%column_local_x)-&
431  v%data(k, current_state%column_local_y-1, current_state%column_local_x))**2)
432 #else
433  ssq22=0.0_default_precision
434 #endif
435 #ifdef W_ACTIVE
436  ssq33=((w%data(k, current_state%column_local_y, current_state%column_local_x)-&
437  w%data(k-1, current_state%column_local_y, current_state%column_local_x))*&
438  current_state%global_grid%configuration%vertical%rdz(k))**2 +&
439  ((w%data(k+1, current_state%column_local_y, current_state%column_local_x)-&
440  w%data(k, current_state%column_local_y, current_state%column_local_x))*&
441  current_state%global_grid%configuration%vertical%rdz(k+1))**2
442 #else
443  ssq33=0.0_default_precision
444 #endif
445 #if defined(U_ACTIVE) && defined(W_ACTIVE)
446  ! Average over 2 points U and W
447  ssq13=(((u%data(k+1, current_state%column_local_y, current_state%column_local_x)-&
448  u%data(k, current_state%column_local_y, current_state%column_local_x))*&
449  current_state%global_grid%configuration%vertical%rdzn(k+1)+&
450  (w%data(k, current_state%column_local_y, current_state%column_local_x+1)-&
451  w%data(k, current_state%column_local_y, current_state%column_local_x))*&
452  current_state%global_grid%configuration%horizontal%cx)**2+&
453  ((u%data(k+1, current_state%column_local_y, current_state%column_local_x-1)-&
454  u%data(k, current_state%column_local_y, current_state%column_local_x-1))*&
455  current_state%global_grid%configuration%vertical%rdzn(k+1)+&
456  (w%data(k, current_state%column_local_y, current_state%column_local_x)-&
457  w%data(k, current_state%column_local_y, current_state%column_local_x-1))*&
458  current_state%global_grid%configuration%horizontal%cx)**2)*0.5_default_precision
459 #elif defined(U_ACTIVE) && !defined(W_ACTIVE)
460  ssq13=(((u%data(k+1, current_state%column_local_y, current_state%column_local_x)-&
461  u%data(k, current_state%column_local_y, current_state%column_local_x))*&
462  current_state%global_grid%configuration%vertical%rdzn(k+1))**2+&
463  ((u%data(k+1, current_state%column_local_y, current_state%column_local_x-1)-&
464  u%data(k, current_state%column_local_y, current_state%column_local_x-1))*&
465  current_state%global_grid%configuration%vertical%rdzn(k+1))**2)*0.5_default_precision
466 #elif !defined(U_ACTIVE) && defined(W_ACTIVE)
467  ssq13=(((w%data(k, current_state%column_local_y, current_state%column_local_x+1)-&
468  w%data(k, current_state%column_local_y, current_state%column_local_x))*&
469  current_state%global_grid%configuration%horizontal%cx)**2+&
470  ((w%data(k, current_state%column_local_y, current_state%column_local_x)-&
471  w%data(k, current_state%column_local_y, current_state%column_local_x-1))*&
472  current_state%global_grid%configuration%horizontal%cx)**2)*0.5_default_precision
473 #else
474  ssq13=0.0_default_precision
475 #endif
476 #if defined(W_ACTIVE) && defined(V_ACTIVE)
477  ! Average over 2 points W and V
478  ssq23=(((w%data(k, current_state%column_local_y, current_state%column_local_x)-&
479  w%data(k, current_state%column_local_y-1, current_state%column_local_x))*&
480  current_state%global_grid%configuration%horizontal%cy+&
481  (v%data(k+1, current_state%column_local_y-1, current_state%column_local_x)-&
482  v%data(k, current_state%column_local_y-1, current_state%column_local_x))*&
483  current_state%global_grid%configuration%vertical%rdzn(k+1))**2+&
484  ((w%data(k, current_state%column_local_y+1, current_state%column_local_x)-&
485  w%data(k, current_state%column_local_y, current_state%column_local_x))*&
486  current_state%global_grid%configuration%horizontal%cy+&
487  (v%data(k+1, current_state%column_local_y, current_state%column_local_x)-&
488  v%data(k, current_state%column_local_y, current_state%column_local_x))*&
489  current_state%global_grid%configuration%vertical%rdzn(k+1))**2)*0.5_default_precision
490 #elif defined(W_ACTIVE) && !defined(V_ACTIVE)
491  ssq23=(((w%data(k, current_state%column_local_y, current_state%column_local_x)-&
492  w%data(k, current_state%column_local_y-1, current_state%column_local_x))*&
493  current_state%global_grid%configuration%horizontal%cy)**2+&
494  ((w%data(k, current_state%column_local_y+1, current_state%column_local_x)-&
495  w%data(k, current_state%column_local_y, current_state%column_local_x))*&
496  current_state%global_grid%configuration%horizontal%cy)**2)*0.5_default_precision
497 #elif !defined(W_ACTIVE) && defined(V_ACTIVE)
498  ssq23=(((v%data(k+1, current_state%column_local_y-1, current_state%column_local_x)-&
499  v%data(k, current_state%column_local_y-1, current_state%column_local_x))*&
500  current_state%global_grid%configuration%vertical%rdzn(k+1))**2+&
501  ((v%data(k+1, current_state%column_local_y, current_state%column_local_x)-&
502  v%data(k, current_state%column_local_y, current_state%column_local_x))*&
503  current_state%global_grid%configuration%vertical%rdzn(k+1))**2)*0.5_default_precision
504 #else
505  ssq23=0.0_default_precision
506 #endif
507 
508 #if defined(U_ACTIVE) && defined(V_ACTIVE)
509  ! Average over 8 points from U and V
510  ssq12=(((((u%data(k, current_state%column_local_y, current_state%column_local_x-1)-&
511  u%data(k, current_state%column_local_y-1, current_state%column_local_x-1))*&
512  current_state%global_grid%configuration%horizontal%cy+&
513  (v%data(k, current_state%column_local_y-1, current_state%column_local_x)-&
514  v%data(k, current_state%column_local_y-1, current_state%column_local_x-1))*&
515  current_state%global_grid%configuration%horizontal%cx)**2 +&
516  ((u%data(k, current_state%column_local_y+1, current_state%column_local_x-1)-&
517  u%data(k, current_state%column_local_y, current_state%column_local_x-1))*&
518  current_state%global_grid%configuration%horizontal%cy+&
519  (v%data(k, current_state%column_local_y, current_state%column_local_x)-&
520  v%data(k, current_state%column_local_y, current_state%column_local_x-1))*&
521  current_state%global_grid%configuration%horizontal%cx)**2) +(&
522  ((u%data(k, current_state%column_local_y, current_state%column_local_x)-&
523  u%data(k, current_state%column_local_y-1, current_state%column_local_x))*&
524  current_state%global_grid%configuration%horizontal%cy+&
525  (v%data(k, current_state%column_local_y-1, current_state%column_local_x+1)-&
526  v%data(k, current_state%column_local_y-1, current_state%column_local_x))*&
527  current_state%global_grid%configuration%horizontal%cx)**2 +&
528  ((u%data(k, current_state%column_local_y+1, current_state%column_local_x)-&
529  u%data(k, current_state%column_local_y, current_state%column_local_x))*&
530  current_state%global_grid%configuration%horizontal%cy+&
531  (v%data(k, current_state%column_local_y, current_state%column_local_x+1)-&
532  v%data(k, current_state%column_local_y, current_state%column_local_x))*&
533  current_state%global_grid%configuration%horizontal%cx)**2))+((&
534  ((u%data(k+1, current_state%column_local_y, current_state%column_local_x-1)-&
535  u%data(k+1, current_state%column_local_y-1, current_state%column_local_x-1))*&
536  current_state%global_grid%configuration%horizontal%cy+&
537  (v%data(k+1, current_state%column_local_y-1, current_state%column_local_x)-&
538  v%data(k+1, current_state%column_local_y-1, current_state%column_local_x-1))*&
539  current_state%global_grid%configuration%horizontal%cx)**2+&
540  ((u%data(k+1, current_state%column_local_y+1, current_state%column_local_x-1)-&
541  u%data(k+1, current_state%column_local_y, current_state%column_local_x-1))*&
542  current_state%global_grid%configuration%horizontal%cy+&
543  (v%data(k+1, current_state%column_local_y, current_state%column_local_x)-&
544  v%data(k+1, current_state%column_local_y, current_state%column_local_x-1))*&
545  current_state%global_grid%configuration%horizontal%cx)**2)+(&
546  ((u%data(k+1, current_state%column_local_y, current_state%column_local_x)-&
547  u%data(k+1, current_state%column_local_y-1, current_state%column_local_x))*&
548  current_state%global_grid%configuration%horizontal%cy+&
549  (v%data(k+1, current_state%column_local_y-1, current_state%column_local_x+1)-&
550  v%data(k+1, current_state%column_local_y-1, current_state%column_local_x))*&
551  current_state%global_grid%configuration%horizontal%cx)**2+&
552  ((u%data(k+1, current_state%column_local_y+1, current_state%column_local_x)-&
553  u%data(k+1, current_state%column_local_y, current_state%column_local_x))*&
554  current_state%global_grid%configuration%horizontal%cy+&
555  (v%data(k+1, current_state%column_local_y, current_state%column_local_x+1)-&
556  v%data(k+1, current_state%column_local_y, current_state%column_local_x))*&
557  current_state%global_grid%configuration%horizontal%cx)**2)))*0.125_default_precision
558 
559 #elif defined(U_ACTIVE) && !defined(V_ACTIVE)
560 
561  ssq12=(((((u%data(k, current_state%column_local_y, current_state%column_local_x-1)-&
562  u%data(k, current_state%column_local_y-1, current_state%column_local_x-1))*&
563  current_state%global_grid%configuration%horizontal%cy)**2 +&
564  ((u%data(k, current_state%column_local_y+1, current_state%column_local_x-1)-&
565  u%data(k, current_state%column_local_y, current_state%column_local_x-1))*&
566  current_state%global_grid%configuration%horizontal%cy)**2) +(&
567  ((u%data(k, current_state%column_local_y, current_state%column_local_x)-&
568  u%data(k, current_state%column_local_y-1, current_state%column_local_x))*&
569  current_state%global_grid%configuration%horizontal%cy)**2 +&
570  ((u%data(k, current_state%column_local_y+1, current_state%column_local_x)-&
571  u%data(k, current_state%column_local_y, current_state%column_local_x))*&
572  current_state%global_grid%configuration%horizontal%cy)**2))+((&
573  ((u%data(k+1, current_state%column_local_y, current_state%column_local_x-1)-&
574  u%data(k+1, current_state%column_local_y-1, current_state%column_local_x-1))*&
575  current_state%global_grid%configuration%horizontal%cy)**2+&
576  ((u%data(k+1, current_state%column_local_y+1, current_state%column_local_x-1)-&
577  u%data(k+1, current_state%column_local_y, current_state%column_local_x-1))*&
578  current_state%global_grid%configuration%horizontal%cy)**2)+(&
579  ((u%data(k+1, current_state%column_local_y, current_state%column_local_x)-&
580  u%data(k+1, current_state%column_local_y-1, current_state%column_local_x))*&
581  current_state%global_grid%configuration%horizontal%cy)**2+&
582  ((u%data(k+1, current_state%column_local_y+1, current_state%column_local_x)-&
583  u%data(k+1, current_state%column_local_y, current_state%column_local_x))*&
584  current_state%global_grid%configuration%horizontal%cy)**2)))*0.125_default_precision
585 
586 #elif !defined(U_ACTIVE) && defined(V_ACTIVE)
587 
588  ssq12=(((((v%data(k, current_state%column_local_y-1, current_state%column_local_x)-&
589  v%data(k, current_state%column_local_y-1, current_state%column_local_x-1))*&
590  current_state%global_grid%configuration%horizontal%cx)**2 +&
591  ((v%data(k, current_state%column_local_y, current_state%column_local_x)-&
592  v%data(k, current_state%column_local_y, current_state%column_local_x-1))*&
593  current_state%global_grid%configuration%horizontal%cx)**2) +&
594  (((v%data(k, current_state%column_local_y-1, current_state%column_local_x+1)-&
595  v%data(k, current_state%column_local_y-1, current_state%column_local_x))*&
596  current_state%global_grid%configuration%horizontal%cx)**2 +&
597  ((v%data(k, current_state%column_local_y, current_state%column_local_x+1)-&
598  v%data(k, current_state%column_local_y, current_state%column_local_x))*&
599  current_state%global_grid%configuration%horizontal%cx)**2))+((&
600  ((v%data(k+1, current_state%column_local_y-1, current_state%column_local_x)-&
601  v%data(k+1, current_state%column_local_y-1, current_state%column_local_x-1))*&
602  current_state%global_grid%configuration%horizontal%cx)**2+&
603  ((v%data(k+1, current_state%column_local_y, current_state%column_local_x)-&
604  v%data(k+1, current_state%column_local_y, current_state%column_local_x-1))*&
605  current_state%global_grid%configuration%horizontal%cx)**2)+(&
606  ((v%data(k+1, current_state%column_local_y-1, current_state%column_local_x+1)-&
607  v%data(k+1, current_state%column_local_y-1, current_state%column_local_x))*&
608  current_state%global_grid%configuration%horizontal%cx)**2+&
609  ((v%data(k+1, current_state%column_local_y, current_state%column_local_x+1)-&
610  v%data(k+1, current_state%column_local_y, current_state%column_local_x))*&
611  current_state%global_grid%configuration%horizontal%cx)**2)))*0.125_default_precision
612 #else
613  ssq12=0.0_default_precision
614 #endif
615  calculate_half_squared_strain_rate(k)=ssq11+ssq22+ssq33+ssq13+ssq23+ssq12+smallp
616  end do
Here is the caller graph for this function:

◆ calculate_richardson_number()

real(kind=default_precision) function, dimension(current_state%local_grid%size(z_index)), public smagorinsky_mod::calculate_richardson_number ( type(model_state_type), intent(inout), target  current_state,
real(kind=default_precision), dimension(:), intent(in)  ssq,
type(prognostic_field_type), intent(inout)  th,
type(prognostic_field_type), dimension(:), intent(inout)  q 
)

Calculates the richardson number depending upon the setup of the model and the method selected.

Parameters
current_stateThe current model state
ssqThe half squared strain rate
Returns
The Richardson number for a specific column

Definition at line 203 of file smagorinsky.F90.

203  type(model_state_type), target, intent(inout) :: current_state
204  real(kind=DEFAULT_PRECISION), dimension(:), intent(in) :: ssq
205  type(prognostic_field_type), intent(inout) :: th
206  type(prognostic_field_type), dimension(:), intent(inout) :: q
207  real(kind=DEFAULT_PRECISION), dimension(current_state%local_grid%size(Z_INDEX)) :: calculate_richardson_number
208 
209  integer :: k
210 
211  if (.not. current_state%passive_th) then
212  if (.not. current_state%passive_q) then
213  if (richardson_number_calculation .eq. 2) then
214  calculate_richardson_number=moist_ri_2(current_state, ssq, th, q)
215  else
216  calculate_richardson_number=moist_ri_1(current_state, ssq, th, q)
217  end if
218  else
219  do k=2, current_state%local_grid%size(z_index)-1
220  calculate_richardson_number(k) = (current_state%global_grid%configuration%vertical%dthref(k) + &
221  current_state%th%data(k+1, current_state%column_local_y, current_state%column_local_x) -&
222  current_state%th%data(k, current_state%column_local_y, current_state%column_local_x))* &
223  current_state%global_grid%configuration%vertical%rdzn(k+1)*&
224  current_state%global_grid%configuration%vertical%buoy_co(k) / ssq(k)
225  end do
226  end if
227  else
228  calculate_richardson_number=0.0_default_precision
229  end if
Here is the call graph for this function:
Here is the caller graph for this function:

◆ calculate_thermal_dissipation_rate()

real(kind=default_precision) function, dimension(current_state%local_grid%size(z_index)), public smagorinsky_mod::calculate_thermal_dissipation_rate ( type(model_state_type), intent(inout), target  current_state,
type(prognostic_field_type), intent(inout)  th 
)

Definition at line 623 of file smagorinsky.F90.

623  type(model_state_type), target, intent(inout) :: current_state
624  type(prognostic_field_type), intent(inout) :: th
625  real(kind=DEFAULT_PRECISION) :: calculate_thermal_dissipation_rate(current_state%local_grid%size(z_index))
626 
627  integer :: k
628 
629  do k=2,current_state%local_grid%size(z_index)-1
630  calculate_thermal_dissipation_rate(k) = &
631  current_state%diff_coefficient%data(k, current_state%column_local_y, current_state%column_local_x) * ( &
632  (current_state%global_grid%configuration%vertical%rdzn(k+1) * &
633  (current_state%th%data(k+1, current_state%column_local_y, current_state%column_local_x) - &
634  current_state%th%data(k, current_state%column_local_y, current_state%column_local_x) - &
635  current_state%global_grid%configuration%vertical%dthref(k) ) ) ** 2 + &
636  0.25 * current_state%global_grid%configuration%horizontal%cx2 * &
637  ( (current_state%th%data(k, current_state%column_local_y, current_state%column_local_x+1) - &
638  current_state%th%data(k, current_state%column_local_y, current_state%column_local_x) ) ** 2 + &
639  (current_state%th%data(k, current_state%column_local_y, current_state%column_local_x) - &
640  current_state%th%data(k, current_state%column_local_y, current_state%column_local_x-1) ) ** 2 + &
641  (current_state%th%data(k+1, current_state%column_local_y, current_state%column_local_x+1) - &
642  current_state%th%data(k+1, current_state%column_local_y, current_state%column_local_x) ) **2 + &
643  (current_state%th%data(k+1, current_state%column_local_y, current_state%column_local_x) - &
644  current_state%th%data(k+1, current_state%column_local_y, current_state%column_local_x-1) ) ** 2 ) + &
645  0.25 * current_state%global_grid%configuration%horizontal%cy2 * &
646  ( (current_state%th%data(k, current_state%column_local_y+1, current_state%column_local_x) - &
647  current_state%th%data(k, current_state%column_local_y, current_state%column_local_x) ) ** 2 + &
648  (current_state%th%data(k, current_state%column_local_y, current_state%column_local_x) - &
649  current_state%th%data(k, current_state%column_local_y-1, current_state%column_local_x) ) ** 2 + &
650  (current_state%th%data(k+1, current_state%column_local_y+1, current_state%column_local_x) - &
651  current_state%th%data(k+1, current_state%column_local_y, current_state%column_local_x) ) **2 + &
652  (current_state%th%data(k+1, current_state%column_local_y, current_state%column_local_x) - &
653  current_state%th%data(k+1, current_state%column_local_y-1, current_state%column_local_x) ) ** 2 ) )
654  end do
655 

◆ copy_diff_corners_to_halo_buffer()

subroutine smagorinsky_mod::copy_diff_corners_to_halo_buffer ( type(model_state_type), intent(inout)  current_state,
type(neighbour_description_type), intent(inout)  neighbour_description,
integer, intent(in)  corner_loc,
integer, intent(in)  x_source_index,
integer, intent(in)  y_source_index,
integer, intent(in)  pid_location,
integer, dimension(:), intent(inout)  current_page,
type(field_data_wrapper_type), dimension(:), intent(in), optional  source_data 
)
private

Copies the diffusion field corner data to halo buffers for a specific process.

Parameters
current_stateThe current model state
neighbour_descriptionsDescription of the neighbour halo swapping status
corner_locThe corner location
x_source_indexThe X source index of the dimension we are reading from in the prognostic field
y_source_indexThe Y source index of the dimension we are reading from in the prognostic field
pid_locationLocation of the neighbouring process in the local stored data structures
current_pageThe current (next) buffer page to copy into
source_dataOptional source data which is read from

Definition at line 738 of file smagorinsky.F90.

738  type(model_state_type), intent(inout) :: current_state
739  integer, intent(in) :: corner_loc, pid_location, x_source_index, y_source_index
740  integer, intent(inout) :: current_page(:)
741  type(neighbour_description_type), intent(inout) :: neighbour_description
742  type(field_data_wrapper_type), dimension(:), intent(in), optional :: source_data
743 
744  call copy_corner_to_buffer(current_state%local_grid, neighbour_description%send_corner_buffer, &
745  current_state%diff_coefficient%data, corner_loc, x_source_index, y_source_index, current_page(pid_location))
746 
747  current_page(pid_location)=current_page(pid_location)+1
Here is the caller graph for this function:

◆ copy_diff_to_halo_buffer()

subroutine smagorinsky_mod::copy_diff_to_halo_buffer ( type(model_state_type), intent(inout)  current_state,
type(neighbour_description_type), intent(inout)  neighbour_description,
integer, intent(in)  dim,
integer, intent(in)  source_index,
integer, intent(in)  pid_location,
integer, dimension(:), intent(inout)  current_page,
type(field_data_wrapper_type), dimension(:), intent(in), optional  source_data 
)
private

Copies the diffusion field data to halo buffers for a specific process in a dimension and halo cell.

Parameters
current_stateThe current model state
neighbour_descriptionsDescription of the neighbour halo swapping status
dimDimension to copy from
source_indexThe source index of the dimension we are reading from in the prognostic field
pid_locationLocation of the neighbouring process in the local stored data structures
current_pageThe current (next) buffer page to copy into
source_dataOptional source data which is read from

Definition at line 715 of file smagorinsky.F90.

715  type(model_state_type), intent(inout) :: current_state
716  integer, intent(in) :: dim, pid_location, source_index
717  integer, intent(inout) :: current_page(:)
718  type(neighbour_description_type), intent(inout) :: neighbour_description
719  type(field_data_wrapper_type), dimension(:), intent(in), optional :: source_data
720 
721  call copy_field_to_buffer(current_state%local_grid, neighbour_description%send_halo_buffer, &
722  current_state%diff_coefficient%data, dim, source_index, current_page(pid_location))
723 
724  current_page(pid_location)=current_page(pid_location)+1
Here is the caller graph for this function:

◆ copy_vis_corners_to_halo_buffer()

subroutine smagorinsky_mod::copy_vis_corners_to_halo_buffer ( type(model_state_type), intent(inout)  current_state,
type(neighbour_description_type), intent(inout)  neighbour_description,
integer, intent(in)  corner_loc,
integer, intent(in)  x_source_index,
integer, intent(in)  y_source_index,
integer, intent(in)  pid_location,
integer, dimension(:), intent(inout)  current_page,
type(field_data_wrapper_type), dimension(:), intent(in), optional  source_data 
)
private

Copies the viscosity field corner data to halo buffers for a specific process.

Parameters
current_stateThe current model state
neighbour_descriptionsDescription of the neighbour halo swapping status
corner_locThe corner location
x_source_indexThe X source index of the dimension we are reading from in the prognostic field
y_source_indexThe Y source index of the dimension we are reading from in the prognostic field
pid_locationLocation of the neighbouring process in the local stored data structures
current_pageThe current (next) buffer page to copy into
source_dataOptional source data which is read from

Definition at line 693 of file smagorinsky.F90.

693  type(model_state_type), intent(inout) :: current_state
694  integer, intent(in) :: corner_loc, pid_location, x_source_index, y_source_index
695  integer, intent(inout) :: current_page(:)
696  type(neighbour_description_type), intent(inout) :: neighbour_description
697  type(field_data_wrapper_type), dimension(:), intent(in), optional :: source_data
698 
699  call copy_corner_to_buffer(current_state%local_grid, neighbour_description%send_corner_buffer, &
700  current_state%vis_coefficient%data, corner_loc, x_source_index, y_source_index, current_page(pid_location))
701 
702  current_page(pid_location)=current_page(pid_location)+1
Here is the caller graph for this function:

◆ copy_vis_to_halo_buffer()

subroutine smagorinsky_mod::copy_vis_to_halo_buffer ( type(model_state_type), intent(inout)  current_state,
type(neighbour_description_type), intent(inout)  neighbour_description,
integer, intent(in)  dim,
integer, intent(in)  source_index,
integer, intent(in)  pid_location,
integer, dimension(:), intent(inout)  current_page,
type(field_data_wrapper_type), dimension(:), intent(in), optional  source_data 
)
private

Copies the viscosity field data to halo buffers for a specific process in a dimension and halo cell.

Parameters
current_stateThe current model state
neighbour_descriptionsDescription of the neighbour halo swapping status
dimDimension to copy from
source_indexThe source index of the dimension we are reading from in the prognostic field
pid_locationLocation of the neighbouring process in the local stored data structures
current_pageThe current (next) buffer page to copy into
source_dataOptional source data which is read from

Definition at line 670 of file smagorinsky.F90.

670  type(model_state_type), intent(inout) :: current_state
671  integer, intent(in) :: dim, pid_location, source_index
672  integer, intent(inout) :: current_page(:)
673  type(neighbour_description_type), intent(inout) :: neighbour_description
674  type(field_data_wrapper_type), dimension(:), intent(in), optional :: source_data
675 
676  call copy_field_to_buffer(current_state%local_grid, neighbour_description%send_halo_buffer, &
677  current_state%vis_coefficient%data, dim, source_index, current_page(pid_location))
678 
679  current_page(pid_location)=current_page(pid_location)+1
Here is the caller graph for this function:

◆ finalisation_callback()

subroutine smagorinsky_mod::finalisation_callback ( type(model_state_type), intent(inout), target  current_state)
private

Called when the model is finishing up, will finalise the halo communications represented by the state.

Parameters
current_stateThe current model state

Definition at line 133 of file smagorinsky.F90.

133  type(model_state_type), target, intent(inout) :: current_state
134 
135  call finalise_halo_communication(current_state%viscosity_halo_swap_state)
136  call finalise_halo_communication(current_state%diffusion_halo_swap_state)
Here is the caller graph for this function:

◆ initialisation_callback()

subroutine smagorinsky_mod::initialisation_callback ( type(model_state_type), intent(inout), target  current_state)
private

Initialisation call back which will read in the coriolis configuration and set up the geostrophic winds.

Parameters
current_stateThe current model state

Definition at line 45 of file smagorinsky.F90.

45  type(model_state_type), target, intent(inout) :: current_state
46 
47 
48  if (.not. is_component_enabled(current_state%options_database, "diffusion")) then
49  call log_master_log(log_error, "Smagorinsky requires the diffusion component to be enabled")
50  end if
51  if (.not. is_component_enabled(current_state%options_database, "viscosity")) then
52  call log_master_log(log_error, "Smagorinsky requires the viscosity component to be enabled")
53  end if
54  if (.not. current_state%use_viscosity_and_diffusion) then
55  call log_master_log(log_error, "Smagorinsky requires use_viscosity_and_diffusion=.true. or monc will fail")
56  endif
57 
58  if (.not. is_component_enabled(current_state%options_database, "lower_bc")) then
59  call log_master_log(log_info, "LOWERBC is disabled, zero diff and vis_coeff on level 1")
60  current_state%vis_coefficient%data(1,:,:)=0.0_default_precision
61  current_state%diff_coefficient%data(1,:,:)=0.0_default_precision
62  endif
63 
64  eps=0.01_default_precision
65  repsh=0.5_default_precision/eps
66  thcona=ratio_mol_wts*current_state%thref0
67  thconb=thcona-current_state%thref0
68  thconap1=thcona
69 
70  subb=options_get_real(current_state%options_database, "smag-subb")
71  subc=options_get_real(current_state%options_database, "smag-subc")
72 
73  subg=1.2_default_precision
74  subh=0.0_default_precision
75  subr=4.0_default_precision
76  pr_n=0.7_default_precision
77  suba=1.0_default_precision/pr_n
78  ric=0.25_default_precision
79  ricinv=1.0_default_precision/ric
80  if (current_state%use_viscosity_and_diffusion) then
81  call init_halo_communication(current_state, get_single_field_per_halo_cell, &
82  current_state%viscosity_halo_swap_state, 2, .true.)
83  call init_halo_communication(current_state, get_single_field_per_halo_cell, &
84  current_state%diffusion_halo_swap_state, 2, .true.)
85  end if
86 
87  iqv = get_q_index(standard_q_names%VAPOUR, 'smagorinsky')
88  current_state%water_vapour_mixing_ratio_index=iqv
89 
90  iql = get_q_index(standard_q_names%CLOUD_LIQUID_MASS, 'smagorinsky')
91  current_state%liquid_water_mixing_ratio_index=iql
92 
Here is the caller graph for this function:

◆ moist_ri_1()

real(kind=default_precision) function, dimension(current_state%local_grid%size(z_index)) smagorinsky_mod::moist_ri_1 ( type(model_state_type), intent(inout), target  current_state,
real(kind=default_precision), dimension(:), intent(in)  ssq,
type(prognostic_field_type), intent(inout)  th,
type(prognostic_field_type), dimension(:), intent(inout)  q 
)
private

Calculates numerator of Richarson number as the difference in buoyancy between a parcel at K+1 and a parcel lifted from K to K+1. Condensation is calculated by assuming that TL=T-(Lv/CP)*QL+gz/CP and QT = QV+QL are conserved on lifting. Microphysics conversions are assumed to happen slowly compared to turbulent processes so are, therefore, neglected in the calculation.

Parameters
current_stateThe current model state
ssqThe half squared strain rate
Returns
The Richardson number calculated here for a specific column

Definition at line 333 of file smagorinsky.F90.

333  type(model_state_type), target, intent(inout) :: current_state
334  real(kind=DEFAULT_PRECISION), dimension(:), intent(in) :: ssq
335  type(prognostic_field_type), intent(inout) :: th
336  type(prognostic_field_type), dimension(:), intent(inout) :: q
337  real(kind=DEFAULT_PRECISION), dimension(current_state%local_grid%size(Z_INDEX)) :: moist_ri_1
338 
339  integer :: k
340  real(kind=DEFAULT_PRECISION) :: liquid_water_at_kp1, temperature_at_kp1, dtlref, tmpinv, &
341  liquid_water_static_energy_temperature, total_water_substance, total_water_substance_p1, scltmp, qloading
342 
343  ! Set up the liquid water static energy & total water
344  do k=2, current_state%local_grid%size(z_index)
345  liquid_water_static_energy_temperature = &
346  th%data(k, current_state%column_local_y, current_state%column_local_x)*&
347  current_state%global_grid%configuration%vertical%rprefrcp(k) - rlvap_over_cp*&
348  q(current_state%liquid_water_mixing_ratio_index)%data(k, current_state%column_local_y, current_state%column_local_x)
349  if (k==2) then
350  total_water_substance = q(current_state%water_vapour_mixing_ratio_index)%data(&
351  k, current_state%column_local_y, current_state%column_local_x) + &
352  q(current_state%liquid_water_mixing_ratio_index)%data(k, current_state%column_local_y, current_state%column_local_x)
353  else
354  total_water_substance=total_water_substance_p1
355  end if
356  total_water_substance_p1 = q(current_state%water_vapour_mixing_ratio_index)%data(&
357  k+1, current_state%column_local_y, current_state%column_local_x) + &
358  q(current_state%liquid_water_mixing_ratio_index)%data(k+1, current_state%column_local_y, current_state%column_local_x)
359  if (k .le. current_state%local_grid%size(z_index)-1) then
360  tmpinv = 1.0_default_precision/ssq(k)
361  dtlref = current_state%global_grid%configuration%vertical%tref(k+1)-&
362  current_state%global_grid%configuration%vertical%tref(k)+(g/cp)*&
363  current_state%global_grid%configuration%vertical%dzn(k+1)
364  temperature_at_kp1=current_state%global_grid%configuration%vertical%tstarpr(k+1)+&
365  current_state%global_grid%configuration%vertical%tref(k+1)+&
366  (liquid_water_static_energy_temperature+rlvap_over_cp*(total_water_substance-&
367  current_state%global_grid%configuration%vertical%qsat(k+1)) - &
368  dtlref - current_state%global_grid%configuration%vertical%tstarpr(k+1) )*&
369  current_state%global_grid%configuration%vertical%qsatfac(k+1)
370  liquid_water_at_kp1=total_water_substance -( current_state%global_grid%configuration%vertical%qsat(k+1)+&
371  current_state%global_grid%configuration%vertical%dqsatdt(k+1)*(temperature_at_kp1- &
372  current_state%global_grid%configuration%vertical%tstarpr(k+1)-&
373  current_state%global_grid%configuration%vertical%tref(k+1)) )
374  if (liquid_water_at_kp1 .le. 0.0_default_precision) then
375  liquid_water_at_kp1=0.0_default_precision
376  temperature_at_kp1=current_state%global_grid%configuration%vertical%tref(k+1)+&
377  (liquid_water_static_energy_temperature-dtlref)
378  end if
379  scltmp=current_state%global_grid%configuration%vertical%rdzn(k+1)*&
380  current_state%global_grid%configuration%vertical%buoy_co(k)
381  ! Calculate the Richardson number taking into account the water
382  ! loading terms from the other hydrometeors in the buoyancy
383  ! calculation
384  qloading = current_state%cq(current_state%water_vapour_mixing_ratio_index)*&
385  (total_water_substance_p1-total_water_substance)-ratio_mol_wts*&
386  (q(current_state%liquid_water_mixing_ratio_index)%data(&
387  k+1, current_state%column_local_y, current_state%column_local_x)-liquid_water_at_kp1)
388  !if(IRAINP.gt.0) qloading = qloading + current_state%cq(IQR)*(l_q(J,K+1,IQR)-l_q(J,K,IQR))
389  !if(ISNOWP.gt.0) qloading = qloading + current_state%cq(IQS)*(l_q(J,K+1,IQS)-l_q(J,K,IQS))
390  !if(ICECLP.gt.0) qloading = qloading + current_state%cq(IQI)*(l_q(J,K+1,IQI)-l_q(J,K,IQI))
391  !if(IGRAUP.gt.0) qloading = qloading + currnet_state%cq(IQG)*(l_q(J,K+1,IQG)-l_q(J,K,IQG))
392 
393  moist_ri_1(k) = ((th%data(k+1, current_state%column_local_y, current_state%column_local_x)+&
394  current_state%global_grid%configuration%vertical%thref(k+1) - temperature_at_kp1*&
395  current_state%global_grid%configuration%vertical%prefrcp(k+1))&
396  +current_state%global_grid%configuration%vertical%thref(k+1)*qloading)*scltmp*tmpinv
397  end if
398  end do
Here is the caller graph for this function:

◆ moist_ri_2()

real(kind=default_precision) function, dimension(current_state%local_grid%size(z_index)) smagorinsky_mod::moist_ri_2 ( type(model_state_type), intent(inout), target  current_state,
real(kind=default_precision), dimension(:), intent(in)  ssq,
type(prognostic_field_type), intent(inout)  th,
type(prognostic_field_type), dimension(:), intent(inout)  q 
)
private

Calculates another "moist version" of the Richardson number based on "change in subgrid buoyancy flux when a fraction EPS is exchanged".

Parameters
current_stateThe current model state
ssqThe half squared strain rate
Returns
The Richardson number calculated here for a specific column

Definition at line 238 of file smagorinsky.F90.

238  type(model_state_type), target, intent(inout) :: current_state
239  real(kind=DEFAULT_PRECISION), dimension(:), intent(in) :: ssq
240  type(prognostic_field_type), intent(inout) :: th
241  type(prognostic_field_type), dimension(:), intent(inout) :: q
242  real(kind=DEFAULT_PRECISION), dimension(current_state%local_grid%size(Z_INDEX)) :: moist_ri_2
243 
244  integer :: k
245  real(kind=DEFAULT_PRECISION) :: tmpinv, scltmp, qt_k, qt_kp1, thlpr_k, thlpr_kp1, qlli, qlui, &
246  thvli, thvui, thlpr_un, thlpr_ln, qtun, qtln, thvln, thvun, qlln, qlun
247 
248  do k=2,current_state%local_grid%size(z_index)-1
249  tmpinv = 1.0_default_precision/ssq(k)
250  scltmp=current_state%global_grid%configuration%vertical%rdzn(k+1)*&
251  current_state%global_grid%configuration%vertical%buoy_co(k)
252  if (current_state%use_anelastic_equations) then
253  thcona=ratio_mol_wts*current_state%global_grid%configuration%vertical%thref(k)
254  thconb=0.5_default_precision*(ratio_mol_wts-1.0_default_precision)*&
255  (current_state%global_grid%configuration%vertical%thref(k)+&
256  current_state%global_grid%configuration%vertical%thref(k+1))
257  thconap1=ratio_mol_wts*current_state%global_grid%configuration%vertical%thref(k+1)
258  end if
259 
260  qt_k = q(current_state%water_vapour_mixing_ratio_index)%data(&
261  k, current_state%column_local_y, current_state%column_local_x)+&
262  q(current_state%liquid_water_mixing_ratio_index)%data(k, current_state%column_local_y, current_state%column_local_x)
263  qt_kp1 = q(current_state%water_vapour_mixing_ratio_index)%data(&
264  k+1, current_state%column_local_y, current_state%column_local_x)+&
265  q(current_state%liquid_water_mixing_ratio_index)%data(k+1, current_state%column_local_y, current_state%column_local_x)
266  thlpr_k = th%data(k, current_state%column_local_y, current_state%column_local_x) -&
267  rlvap_over_cp*q(current_state%liquid_water_mixing_ratio_index)%data(&
268  k, current_state%column_local_y, current_state%column_local_x)*&
269  current_state%global_grid%configuration%vertical%prefrcp(k)
270  ! calculate theta_l at constant (level K) pressure
271  thlpr_kp1 = th%data(k+1, current_state%column_local_y, current_state%column_local_x) - &
272  rlvap_over_cp*q(current_state%liquid_water_mixing_ratio_index)%data(&
273  k+1, current_state%column_local_y, current_state%column_local_x)*&
274  current_state%global_grid%configuration%vertical%prefrcp(k)
275 
276  !.......CALCULATE QL AND BUOYANCY OF LEVELS K AND K+1 (LOWER AND UPPER)
277  !.......QL CALC AT CONSTANT PRESSURE
278  qlli = max(0.0_default_precision, (qt_k-(current_state%global_grid%configuration%vertical%qsat(k) +&
279  current_state%global_grid%configuration%vertical%dqsatdt(k)*&
280  (current_state%global_grid%configuration%vertical%rprefrcp(k)*thlpr_k-&
281  current_state%global_grid%configuration%vertical%tstarpr(k))))*&
282  current_state%global_grid%configuration%vertical%qsatfac(k))
283  qlui = max(0.0_default_precision, (qt_kp1-(current_state%global_grid%configuration%vertical%qsat(k) +&
284  current_state%global_grid%configuration%vertical%dqsatdt(k)*&
285  (current_state%global_grid%configuration%vertical%rprefrcp(k)*(thlpr_kp1+&
286  current_state%global_grid%configuration%vertical%thref(k+1))-&
287  current_state%global_grid%configuration%vertical%tstarpr(k)-&
288  current_state%global_grid%configuration%vertical%tref(k))))*&
289  current_state%global_grid%configuration%vertical%qsatfac(k))
290  thvli = thlpr_k+current_state%global_grid%configuration%vertical%thref(k) +&
291  ( rlvap_over_cp*current_state%global_grid%configuration%vertical%prefrcp(k)-thcona ) * qlli
292  thvui = thlpr_kp1+current_state%global_grid%configuration%vertical%thref(k+1) +&
293  ( rlvap_over_cp*current_state%global_grid%configuration%vertical%prefrcp(k)-thconap1 ) * qlui
294 
295  !.......MIX CONSERVATIVE QUANTITIES AT CONSTANT PRESSURE
296  !.......NOTE THAT FOR COMPUTATIONAL CONVENIENCE THE RELEVANT LEVEL'S
297  !.......THREF HAS BEEN OMITTED FROM ALL THV'S, THLPR_LN AND THLPR_UN AS
298  !.......ONLY THE DIFFERENCE BETWEEN THV'S AT THE SAME LEVEL IS REQUIRED.
299  thlpr_un = thlpr_kp1 - eps*(thlpr_kp1 - thlpr_k + current_state%global_grid%configuration%vertical%dthref(k))
300  thlpr_ln = thlpr_k + eps*(thlpr_kp1 - thlpr_k + current_state%global_grid%configuration%vertical%dthref(k))
301  qtun = qt_kp1 - eps*(qt_kp1 - qt_k)
302  qtln = qt_k + eps*(qt_kp1 - qt_k)
303  !.......CALCULATE LIQUID WATER AND BUOYANCY OF MIXED AIR
304  qlln = max(0.0_default_precision,( qtln - ( current_state%global_grid%configuration%vertical%qsat(k) +&
305  current_state%global_grid%configuration%vertical%dqsatdt(k)*&
306  (current_state%global_grid%configuration%vertical%rprefrcp(k)*&
307  thlpr_ln-current_state%global_grid%configuration%vertical%tstarpr(k)))&
308  )*current_state%global_grid%configuration%vertical%qsatfac(k))
309  qlun = max(0.0_default_precision,( qtun - ( current_state%global_grid%configuration%vertical%qsat(k) + &
310  current_state%global_grid%configuration%vertical%dqsatdt(k)*&
311  (current_state%global_grid%configuration%vertical%rprefrcp(k)*&
312  (thlpr_un+current_state%global_grid%configuration%vertical%thref(k+1))-&
313  current_state%global_grid%configuration%vertical%tstarpr(k)-current_state%global_grid%configuration%vertical%tref(k)))&
314  )*current_state%global_grid%configuration%vertical%qsatfac(k))
315  thvln = thlpr_ln+current_state%global_grid%configuration%vertical%thref(k) +&
316  ( rlvap_over_cp*current_state%global_grid%configuration%vertical%prefrcp(k)-thcona) * qlln
317  thvun = thlpr_un+current_state%global_grid%configuration%vertical%thref(k+1) +&
318  ( rlvap_over_cp*current_state%global_grid%configuration%vertical%prefrcp(k)-thconap1) * qlun
319 
320  !.....CHANGE IN BUOYANCY DUE TO FRAC MIXING IS USED TO DETERMINE l_ri
321  moist_ri_2(k) = scltmp*tmpinv*(thconb*(qt_kp1-qt_k) +((thvln-thvli)-(thvun-thvui))*repsh)
322  end do
Here is the caller graph for this function:

◆ setfri()

subroutine smagorinsky_mod::setfri ( type(model_state_type), intent(inout), target  current_state,
real(kind=default_precision), dimension(:), intent(in)  richardson_number,
real(kind=default_precision), dimension(:), intent(in)  ssq 
)
private

Calculates the eddy viscosity (VIS) and diffusivity (DIFF) depending on the Richardson Number (RI) and half squared strain rate.

Parameters
current_stateThe current model state
richardson_numberRichardson number
ssqThe half squared strain rate

Definition at line 163 of file smagorinsky.F90.

163  type(model_state_type), target, intent(inout) :: current_state
164  real(kind=DEFAULT_PRECISION), dimension(:), intent(in) :: richardson_number, ssq
165 
166  integer :: k
167  real(kind=DEFAULT_PRECISION) :: rifac, sctmp
168 
169  do k=2,current_state%local_grid%size(z_index)-1
170  if (richardson_number(k) .le. 0.0_default_precision) then
171  current_state%vis_coefficient%data(k, current_state%column_local_y, current_state%column_local_x)=&
172  sqrt(1.0_default_precision-subc*richardson_number(k))
173  current_state%diff_coefficient%data(k, current_state%column_local_y, current_state%column_local_x)=&
174  suba*sqrt(1.-subb*richardson_number(k))
175  else if ((richardson_number(k) .gt. 0.0_default_precision) .and. (richardson_number(k) .lt. ric)) then
176  rifac=(1.-richardson_number(k)*ricinv)**4
177  current_state%vis_coefficient%data(k, current_state%column_local_y, current_state%column_local_x)=&
178  rifac*(1.-subh*richardson_number(k))
179  current_state%diff_coefficient%data(k, current_state%column_local_y, current_state%column_local_x)=&
180  rifac*suba*(1.-subg*richardson_number(k))
181  else
182  current_state%vis_coefficient%data(k, current_state%column_local_y, current_state%column_local_x)=0.0_default_precision
183  current_state%diff_coefficient%data(k, current_state%column_local_y, current_state%column_local_x)=0.0_default_precision
184  end if
185  sctmp=current_state%global_grid%configuration%vertical%rneutml_sq(k)*sqrt(ssq(k))
186  current_state%vis_coefficient%data(k, current_state%column_local_y, current_state%column_local_x)=&
187  current_state%vis_coefficient%data(k, current_state%column_local_y, current_state%column_local_x)*sctmp
188  current_state%diff_coefficient%data(k, current_state%column_local_y, current_state%column_local_x)=&
189  current_state%diff_coefficient%data(k, current_state%column_local_y, current_state%column_local_x)*sctmp
190  end do
191 
192  current_state%vis_coefficient%data(current_state%local_grid%size(z_index), &
193  current_state%column_local_y, current_state%column_local_x)= 0.0_default_precision
194  current_state%diff_coefficient%data(current_state%local_grid%size(z_index), &
195  current_state%column_local_y, current_state%column_local_x)= 0.0_default_precision
Here is the caller graph for this function:

◆ smagorinsky_get_descriptor()

type(component_descriptor_type) function, public smagorinsky_mod::smagorinsky_get_descriptor ( )

Provides the descriptor back to the caller and is used in component registration.

Returns
The termination check component descriptor

Definition at line 35 of file smagorinsky.F90.

35  smagorinsky_get_descriptor%name="smagorinsky"
36  smagorinsky_get_descriptor%version=0.1
37  smagorinsky_get_descriptor%initialisation=>initialisation_callback
38  smagorinsky_get_descriptor%timestep=>timestep_callback
39  smagorinsky_get_descriptor%finalisation=>finalisation_callback
Here is the call graph for this function:

◆ timestep_callback()

subroutine smagorinsky_mod::timestep_callback ( type(model_state_type), intent(inout), target  current_state)
private

For each none halo cell this will calculate the subgrid terms for viscosity and diffusion.

Parameters
current_stateThe current model state

Definition at line 98 of file smagorinsky.F90.

98  type(model_state_type), target, intent(inout) :: current_state
99 
100  real(kind=DEFAULT_PRECISION), dimension(current_state%local_grid%size(Z_INDEX)) :: richardson_number, ssq
101 
102  if (.not. current_state%halo_column) then
103  if (.not. current_state%use_viscosity_and_diffusion) then
104  current_state%vis_coefficient%data(2:,:,:)=0.0_default_precision
105  current_state%diff_coefficient%data(2:,:,:)=0.0_default_precision
106  current_state%cvis=0.0_default_precision
107  else
108  if (current_state%field_stepping == forward_stepping) then
109  ssq=calculate_half_squared_strain_rate(current_state, current_state%u, current_state%v, current_state%w)
110  richardson_number=calculate_richardson_number(current_state, ssq, current_state%th, current_state%q)
111  else
112  ssq=calculate_half_squared_strain_rate(current_state, current_state%zu, current_state%zv, current_state%zw)
113  richardson_number=calculate_richardson_number(current_state, ssq, current_state%zth, current_state%zq)
114  end if
115  call setfri(current_state, richardson_number, ssq)
116  if (is_component_enabled(current_state%options_database, "cfltest")) then
117  call update_viscous_number(current_state)
118  endif
119  end if
120  end if
121  if (current_state%last_timestep_column) then
122  ! need to ensure not already in progress
123  call initiate_nonblocking_halo_swap(current_state, current_state%diffusion_halo_swap_state, &
124  copy_diff_to_halo_buffer, copy_diff_corners_to_halo_buffer)
125  call initiate_nonblocking_halo_swap(current_state, current_state%viscosity_halo_swap_state, &
126  copy_vis_to_halo_buffer, copy_vis_corners_to_halo_buffer)
127  end if
Here is the call graph for this function:
Here is the caller graph for this function:

◆ update_viscous_number()

subroutine smagorinsky_mod::update_viscous_number ( type(model_state_type), intent(inout), target  current_state)
private

Update viscous number based upon the viscosity and diffusivity coefficients.

Parameters
current_stateThe current model state

Definition at line 142 of file smagorinsky.F90.

142  type(model_state_type), target, intent(inout) :: current_state
143 
144  integer :: k
145 
146  if (mod(current_state%timestep, current_state%cfl_frequency) == 1 .or. &
147  current_state%timestep-current_state%start_timestep .le. current_state%cfl_frequency) then
148  do k=2, current_state%local_grid%size(z_index)-1
149  current_state%cvis=max(current_state%cvis, max(current_state%vis_coefficient%data(k, current_state%column_local_y, &
150  current_state%column_local_x),current_state%diff_coefficient%data(k, current_state%column_local_y, &
151  current_state%column_local_x))*(current_state%global_grid%configuration%vertical%rdzn(k+1)**2+&
152  current_state%global_grid%configuration%horizontal%cx2+current_state%global_grid%configuration%horizontal%cy2))
153  end do
154  end if
Here is the caller graph for this function:

Variable Documentation

◆ eps

real(kind=default_precision) smagorinsky_mod::eps
private

Definition at line 23 of file smagorinsky.F90.

23  real(kind=DEFAULT_PRECISION) :: eps, repsh, thcona, thconb, thconap1, suba, subb, subc, subg, subh, subr, pr_n, ric, ricinv

◆ iql

integer smagorinsky_mod::iql
private

Definition at line 28 of file smagorinsky.F90.

◆ iqv

integer smagorinsky_mod::iqv
private

Definition at line 28 of file smagorinsky.F90.

28  integer :: iqv, iql ! index for water vapour and liquid

◆ pr_n

real(kind=default_precision) smagorinsky_mod::pr_n
private

Definition at line 23 of file smagorinsky.F90.

◆ repsh

real(kind=default_precision) smagorinsky_mod::repsh
private

Definition at line 23 of file smagorinsky.F90.

◆ ric

real(kind=default_precision) smagorinsky_mod::ric
private

Definition at line 23 of file smagorinsky.F90.

◆ richardson_number_calculation

integer, parameter smagorinsky_mod::richardson_number_calculation =2
private

Definition at line 22 of file smagorinsky.F90.

22  integer, parameter :: richardson_number_calculation=2

◆ ricinv

real(kind=default_precision) smagorinsky_mod::ricinv
private

Definition at line 23 of file smagorinsky.F90.

◆ suba

real(kind=default_precision) smagorinsky_mod::suba
private

Definition at line 23 of file smagorinsky.F90.

◆ subb

real(kind=default_precision) smagorinsky_mod::subb
private

Definition at line 23 of file smagorinsky.F90.

◆ subc

real(kind=default_precision) smagorinsky_mod::subc
private

Definition at line 23 of file smagorinsky.F90.

◆ subg

real(kind=default_precision) smagorinsky_mod::subg
private

Definition at line 23 of file smagorinsky.F90.

◆ subh

real(kind=default_precision) smagorinsky_mod::subh
private

Definition at line 23 of file smagorinsky.F90.

◆ subr

real(kind=default_precision) smagorinsky_mod::subr
private

Definition at line 23 of file smagorinsky.F90.

◆ thcona

real(kind=default_precision) smagorinsky_mod::thcona
private

Definition at line 23 of file smagorinsky.F90.

◆ thconap1

real(kind=default_precision) smagorinsky_mod::thconap1
private

Definition at line 23 of file smagorinsky.F90.

◆ thconb

real(kind=default_precision) smagorinsky_mod::thconb
private

Definition at line 23 of file smagorinsky.F90.