MONC
Functions/Subroutines | Variables
forcing_mod Module Reference

Forcing, both subsidence and large scale. More...

Functions/Subroutines

type(component_descriptor_type) function, public forcing_get_descriptor ()
 Provides the component descriptor for the core to register. More...
 
subroutine field_information_retrieval_callback (current_state, name, field_information)
 Field information retrieval callback, this returns information for a specific components published field. More...
 
subroutine field_value_retrieval_callback (current_state, name, field_value)
 Field value retrieval callback, this returns the value of a specific published field. More...
 
subroutine init_callback (current_state)
 Initialises the forcing data structures. More...
 
subroutine timestep_callback (current_state)
 Called for each data column and will determine the forcing values in x and y which are then applied to the field source terms. More...
 
subroutine apply_subsidence_to_flow_fields (current_state)
 
subroutine apply_subsidence_to_theta (current_state)
 
subroutine apply_subsidence_to_q_fields (current_state)
 
subroutine apply_time_independent_forcing_to_theta (current_state)
 
subroutine apply_time_independent_forcing_to_q (current_state)
 
subroutine apply_time_independent_forcing_to_u (current_state)
 
subroutine apply_time_independent_forcing_to_v (current_state)
 
subroutine finalisation_callback (current_state)
 Finalises the component @current_state Current model state. More...
 
real(kind=default_precision) function, dimension(size(diagnostics_summed)) get_averaged_diagnostics (current_state, diagnostics_summed)
 Averages some diagnostic values across all local horizontal points. More...
 
subroutine save_precomponent_tendencies (current_state, cxn, cyn, txn, tyn)
 Save the 3d tendencies coming into the component. More...
 
subroutine compute_component_tendencies (current_state, cxn, cyn, txn, tyn)
 Computation of component tendencies. More...
 
subroutine set_published_field_value (field_value, real_1d_field, real_2d_field, real_3d_field)
 Sets the published field value from the temporary diagnostic values held by this component. More...
 
subroutine check_forcing_status (status)
 Will check a NetCDF status and write to log_log error any decoded statuses. More...
 
subroutine read_2d_forcing_dimensions (ncid, time_dim, z_dim)
 Reads the dimensions for forcing from the NetCDF file. This routine assumes the forcing uses only time and height. More...
 
subroutine read_2d_forcing_variables (filename, ncid, time_dim, time, z_dim, z_profile, force_2d_key, force_2d_var)
 Reads the variables from the NetCDF forcing file. The 2d variables are assumed to be time and height. More...
 
subroutine read_single_forcing_variable (ncid, key, data1d, data2d)
 Reads a single variable out of a NetCDF file. More...
 

Variables

character(len= *), parameter time_key = "time"
 NetCDF data time key. More...
 
character(len= *), parameter z_key = "z"
 NetCDF data height(z) key. More...
 
character(len= *), parameter wsubs_key = "wsubs"
 NetCDF data subsidence velocity. More...
 
integer, parameter max_file_len =200
 Maximum length of surface condition input filename. More...
 
character(max_file_leninput_file
 
integer, parameter divergence =0
 
integer, parameter subsidence =1
 
integer, parameter tendency =0
 
integer, parameter relaxation =1
 
integer, parameter increment =2
 
real(kind=default_precision), dimension(:), allocatable theta_profile
 
real(kind=default_precision), dimension(:), allocatable q_profile
 
real(kind=default_precision), dimension(:), allocatable u_profile
 
real(kind=default_precision), dimension(:), allocatable v_profile
 
real(kind=default_precision), dimension(:), allocatable dtheta_profile
 
real(kind=default_precision), dimension(:), allocatable dq_profile
 
real(kind=default_precision), dimension(:), allocatable du_profile
 
real(kind=default_precision), dimension(:), allocatable dv_profile
 
real(kind=default_precision), dimension(:), allocatable du_profile_diag
 
real(kind=default_precision), dimension(:), allocatable dv_profile_diag
 
real(kind=default_precision), dimension(:), allocatable dtheta_profile_diag
 
real(kind=default_precision), dimension(:,:), allocatable dq_profile_diag
 
real(kind=default_precision), dimension(:), allocatable du_subs_profile_diag
 
real(kind=default_precision), dimension(:), allocatable dv_subs_profile_diag
 
real(kind=default_precision), dimension(:), allocatable dtheta_subs_profile_diag
 
real(kind=default_precision), dimension(:,:), allocatable dq_subs_profile_diag
 
real(kind=default_precision), dimension(:,:), allocatable w_subs_varies_with_time
 
real(kind=default_precision), dimension(:), allocatable forcing_input_times
 
real(kind=default_precision) forcing_timescale_theta
 
real(kind=default_precision) forcing_timescale_q
 
real(kind=default_precision) forcing_timescale_u
 
real(kind=default_precision) forcing_timescale_v
 
logical l_constant_forcing_theta
 
logical l_constant_forcing_q
 
logical l_constant_forcing_u
 
logical l_constant_forcing_v
 
integer constant_forcing_type_theta =TENDENCY
 
integer constant_forcing_type_q =TENDENCY
 
integer constant_forcing_type_u =RELAXATION
 
integer constant_forcing_type_v =RELAXATION
 
logical l_constant_forcing_theta_z2pressure
 
logical relax_to_initial_u_profile
 
logical relax_to_initial_v_profile
 
logical relax_to_initial_theta_profile
 
logical use_time_varying_subsidence
 
logical l_subs_pl_theta
 
logical l_subs_pl_q
 
logical l_subs_local_theta
 
logical l_subs_local_q
 
character(len=string_length), dimension(:), allocatable names_force_pl_q
 
real(kind=default_precision), dimension(:,:,:), allocatable tend_3d_u
 
real(kind=default_precision), dimension(:,:,:), allocatable tend_3d_v
 
real(kind=default_precision), dimension(:,:,:), allocatable tend_3d_th
 
real(kind=default_precision), dimension(:,:,:), allocatable tend_3d_qv
 
real(kind=default_precision), dimension(:,:,:), allocatable tend_3d_ql
 
real(kind=default_precision), dimension(:,:,:), allocatable tend_3d_qi
 
real(kind=default_precision), dimension(:,:,:), allocatable tend_3d_qr
 
real(kind=default_precision), dimension(:,:,:), allocatable tend_3d_qs
 
real(kind=default_precision), dimension(:,:,:), allocatable tend_3d_qg
 
real(kind=default_precision), dimension(:,:,:), allocatable tend_3d_tabs
 
logical l_tend_3d_u
 
logical l_tend_3d_v
 
logical l_tend_3d_th
 
logical l_tend_3d_qv
 
logical l_tend_3d_ql
 
logical l_tend_3d_qi
 
logical l_tend_3d_qr
 
logical l_tend_3d_qs
 
logical l_tend_3d_qg
 
logical l_tend_3d_tabs
 
real(kind=default_precision), dimension(:), allocatable tend_pr_tot_u
 
real(kind=default_precision), dimension(:), allocatable tend_pr_tot_v
 
real(kind=default_precision), dimension(:), allocatable tend_pr_tot_th
 
real(kind=default_precision), dimension(:), allocatable tend_pr_tot_qv
 
real(kind=default_precision), dimension(:), allocatable tend_pr_tot_ql
 
real(kind=default_precision), dimension(:), allocatable tend_pr_tot_qi
 
real(kind=default_precision), dimension(:), allocatable tend_pr_tot_qr
 
real(kind=default_precision), dimension(:), allocatable tend_pr_tot_qs
 
real(kind=default_precision), dimension(:), allocatable tend_pr_tot_qg
 
real(kind=default_precision), dimension(:), allocatable tend_pr_tot_tabs
 
logical l_tend_pr_tot_u
 
logical l_tend_pr_tot_v
 
logical l_tend_pr_tot_th
 
logical l_tend_pr_tot_qv
 
logical l_tend_pr_tot_ql
 
logical l_tend_pr_tot_qi
 
logical l_tend_pr_tot_qr
 
logical l_tend_pr_tot_qs
 
logical l_tend_pr_tot_qg
 
logical l_tend_pr_tot_tabs
 
integer iqv =0
 
integer iql =0
 
integer iqr =0
 
integer iqi =0
 
integer iqs =0
 
integer iqg =0
 
integer diagnostic_generation_frequency
 

Detailed Description

Forcing, both subsidence and large scale.

Function/Subroutine Documentation

◆ apply_subsidence_to_flow_fields()

subroutine forcing_mod::apply_subsidence_to_flow_fields ( type(model_state_type), intent(inout)  current_state)
private

Definition at line 959 of file forcing.F90.

960  type(model_state_type), intent(inout) :: current_state
961 
962  integer :: k
963  real(kind=default_precision) :: usub, vsub
964 
965 
966  if (l_subs_local_theta)then ! Use local gradients not global means
967  u_profile(:)=current_state%zu%data(:,current_state%column_local_y,current_state%column_local_x)
968  v_profile(:)=current_state%zv%data(:,current_state%column_local_y,current_state%column_local_x)
969  else
970  u_profile(:)=current_state%global_grid%configuration%vertical%olzubar(:)
971  v_profile(:)=current_state%global_grid%configuration%vertical%olzvbar(:)
972  end if
973 
974  do k=2,current_state%local_grid%size(z_index)-1
975 #ifdef U_ACTIVE
976  usub = 2.0 * (current_state%global_grid%configuration%vertical%w_subs(k-1)* &
977  current_state%global_grid%configuration%vertical%tzc1(k)*(u_profile(k)-u_profile(k-1)) &
978  + current_state%global_grid%configuration%vertical%w_subs(k)* &
979  current_state%global_grid%configuration%vertical%tzc2(k)* &
980  (u_profile(k+1) - u_profile(k)))
981  current_state%su%data(k,current_state%column_local_y,current_state%column_local_x) = &
982  current_state%su%data(k,current_state%column_local_y,current_state%column_local_x) - usub
983  du_subs_profile_diag(k) = du_subs_profile_diag(k) - usub
984 #endif
985 #ifdef V_ACTIVE
986  vsub = 2.0 * (current_state%global_grid%configuration%vertical%w_subs(k-1)* &
987  current_state%global_grid%configuration%vertical%tzc1(k)*(v_profile(k)-v_profile(k-1)) &
988  + current_state%global_grid%configuration%vertical%w_subs(k)* &
989  current_state%global_grid%configuration%vertical%tzc2(k)* &
990  (v_profile(k+1) - v_profile(k)))
991  current_state%sv%data(k,current_state%column_local_y,current_state%column_local_x) = &
992  current_state%sv%data(k,current_state%column_local_y,current_state%column_local_x) - vsub
993  dv_subs_profile_diag(k) = dv_subs_profile_diag(k) - vsub
994 #endif
995  end do
996  k=current_state%local_grid%size(z_index)
997 #ifdef U_ACTIVE
998  usub = 2.0 * (current_state%global_grid%configuration%vertical%w_subs(k-1)* &
999  current_state%global_grid%configuration%vertical%tzc1(k)*(u_profile(k)-u_profile(k-1)))
1000  current_state%su%data(k,current_state%column_local_y,current_state%column_local_x) = &
1001  current_state%su%data(k,current_state%column_local_y,current_state%column_local_x) - usub
1002  du_subs_profile_diag(k) = du_subs_profile_diag(k) - usub
1003 #endif
1004 #ifdef V_ACTIVE
1005  vsub = 2.0 * (current_state%global_grid%configuration%vertical%w_subs(k-1)* &
1006  current_state%global_grid%configuration%vertical%tzc1(k)*(v_profile(k)-v_profile(k-1)))
1007  current_state%sv%data(k,current_state%column_local_y,current_state%column_local_x) = &
1008  current_state%sv%data(k,current_state%column_local_y,current_state%column_local_x) - vsub
1009  dv_subs_profile_diag(k) = dv_subs_profile_diag(k) - vsub
1010 #endif
Here is the caller graph for this function:

◆ apply_subsidence_to_q_fields()

subroutine forcing_mod::apply_subsidence_to_q_fields ( type(model_state_type), intent(inout)  current_state)
private

Definition at line 1046 of file forcing.F90.

1047  type(model_state_type), intent(inout) :: current_state
1048 
1049  integer :: k, n
1050  real(kind=default_precision) :: qsub
1051 
1052 
1053  do n=1,current_state%number_q_fields
1054  if (l_subs_local_q)then ! Use local gradients not global means
1055  q_profile(:)=current_state%zq(n)%data(:,current_state%column_local_y,current_state%column_local_x)
1056  else
1057  q_profile(:)=current_state%global_grid%configuration%vertical%olzqbar(:,n)
1058  end if
1059  do k=2,current_state%local_grid%size(z_index)-1
1060  qsub = current_state%global_grid%configuration%vertical%w_subs(k)* &
1061  (q_profile(k+1) - q_profile(k))* &
1062  current_state%global_grid%configuration%vertical%rdzn(k+1)
1063  current_state%sq(n)%data(k,current_state%column_local_y,current_state%column_local_x) = &
1064  current_state%sq(n)%data(k,current_state%column_local_y,current_state%column_local_x) - qsub
1065  dq_subs_profile_diag(k,n) = dq_subs_profile_diag(k,n) - qsub
1066  end do
1067  k=current_state%local_grid%size(z_index)
1068  qsub = current_state%global_grid%configuration%vertical%w_subs(k)* &
1069  (q_profile(k) - q_profile(k-1))* &
1070  current_state%global_grid%configuration%vertical%rdzn(k)
1071  current_state%sq(n)%data(k,current_state%column_local_y,current_state%column_local_x) = &
1072  current_state%sq(n)%data(k,current_state%column_local_y,current_state%column_local_x) - qsub
1073  dq_subs_profile_diag(k,n) = dq_subs_profile_diag(k,n) - qsub
1074  end do
Here is the caller graph for this function:

◆ apply_subsidence_to_theta()

subroutine forcing_mod::apply_subsidence_to_theta ( type(model_state_type), intent(inout)  current_state)
private

Definition at line 1013 of file forcing.F90.

1014  type(model_state_type), intent(inout) :: current_state
1015 
1016  integer :: k
1017  real(kind=default_precision) :: thsub
1018 
1019  if (l_subs_local_theta)then ! Use local gradients not global means
1020  theta_profile(:)=current_state%zth%data(:,current_state%column_local_y,current_state%column_local_x) &
1021  + current_state%global_grid%configuration%vertical%thref(:)
1022  else
1023  theta_profile(:)=current_state%global_grid%configuration%vertical%olzthbar(:) &
1024  + current_state%global_grid%configuration%vertical%thref(:)
1025  end if
1026 
1027  do k=2,current_state%local_grid%size(z_index)-1
1028  thsub = current_state%global_grid%configuration%vertical%w_subs(k)* &
1029  (theta_profile(k+1) - theta_profile(k))* &
1030  current_state%global_grid%configuration%vertical%rdzn(k+1)
1031 
1032  current_state%sth%data(k,current_state%column_local_y,current_state%column_local_x) = &
1033  current_state%sth%data(k,current_state%column_local_y,current_state%column_local_x) - thsub
1034  dtheta_subs_profile_diag(k) = dtheta_subs_profile_diag(k) - thsub
1035  end do
1036  k=current_state%local_grid%size(z_index)
1037  thsub = current_state%global_grid%configuration%vertical%w_subs(k)* &
1038  (theta_profile(k) - theta_profile(k-1))* &
1039  current_state%global_grid%configuration%vertical%rdzn(k)
1040 
1041  current_state%sth%data(k,current_state%column_local_y,current_state%column_local_x) = &
1042  current_state%sth%data(k,current_state%column_local_y,current_state%column_local_x) - thsub
1043  dtheta_subs_profile_diag(k) = dtheta_subs_profile_diag(k) - thsub
Here is the caller graph for this function:

◆ apply_time_independent_forcing_to_q()

subroutine forcing_mod::apply_time_independent_forcing_to_q ( type(model_state_type), intent(inout)  current_state)
private

Definition at line 1108 of file forcing.F90.

1109  type(model_state_type), intent(inout) :: current_state
1110 
1111  integer :: n, k
1112  real(kind=default_precision) :: dtm_scale
1113 
1114  do n=1,current_state%number_q_fields
1115  if (current_state%l_forceq(n))then
1116  if (constant_forcing_type_q==tendency)then
1117  dtm_scale=1.0_default_precision
1118  else ! constant_forcing_type_q==(RELAXATION or INCREMENT)
1119  dtm_scale=1.0_default_precision/forcing_timescale_q
1120  end if
1121 
1122  if (constant_forcing_type_q==relaxation)then
1123  dq_profile(:)=dtm_scale * (current_state%global_grid%configuration%vertical%q_force(:,n) - &
1124  current_state%zq(n)%data(:,current_state%column_local_y,current_state%column_local_x))
1125  else ! constant_forcing_type_q==(TENDENCY or INCREMENT)
1126  dq_profile(:)=dtm_scale * current_state%global_grid%configuration%vertical%q_force(:,n)
1127  end if
1128 
1129  dq_profile_diag(:,n)=dq_profile_diag(:,n)+dq_profile
1130 
1131  do k=2,current_state%local_grid%size(z_index)-1
1132  current_state%sq(n)%data(k,current_state%column_local_y,current_state%column_local_x) = &
1133  current_state%sq(n)%data(k,current_state%column_local_y,current_state%column_local_x) &
1134  + dq_profile(k)
1135  end do
1136  end if
1137  end do
1138 
Here is the caller graph for this function:

◆ apply_time_independent_forcing_to_theta()

subroutine forcing_mod::apply_time_independent_forcing_to_theta ( type(model_state_type), intent(inout)  current_state)
private

Definition at line 1077 of file forcing.F90.

1078  type(model_state_type), intent(inout) :: current_state
1079 
1080  integer :: k
1081  real(kind=default_precision) :: dtm_scale
1082 
1083  if (constant_forcing_type_theta==tendency)then
1084  dtm_scale=1.0_default_precision
1085  else ! constant_forcing_type_theta==(RELAXATION or INCREMENT)
1086  dtm_scale=1.0_default_precision/forcing_timescale_theta
1087  end if
1088 
1089  if (constant_forcing_type_theta==relaxation)then
1090  dtheta_profile(:)=dtm_scale * (current_state%global_grid%configuration%vertical%theta_force(:) - &
1091  current_state%zth%data(:,current_state%column_local_y,current_state%column_local_x) - &
1092  current_state%global_grid%configuration%vertical%thref(:))
1093  else ! constant_forcing_type_theta==(TENDENCY or INCREMENT)
1094  dtheta_profile(:)=dtm_scale * current_state%global_grid%configuration%vertical%theta_force(:)
1095  end if
1096 
1097 
1098  dtheta_profile_diag=dtheta_profile_diag+(dtheta_profile)
1099 
1100  do k=2,current_state%local_grid%size(z_index)-1
1101  current_state%sth%data(k,current_state%column_local_y,current_state%column_local_x) = &
1102  current_state%sth%data(k,current_state%column_local_y,current_state%column_local_x) &
1103  + dtheta_profile(k)
1104  end do
1105 
Here is the caller graph for this function:

◆ apply_time_independent_forcing_to_u()

subroutine forcing_mod::apply_time_independent_forcing_to_u ( type(model_state_type), intent(inout)  current_state)
private

Definition at line 1141 of file forcing.F90.

1142  type(model_state_type), intent(inout) :: current_state
1143 
1144  integer :: k
1145  real(kind=default_precision) :: dtm_scale
1146 
1147  if (constant_forcing_type_u==tendency)then
1148  dtm_scale=1.0_default_precision
1149  else ! constant_forcing_type_u==(RELAXATION or INCREMENT)
1150  dtm_scale=1.0_default_precision/forcing_timescale_u
1151  end if
1152 
1153  if (constant_forcing_type_u==relaxation)then
1154  du_profile(:)=dtm_scale * (current_state%global_grid%configuration%vertical%u_force(:) - &
1155  current_state%global_grid%configuration%vertical%olzubar(:))
1156  else ! constant_forcing_type_u==(TENDENCY or INCREMENT)
1157  du_profile(:)=dtm_scale * current_state%global_grid%configuration%vertical%u_force(:)
1158  end if
1159 
1160  du_profile_diag=du_profile_diag+du_profile
1161 
1162  do k=2,current_state%local_grid%size(z_index)-1
1163  current_state%su%data(k,current_state%column_local_y,current_state%column_local_x) = &
1164  current_state%su%data(k,current_state%column_local_y,current_state%column_local_x) &
1165  + du_profile(k)
1166  end do
1167 
Here is the caller graph for this function:

◆ apply_time_independent_forcing_to_v()

subroutine forcing_mod::apply_time_independent_forcing_to_v ( type(model_state_type), intent(inout)  current_state)
private

Definition at line 1170 of file forcing.F90.

1171  type(model_state_type), intent(inout) :: current_state
1172 
1173  integer :: k
1174  real(kind=default_precision) :: dtm_scale
1175 
1176  if (constant_forcing_type_v==tendency)then
1177  dtm_scale=1.0_default_precision
1178  else ! constant_forcing_type_v==(RELAXATION or INCREMENT)
1179  dtm_scale=1.0_default_precision/forcing_timescale_v
1180  end if
1181 
1182  if (constant_forcing_type_v==relaxation)then
1183  dv_profile(:)=dtm_scale * (current_state%global_grid%configuration%vertical%v_force(:) - &
1184  current_state%global_grid%configuration%vertical%olzvbar(:) )
1185  else ! constant_forcing_type_v==(TENDENCY or INCREMENT)
1186  dv_profile(:)=dtm_scale * current_state%global_grid%configuration%vertical%v_force(:)
1187  end if
1188 
1189  dv_profile_diag=dv_profile_diag+dv_profile
1190 
1191  do k=2,current_state%local_grid%size(z_index)-1
1192  current_state%sv%data(k,current_state%column_local_y,current_state%column_local_x) = &
1193  current_state%sv%data(k,current_state%column_local_y,current_state%column_local_x) &
1194  + dv_profile(k)
1195  end do
Here is the caller graph for this function:

◆ check_forcing_status()

subroutine forcing_mod::check_forcing_status ( integer, intent(in)  status)
private

Will check a NetCDF status and write to log_log error any decoded statuses.

Parameters
statusThe NetCDF status flag. This is a copy of the routine called in setfluxlook

Definition at line 1394 of file forcing.F90.

1395  integer, intent(in) :: status
1396 
1397  if (status /= nf90_noerr) then
1398  call log_log(log_error, "NetCDF returned error code of "//trim(nf90_strerror(status)))
1399  end if
Here is the caller graph for this function:

◆ compute_component_tendencies()

subroutine forcing_mod::compute_component_tendencies ( type(model_state_type), intent(inout), target  current_state,
integer, intent(in)  cxn,
integer, intent(in)  cyn,
integer, intent(in)  txn,
integer, intent(in)  tyn 
)
private

Computation of component tendencies.

Parameters
current_stateCurrent model state
cxnThe current slice, x, index
cynThe current column, y, index.
txntarget_x_index
tyntarget_y_index

Definition at line 1299 of file forcing.F90.

1300  type(model_state_type), target, intent(inout) :: current_state
1301  integer, intent(in) :: cxn, cyn, txn, tyn
1302 
1303  ! Calculate change in tendency due to component
1304  if (l_tend_3d_u) then
1305  tend_3d_u(:,tyn,txn)=current_state%su%data(:,cyn,cxn) - tend_3d_u(:,tyn,txn)
1306  endif
1307  if (l_tend_3d_v) then
1308  tend_3d_v(:,tyn,txn)=current_state%sv%data(:,cyn,cxn) - tend_3d_v(:,tyn,txn)
1309  endif
1310  if (l_tend_3d_th) then
1311  tend_3d_th(:,tyn,txn)=current_state%sth%data(:,cyn,cxn) - tend_3d_th(:,tyn,txn)
1312  endif
1313  if (l_tend_3d_qv) then
1314  tend_3d_qv(:,tyn,txn)=current_state%sq(iqv)%data(:,cyn,cxn) - tend_3d_qv(:,tyn,txn)
1315  endif
1316  if (l_tend_3d_ql) then
1317  tend_3d_ql(:,tyn,txn)=current_state%sq(iql)%data(:,cyn,cxn) - tend_3d_ql(:,tyn,txn)
1318  endif
1319  if (l_tend_3d_qi) then
1320  tend_3d_qi(:,tyn,txn)=current_state%sq(iqi)%data(:,cyn,cxn) - tend_3d_qi(:,tyn,txn)
1321  endif
1322  if (l_tend_3d_qr) then
1323  tend_3d_qr(:,tyn,txn)=current_state%sq(iqr)%data(:,cyn,cxn) - tend_3d_qr(:,tyn,txn)
1324  endif
1325  if (l_tend_3d_qs) then
1326  tend_3d_qs(:,tyn,txn)=current_state%sq(iqs)%data(:,cyn,cxn) - tend_3d_qs(:,tyn,txn)
1327  endif
1328  if (l_tend_3d_qg) then
1329  tend_3d_qg(:,tyn,txn)=current_state%sq(iqg)%data(:,cyn,cxn) - tend_3d_qg(:,tyn,txn)
1330  endif
1331  if (l_tend_3d_tabs) then
1332  tend_3d_tabs(:,tyn,txn)= &
1333  current_state%sth%data(:,cyn,cxn) * current_state%global_grid%configuration%vertical%rprefrcp(:) &
1334  - tend_3d_tabs(:,tyn,txn)
1335  endif
1336 
1337  ! Add local tendency fields to the profile total
1338  if (l_tend_pr_tot_u) then
1339  tend_pr_tot_u(:)=tend_pr_tot_u(:) + tend_3d_u(:,tyn,txn)
1340  endif
1341  if (l_tend_pr_tot_v) then
1342  tend_pr_tot_v(:)=tend_pr_tot_v(:) + tend_3d_v(:,tyn,txn)
1343  endif
1344  if (l_tend_pr_tot_th) then
1345  tend_pr_tot_th(:)=tend_pr_tot_th(:) + tend_3d_th(:,tyn,txn)
1346  endif
1347  if (l_tend_pr_tot_qv) then
1348  tend_pr_tot_qv(:)=tend_pr_tot_qv(:) + tend_3d_qv(:,tyn,txn)
1349  endif
1350  if (l_tend_pr_tot_ql) then
1351  tend_pr_tot_ql(:)=tend_pr_tot_ql(:) + tend_3d_ql(:,tyn,txn)
1352  endif
1353  if (l_tend_pr_tot_qi) then
1354  tend_pr_tot_qi(:)=tend_pr_tot_qi(:) + tend_3d_qi(:,tyn,txn)
1355  endif
1356  if (l_tend_pr_tot_qr) then
1357  tend_pr_tot_qr(:)=tend_pr_tot_qr(:) + tend_3d_qr(:,tyn,txn)
1358  endif
1359  if (l_tend_pr_tot_qs) then
1360  tend_pr_tot_qs(:)=tend_pr_tot_qs(:) + tend_3d_qs(:,tyn,txn)
1361  endif
1362  if (l_tend_pr_tot_qg) then
1363  tend_pr_tot_qg(:)=tend_pr_tot_qg(:) + tend_3d_qg(:,tyn,txn)
1364  endif
1365  if (l_tend_pr_tot_tabs) then
1366  tend_pr_tot_tabs(:)=tend_pr_tot_tabs(:) + tend_3d_tabs(:,tyn,txn)
1367  endif
1368 
Here is the caller graph for this function:

◆ field_information_retrieval_callback()

subroutine forcing_mod::field_information_retrieval_callback ( type(model_state_type), intent(inout), target  current_state,
character(len=*), intent(in)  name,
type(component_field_information_type), intent(out)  field_information 
)
private

Field information retrieval callback, this returns information for a specific components published field.

Parameters
current_stateCurrent model state
nameThe name of the field to retrieve information for
field_informationPopulated with information about the field

Definition at line 181 of file forcing.F90.

182  type(model_state_type), target, intent(inout) :: current_state
183  character(len=*), intent(in) :: name
184  type(component_field_information_type), intent(out) :: field_information
185  integer :: strcomp
186 
187  field_information%field_type=component_array_field_type
188  field_information%data_type=component_double_data_type
189  field_information%number_dimensions=1
190  field_information%dimension_sizes(1)=current_state%local_grid%size(z_index)
191 
192  if (name .eq. "u_subsidence") then
193  field_information%enabled=current_state%u%active .and. l_subs_pl_theta .and. &
194  allocated(current_state%global_grid%configuration%vertical%olzubar)
195  else if (name .eq. "v_subsidence") then
196  field_information%enabled=current_state%v%active .and. l_subs_pl_theta .and. &
197  allocated(current_state%global_grid%configuration%vertical%olzvbar)
198  else if (name .eq. "th_subsidence") then
199  field_information%enabled=current_state%th%active .and. l_subs_pl_theta .and. &
200  allocated(current_state%global_grid%configuration%vertical%olzthbar)
201  else if (name .eq. "vapour_mmr_subsidence" .or. name .eq. "vapour_mmr_subsidence" .or. &
202  name .eq. "cloud_mmr_subsidence" .or. name .eq. "cloud_mmr_subsidence" ) then
203  field_information%enabled=.not. current_state%passive_q .and. &
204  current_state%number_q_fields .gt. 0 .and. l_subs_pl_q .and. &
205  allocated(current_state%global_grid%configuration%vertical%olzqbar)
206  else if (name .eq. "rain_mmr_subsidence" ) then
207  field_information%enabled=current_state%rain_water_mixing_ratio_index .gt. 0 .and. &
208  allocated(current_state%global_grid%configuration%vertical%olzqbar)
209  else if (name .eq. "ice_mmr_subsidence" ) then
210  field_information%enabled= current_state%ice_water_mixing_ratio_index .gt. 0 .and. &
211  allocated(current_state%global_grid%configuration%vertical%olzqbar)
212  else if (name .eq. "snow_mmr_subsidence" ) then
213  field_information%enabled= current_state%snow_water_mixing_ratio_index .gt. 0 .and. &
214  allocated(current_state%global_grid%configuration%vertical%olzqbar)
215  else if (name .eq. "graupel_mmr_subsidence" ) then
216  field_information%enabled= current_state%graupel_water_mixing_ratio_index .gt. 0 .and. &
217  allocated(current_state%global_grid%configuration%vertical%olzqbar)
218  else if (name .eq. "u_large_scale") then
219  field_information%enabled=current_state%u%active .and. l_constant_forcing_u
220  else if (name .eq. "v_large_scale") then
221  field_information%enabled=current_state%v%active .and. l_constant_forcing_v
222  else if (name .eq. "th_large_scale") then
223  field_information%enabled=current_state%th%active .and. l_constant_forcing_theta
224  else if (name .eq. "vapour_mmr_large_scale" .or. name .eq. "vapour_mmr_large_scale" .or. &
225  name .eq. "cloud_mmr_large_scale" .or. name .eq. "cloud_mmr_large_scale" ) then
226  field_information%enabled=.not. current_state%passive_q .and. &
227  current_state%number_q_fields .gt. 0 .and. l_subs_pl_q .and. &
228  allocated(current_state%global_grid%configuration%vertical%olzqbar)
229  else if (name .eq. "rain_mmr_large_scale" ) then
230  field_information%enabled=current_state%rain_water_mixing_ratio_index .gt. 0 .and. &
231  allocated(current_state%global_grid%configuration%vertical%olzqbar)
232  else if (name .eq. "ice_mmr_large_scale" ) then
233  field_information%enabled= current_state%ice_water_mixing_ratio_index .gt. 0 .and. &
234  allocated(current_state%global_grid%configuration%vertical%olzqbar)
235  else if (name .eq. "snow_mmr_large_scale" ) then
236  field_information%enabled= current_state%snow_water_mixing_ratio_index .gt. 0 .and. &
237  allocated(current_state%global_grid%configuration%vertical%olzqbar)
238  else if (name .eq. "graupel_mmr_large_scale" ) then
239  field_information%enabled= current_state%graupel_water_mixing_ratio_index .gt. 0 .and. &
240  allocated(current_state%global_grid%configuration%vertical%olzqbar)
241  end if
242 
243  ! Field information for 3d
244  strcomp=index(name, "forcing_3d_local")
245  if (strcomp .ne. 0) then
246  field_information%field_type=component_array_field_type
247  field_information%number_dimensions=3
248  field_information%dimension_sizes(1)=current_state%local_grid%size(z_index)
249  field_information%dimension_sizes(2)=current_state%local_grid%size(y_index)
250  field_information%dimension_sizes(3)=current_state%local_grid%size(x_index)
251  field_information%data_type=component_double_data_type
252 
253  if (name .eq. "tend_u_forcing_3d_local") then
254  field_information%enabled=l_tend_3d_u
255  else if (name .eq. "tend_v_forcing_3d_local") then
256  field_information%enabled=l_tend_3d_v
257  else if (name .eq. "tend_th_forcing_3d_local") then
258  field_information%enabled=l_tend_3d_th
259  else if (name .eq. "tend_qv_forcing_3d_local") then
260  field_information%enabled=l_tend_3d_qv
261  else if (name .eq. "tend_ql_forcing_3d_local") then
262  field_information%enabled=l_tend_3d_ql
263  else if (name .eq. "tend_qi_forcing_3d_local") then
264  field_information%enabled=l_tend_3d_qi
265  else if (name .eq. "tend_qr_forcing_3d_local") then
266  field_information%enabled=l_tend_3d_qr
267  else if (name .eq. "tend_qs_forcing_3d_local") then
268  field_information%enabled=l_tend_3d_qs
269  else if (name .eq. "tend_qg_forcing_3d_local") then
270  field_information%enabled=l_tend_3d_qg
271  else if (name .eq. "tend_tabs_forcing_3d_local") then
272  field_information%enabled=l_tend_3d_tabs
273  else
274  field_information%enabled=.true.
275  end if
276 
277  end if !end 3d check
278 
279  ! Field information for profiles
280  strcomp=index(name, "forcing_profile_total_local")
281  if (strcomp .ne. 0) then
282  field_information%field_type=component_array_field_type
283  field_information%number_dimensions=1
284  field_information%dimension_sizes(1)=current_state%local_grid%size(z_index)
285  field_information%data_type=component_double_data_type
286 
287  if (name .eq. "tend_u_forcing_profile_total_local") then
288  field_information%enabled=l_tend_pr_tot_u
289  else if (name .eq. "tend_v_forcing_profile_total_local") then
290  field_information%enabled=l_tend_pr_tot_v
291  else if (name .eq. "tend_th_forcing_profile_total_local") then
292  field_information%enabled=l_tend_pr_tot_th
293  else if (name .eq. "tend_qv_forcing_profile_total_local") then
294  field_information%enabled=l_tend_pr_tot_qv
295  else if (name .eq. "tend_ql_forcing_profile_total_local") then
296  field_information%enabled=l_tend_pr_tot_ql
297  else if (name .eq. "tend_qi_forcing_profile_total_local") then
298  field_information%enabled=l_tend_pr_tot_qi
299  else if (name .eq. "tend_qr_forcing_profile_total_local") then
300  field_information%enabled=l_tend_pr_tot_qr
301  else if (name .eq. "tend_qs_forcing_profile_total_local") then
302  field_information%enabled=l_tend_pr_tot_qs
303  else if (name .eq. "tend_qg_forcing_profile_total_local") then
304  field_information%enabled=l_tend_pr_tot_qg
305  else if (name .eq. "tend_tabs_forcing_profile_total_local") then
306  field_information%enabled=l_tend_pr_tot_tabs
307  else
308  field_information%enabled=.true.
309  end if
310 
311  end if !end profile check
312 
Here is the caller graph for this function:

◆ field_value_retrieval_callback()

subroutine forcing_mod::field_value_retrieval_callback ( type(model_state_type), intent(inout), target  current_state,
character(len=*), intent(in)  name,
type(component_field_value_type), intent(out)  field_value 
)
private

Field value retrieval callback, this returns the value of a specific published field.

Parameters
current_stateCurrent model state
nameThe name of the field to retrieve the value for
field_valuePopulated with the value of the field

Definition at line 319 of file forcing.F90.

320  type(model_state_type), target, intent(inout) :: current_state
321  character(len=*), intent(in) :: name
322  type(component_field_value_type), intent(out) :: field_value
323 
324  integer :: k, n, column_size
325 
326  column_size=current_state%local_grid%size(z_index)
327 
328  ! subsidence diagnostics
329  if (name .eq. "u_subsidence") then
330  allocate(field_value%real_1d_array(column_size))
331  field_value%real_1d_array(:)=du_subs_profile_diag(:)
332  else if (name .eq. "v_subsidence") then
333  allocate(field_value%real_1d_array(column_size))
334  field_value%real_1d_array(:)=dv_subs_profile_diag(:)
335  else if (name .eq. "th_subsidence") then
336  allocate(field_value%real_1d_array(column_size))
337  field_value%real_1d_array(:)= dtheta_subs_profile_diag(:)
338  else if (name .eq. "vapour_mmr_subsidence") then
339  allocate(field_value%real_1d_array(column_size))
340  field_value%real_1d_array(:)=dq_subs_profile_diag(:,iqv)
341  else if (name .eq. "cloud_mmr_subsidence") then
342  allocate(field_value%real_1d_array(column_size))
343  field_value%real_1d_array(:)=dq_subs_profile_diag(:,iql)
344  else if (name .eq. "rain_mmr_subsidence") then
345  allocate(field_value%real_1d_array(column_size))
346  field_value%real_1d_array(:)=dq_subs_profile_diag(:,iqr)
347  else if (name .eq. "ice_mmr_subsidence") then
348  allocate(field_value%real_1d_array(column_size))
349  field_value%real_1d_array(:)=dq_subs_profile_diag(:,iqi)
350  else if (name .eq. "snow_mmr_subsidence") then
351  allocate(field_value%real_1d_array(column_size))
352  field_value%real_1d_array(:)=dq_subs_profile_diag(:,iqs)
353  else if (name .eq. "graupel_mmr_subsidence") then
354  allocate(field_value%real_1d_array(column_size))
355  field_value%real_1d_array(:)=dq_subs_profile_diag(:,iqg)
356  ! Large-scale forcing diagnostics
357  else if (name .eq. "u_large_scale") then
358  allocate(field_value%real_1d_array(column_size))
359  field_value%real_1d_array=get_averaged_diagnostics(current_state, du_profile_diag)
360  else if (name .eq. "v_large_scale") then
361  allocate(field_value%real_1d_array(column_size))
362  field_value%real_1d_array=get_averaged_diagnostics(current_state, dv_profile_diag)
363  else if (name .eq. "th_large_scale") then
364  allocate(field_value%real_1d_array(column_size))
365  field_value%real_1d_array=dtheta_profile_diag
366  else if (name .eq. "vapour_mmr_large_scale") then
367  allocate(field_value%real_1d_array(column_size))
368  field_value%real_1d_array(:)=get_averaged_diagnostics(current_state, dq_profile_diag(:,iqv))
369  else if (name .eq. "cloud_mmr_large_scale") then
370  allocate(field_value%real_1d_array(column_size))
371  field_value%real_1d_array(:)=get_averaged_diagnostics(current_state, dq_profile_diag(:,iql))
372  else if (name .eq. "rain_mmr_large_scale") then
373  allocate(field_value%real_1d_array(column_size))
374  field_value%real_1d_array(:)=get_averaged_diagnostics(current_state, dq_profile_diag(:,iqr))
375  else if (name .eq. "ice_mmr_large_scale") then
376  allocate(field_value%real_1d_array(column_size))
377  field_value%real_1d_array(:)=get_averaged_diagnostics(current_state, dq_profile_diag(:,iqi))
378  else if (name .eq. "snow_mmr_large_scale") then
379  allocate(field_value%real_1d_array(column_size))
380  field_value%real_1d_array(:)=get_averaged_diagnostics(current_state, dq_profile_diag(:,iqs))
381  else if (name .eq. "graupel_mmr_large_scale") then
382  allocate(field_value%real_1d_array(column_size))
383  field_value%real_1d_array(:)=get_averaged_diagnostics(current_state, dq_profile_diag(:,iqg))
384 
385  ! 3d Tendency Fields
386  else if (name .eq. "tend_u_forcing_3d_local" .and. allocated(tend_3d_u)) then
387  call set_published_field_value(field_value, real_3d_field=tend_3d_u)
388  else if (name .eq. "tend_v_forcing_3d_local" .and. allocated(tend_3d_v)) then
389  call set_published_field_value(field_value, real_3d_field=tend_3d_v)
390  else if (name .eq. "tend_th_forcing_3d_local" .and. allocated(tend_3d_th)) then
391  call set_published_field_value(field_value, real_3d_field=tend_3d_th)
392  else if (name .eq. "tend_qv_forcing_3d_local" .and. allocated(tend_3d_qv)) then
393  call set_published_field_value(field_value, real_3d_field=tend_3d_qv)
394  else if (name .eq. "tend_ql_forcing_3d_local" .and. allocated(tend_3d_ql)) then
395  call set_published_field_value(field_value, real_3d_field=tend_3d_ql)
396  else if (name .eq. "tend_qi_forcing_3d_local" .and. allocated(tend_3d_qi)) then
397  call set_published_field_value(field_value, real_3d_field=tend_3d_qi)
398  else if (name .eq. "tend_qr_forcing_3d_local" .and. allocated(tend_3d_qr)) then
399  call set_published_field_value(field_value, real_3d_field=tend_3d_qr)
400  else if (name .eq. "tend_qs_forcing_3d_local" .and. allocated(tend_3d_qs)) then
401  call set_published_field_value(field_value, real_3d_field=tend_3d_qs)
402  else if (name .eq. "tend_qg_forcing_3d_local" .and. allocated(tend_3d_qg)) then
403  call set_published_field_value(field_value, real_3d_field=tend_3d_qg)
404  else if (name .eq. "tend_tabs_forcing_3d_local" .and. allocated(tend_3d_tabs)) then
405  call set_published_field_value(field_value, real_3d_field=tend_3d_tabs)
406 
407  ! Profile Tendency Fields
408  else if (name .eq. "tend_u_forcing_profile_total_local" .and. allocated(tend_pr_tot_u)) then
409  call set_published_field_value(field_value, real_1d_field=tend_pr_tot_u)
410  else if (name .eq. "tend_v_forcing_profile_total_local" .and. allocated(tend_pr_tot_v)) then
411  call set_published_field_value(field_value, real_1d_field=tend_pr_tot_v)
412  else if (name .eq. "tend_th_forcing_profile_total_local" .and. allocated(tend_pr_tot_th)) then
413  call set_published_field_value(field_value, real_1d_field=tend_pr_tot_th)
414  else if (name .eq. "tend_qv_forcing_profile_total_local" .and. allocated(tend_pr_tot_qv)) then
415  call set_published_field_value(field_value, real_1d_field=tend_pr_tot_qv)
416  else if (name .eq. "tend_ql_forcing_profile_total_local" .and. allocated(tend_pr_tot_ql)) then
417  call set_published_field_value(field_value, real_1d_field=tend_pr_tot_ql)
418  else if (name .eq. "tend_qi_forcing_profile_total_local" .and. allocated(tend_pr_tot_qi)) then
419  call set_published_field_value(field_value, real_1d_field=tend_pr_tot_qi)
420  else if (name .eq. "tend_qr_forcing_profile_total_local" .and. allocated(tend_pr_tot_qr)) then
421  call set_published_field_value(field_value, real_1d_field=tend_pr_tot_qr)
422  else if (name .eq. "tend_qs_forcing_profile_total_local" .and. allocated(tend_pr_tot_qs)) then
423  call set_published_field_value(field_value, real_1d_field=tend_pr_tot_qs)
424  else if (name .eq. "tend_qg_forcing_profile_total_local" .and. allocated(tend_pr_tot_qg)) then
425  call set_published_field_value(field_value, real_1d_field=tend_pr_tot_qg)
426  else if (name .eq. "tend_tabs_forcing_profile_total_local" .and. allocated(tend_pr_tot_tabs)) then
427  call set_published_field_value(field_value, real_1d_field=tend_pr_tot_tabs)
428  end if
429 
Here is the call graph for this function:
Here is the caller graph for this function:

◆ finalisation_callback()

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

Finalises the component @current_state Current model state.

Definition at line 1200 of file forcing.F90.

1201  type(model_state_type), target, intent(inout) :: current_state
1202  deallocate(theta_profile, q_profile, u_profile, v_profile)
1203  deallocate(dtheta_profile, dq_profile, du_profile, dv_profile)
1204  if (allocated(du_profile_diag)) deallocate(du_profile_diag)
1205  if (allocated(dv_profile_diag)) deallocate(dv_profile_diag)
1206  if (allocated(dtheta_profile_diag)) deallocate(dtheta_profile_diag)
1207  if (allocated(dq_profile_diag)) deallocate(dq_profile_diag)
1208 
1209  if (allocated(tend_3d_u)) deallocate(tend_3d_u)
1210  if (allocated(tend_3d_v)) deallocate(tend_3d_v)
1211  if (allocated(tend_3d_th)) deallocate(tend_3d_th)
1212  if (allocated(tend_3d_qv)) deallocate(tend_3d_qv)
1213  if (allocated(tend_3d_ql)) deallocate(tend_3d_ql)
1214  if (allocated(tend_3d_qi)) deallocate(tend_3d_qi)
1215  if (allocated(tend_3d_qr)) deallocate(tend_3d_qr)
1216  if (allocated(tend_3d_qs)) deallocate(tend_3d_qs)
1217  if (allocated(tend_3d_qg)) deallocate(tend_3d_qg)
1218  if (allocated(tend_3d_tabs)) deallocate(tend_3d_tabs)
1219 
1220  if (allocated(tend_pr_tot_u)) deallocate(tend_pr_tot_u)
1221  if (allocated(tend_pr_tot_v)) deallocate(tend_pr_tot_v)
1222  if (allocated(tend_pr_tot_th)) deallocate(tend_pr_tot_th)
1223  if (allocated(tend_pr_tot_qv)) deallocate(tend_pr_tot_qv)
1224  if (allocated(tend_pr_tot_ql)) deallocate(tend_pr_tot_ql)
1225  if (allocated(tend_pr_tot_qi)) deallocate(tend_pr_tot_qi)
1226  if (allocated(tend_pr_tot_qr)) deallocate(tend_pr_tot_qr)
1227  if (allocated(tend_pr_tot_qs)) deallocate(tend_pr_tot_qs)
1228  if (allocated(tend_pr_tot_qg)) deallocate(tend_pr_tot_qg)
1229  if (allocated(tend_pr_tot_tabs)) deallocate(tend_pr_tot_tabs)
1230 
Here is the caller graph for this function:

◆ forcing_get_descriptor()

type(component_descriptor_type) function, public forcing_mod::forcing_get_descriptor

Provides the component descriptor for the core to register.

Returns
The descriptor describing this component

Definition at line 123 of file forcing.F90.

124  forcing_get_descriptor%name="forcing"
125  forcing_get_descriptor%version=0.1
126  forcing_get_descriptor%initialisation=>init_callback
127  forcing_get_descriptor%timestep=>timestep_callback
128  forcing_get_descriptor%finalisation=>finalisation_callback
129 
130  forcing_get_descriptor%field_value_retrieval=>field_value_retrieval_callback
131  forcing_get_descriptor%field_information_retrieval=>field_information_retrieval_callback
132  allocate(forcing_get_descriptor%published_fields(18+10+10))
133 
134  forcing_get_descriptor%published_fields(1)="u_subsidence"
135  forcing_get_descriptor%published_fields(2)="v_subsidence"
136  forcing_get_descriptor%published_fields(3)="th_subsidence"
137  forcing_get_descriptor%published_fields(4)="vapour_mmr_subsidence"
138  forcing_get_descriptor%published_fields(5)="cloud_mmr_subsidence"
139  forcing_get_descriptor%published_fields(6)="rain_mmr_subsidence"
140  forcing_get_descriptor%published_fields(7)="ice_mmr_subsidence"
141  forcing_get_descriptor%published_fields(8)="snow_mmr_subsidence"
142  forcing_get_descriptor%published_fields(9)="graupel_mmr_subsidence"
143  forcing_get_descriptor%published_fields(10)="u_large_scale"
144  forcing_get_descriptor%published_fields(11)="v_large_scale"
145  forcing_get_descriptor%published_fields(12)="th_large_scale"
146  forcing_get_descriptor%published_fields(13)="vapour_mmr_large_scale"
147  forcing_get_descriptor%published_fields(14)="cloud_mmr_large_scale"
148  forcing_get_descriptor%published_fields(15)="rain_mmr_large_scale"
149  forcing_get_descriptor%published_fields(16)="ice_mmr_large_scale"
150  forcing_get_descriptor%published_fields(17)="snow_mmr_large_scale"
151  forcing_get_descriptor%published_fields(18)="graupel_mmr_large_scale"
152 
153  forcing_get_descriptor%published_fields(18+1)= "tend_u_forcing_3d_local"
154  forcing_get_descriptor%published_fields(18+2)= "tend_v_forcing_3d_local"
155  forcing_get_descriptor%published_fields(18+3)= "tend_th_forcing_3d_local"
156  forcing_get_descriptor%published_fields(18+4)= "tend_qv_forcing_3d_local"
157  forcing_get_descriptor%published_fields(18+5)= "tend_ql_forcing_3d_local"
158  forcing_get_descriptor%published_fields(18+6)= "tend_qi_forcing_3d_local"
159  forcing_get_descriptor%published_fields(18+7)= "tend_qr_forcing_3d_local"
160  forcing_get_descriptor%published_fields(18+8)= "tend_qs_forcing_3d_local"
161  forcing_get_descriptor%published_fields(18+9)= "tend_qg_forcing_3d_local"
162  forcing_get_descriptor%published_fields(18+10)="tend_tabs_forcing_3d_local"
163 
164  forcing_get_descriptor%published_fields(18+10+1)= "tend_u_forcing_profile_total_local"
165  forcing_get_descriptor%published_fields(18+10+2)= "tend_v_forcing_profile_total_local"
166  forcing_get_descriptor%published_fields(18+10+3)= "tend_th_forcing_profile_total_local"
167  forcing_get_descriptor%published_fields(18+10+4)= "tend_qv_forcing_profile_total_local"
168  forcing_get_descriptor%published_fields(18+10+5)= "tend_ql_forcing_profile_total_local"
169  forcing_get_descriptor%published_fields(18+10+6)= "tend_qi_forcing_profile_total_local"
170  forcing_get_descriptor%published_fields(18+10+7)= "tend_qr_forcing_profile_total_local"
171  forcing_get_descriptor%published_fields(18+10+8)= "tend_qs_forcing_profile_total_local"
172  forcing_get_descriptor%published_fields(18+10+9)= "tend_qg_forcing_profile_total_local"
173  forcing_get_descriptor%published_fields(18+10+10)="tend_tabs_forcing_profile_total_local"
174 
Here is the call graph for this function:

◆ get_averaged_diagnostics()

real(kind=default_precision) function, dimension(size(diagnostics_summed)) forcing_mod::get_averaged_diagnostics ( type(model_state_type), intent(inout), target  current_state,
real(kind=default_precision), dimension(:), intent(in)  diagnostics_summed 
)
private

Averages some diagnostic values across all local horizontal points.

Parameters
current_stateCurrent model state
diagnostics_summedThe diagnostics, which are summed up values from each local point, which need to be averaged
Returns
The averaged diagnostics across all local horizontal points

Definition at line 1237 of file forcing.F90.

1238  type(model_state_type), target, intent(inout) :: current_state
1239  real(kind=default_precision), dimension(:), intent(in) :: diagnostics_summed
1240  real(kind=default_precision), dimension(size(diagnostics_summed)) :: get_averaged_diagnostics
1241 
1242  integer :: horizontal_points
1243 
1244  horizontal_points=current_state%local_grid%size(x_index) * current_state%local_grid%size(y_index)
1245 
1246  get_averaged_diagnostics(:)=diagnostics_summed(:)/horizontal_points
Here is the caller graph for this function:

◆ init_callback()

subroutine forcing_mod::init_callback ( type(model_state_type), intent(inout), target  current_state)
private

Initialises the forcing data structures.

Definition at line 433 of file forcing.F90.

434 
435  type(model_state_type), target, intent(inout) :: current_state
436 
437  integer :: nq_force ! The number of q fields apply large-scale time-independent forcing
438  integer :: nzq ! The number of input levels for subsidence/divergence profile
439  integer :: i,n ! loop counters
440  integer :: iq ! temporary q varible index
441 
442  integer :: ncid ! id for the netcdf file
443  integer :: time_dim ! number of elements in time variable, read from input file
444  integer :: z_dim ! number of elements in height variable, read from input file
445 
446  ! Input arrays from config (always 1D) - subsidence profile
447  real(kind=default_precision), dimension(:), allocatable :: f_subs_pl ! subsidence node for q variables
448  real(kind=default_precision), dimension(:), allocatable :: z_subs_pl ! subsidence node height values for q variables
449 
450  ! Input arrays from config (always 1D) - time-independent forcing
451  real(kind=default_precision), dimension(:, :), allocatable :: f_force_pl_q ! Forcing values for q variables
452  real(kind=default_precision), dimension(:), allocatable :: z_force_pl_q ! Forcing height values for q variables
453  real(kind=default_precision), dimension(:), allocatable :: f_force_pl_theta ! Forcing values for theta variable
454  real(kind=default_precision), dimension(:), allocatable :: z_force_pl_theta ! Forcing height values for theta variable
455  real(kind=default_precision), dimension(:), allocatable :: f_force_pl_u ! Forcing values for u variable
456  real(kind=default_precision), dimension(:), allocatable :: z_force_pl_u ! Forcing height values for u variable
457  real(kind=default_precision), dimension(:), allocatable :: f_force_pl_v ! Forcing values for v variable
458  real(kind=default_precision), dimension(:), allocatable :: z_force_pl_v ! Forcing height values for v variabl
459 
460  ! Read from netcdf file if used - always 2D
461  real(kind=default_precision), dimension(:,:), allocatable :: f_subs_2d ! subsidence node for q variables
462 
463  integer :: subsidence_input_type=divergence ! Determines if we're reading in a subsidence velocity or divergence
464 
465  real(kind=default_precision), allocatable :: f_force_pl_q_tmp(:) !temporary 1D storage of forcing for q fields
466  real(kind=default_precision), allocatable :: zgrid(:) ! z grid to use in interpolation
467 
468  character(len=STRING_LENGTH), dimension(:), allocatable :: units_q_force ! units of q variable forcing
469  character(len=STRING_LENGTH) :: units_theta_force='unset' ! units of theta variable forcing
470  character(len=STRING_LENGTH) :: units_u_force='unset' ! units of theta variable forcing
471  character(len=STRING_LENGTH) :: units_v_force='unset' ! units of theta variable forcing
472 
473  logical :: convert_input_theta_from_temperature=.false. ! If .true. input forcing data is for temperature and should
474  ! be converted to theta (potential temerature).
475 
476  integer :: k
477  logical :: l_qdiag
478 
479  allocate(u_profile(current_state%local_grid%size(z_index)), &
480  v_profile(current_state%local_grid%size(z_index)), &
481  theta_profile(current_state%local_grid%size(z_index)), &
482  q_profile(current_state%local_grid%size(z_index)))
483 
484  allocate(dtheta_profile(current_state%local_grid%size(z_index)), &
485  dq_profile(current_state%local_grid%size(z_index)), &
486  du_profile(current_state%local_grid%size(z_index)), &
487  dv_profile(current_state%local_grid%size(z_index)))
488 
489  allocate(du_profile_diag(current_state%local_grid%size(z_index)), &
490  dv_profile_diag(current_state%local_grid%size(z_index)), &
491  dtheta_profile_diag(current_state%local_grid%size(z_index)), &
492  dq_profile_diag(current_state%local_grid%size(z_index), current_state%number_q_fields))
493 
494  allocate(du_subs_profile_diag(current_state%local_grid%size(z_index)), &
495  dv_subs_profile_diag(current_state%local_grid%size(z_index)), &
496  dtheta_subs_profile_diag(current_state%local_grid%size(z_index)), &
497  dq_subs_profile_diag(current_state%local_grid%size(z_index), current_state%number_q_fields))
498 
499  allocate(zgrid(current_state%local_grid%size(z_index)))
500 
501  ! assign microphysics indexes, needed for the diagnostic output
502  if (.not. current_state%passive_q .and. current_state%number_q_fields .gt. 0) then
503  iqv=get_q_index(standard_q_names%VAPOUR, 'forcing')
504  iql=get_q_index(standard_q_names%CLOUD_LIQUID_MASS, 'forcing')
505  endif
506  if (current_state%rain_water_mixing_ratio_index > 0) &
507  iqr = current_state%rain_water_mixing_ratio_index
508  if (current_state%ice_water_mixing_ratio_index > 0) &
509  iqi = current_state%ice_water_mixing_ratio_index
510  if (current_state%snow_water_mixing_ratio_index > 0) &
511  iqs = current_state%snow_water_mixing_ratio_index
512  if (current_state%graupel_water_mixing_ratio_index > 0) &
513  iqg = current_state%graupel_water_mixing_ratio_index
514 
515  ! Subsidence forcing initialization
516 
517  use_time_varying_subsidence= &
518  options_get_logical(current_state%options_database, "use_time_varying_subsidence")
519  l_subs_pl_theta=options_get_logical(current_state%options_database, "l_subs_pl_theta")
520  l_subs_pl_q=options_get_logical(current_state%options_database, "l_subs_pl_q")
521  subsidence_input_type=options_get_integer(current_state%options_database, "subsidence_input_type")
522  l_subs_local_theta=options_get_logical(current_state%options_database, "subsidence_local_theta")
523  l_subs_local_q=options_get_logical(current_state%options_database, "subsidence_local_q")
524 
525  input_file=options_get_string(current_state%options_database, "forcing_file")
526 
527  if ((l_subs_pl_theta .and. .not. l_subs_local_theta) .or. &
528  (l_subs_pl_q .and. .not. l_subs_local_q))then
529  if (.not. is_component_enabled(current_state%options_database, "mean_profiles")) then
530  call log_master_log(log_error, "subsidence requires the mean profiles component to be enabled")
531  end if
532  end if
533 
534  if (l_subs_pl_theta .or. l_subs_pl_q)then
535  ! Read in the input_file
536  if (trim(input_file)=='' .or. trim(input_file)=='None') then
537  if (.not. use_time_varying_subsidence) then
538  allocate(z_subs_pl(options_get_array_size(current_state%options_database, "z_subs_pl")), &
539  f_subs_pl(options_get_array_size(current_state%options_database, "f_subs_pl")))
540  call options_get_real_array(current_state%options_database, "z_subs_pl", z_subs_pl)
541  call options_get_real_array(current_state%options_database, "f_subs_pl", f_subs_pl)
542  ! Get profiles
543  zgrid=current_state%global_grid%configuration%vertical%z(:)
544  call piecewise_linear_1d(z_subs_pl(1:size(z_subs_pl)), f_subs_pl(1:size(f_subs_pl)), zgrid, &
545  current_state%global_grid%configuration%vertical%w_subs)
546  if (subsidence_input_type==divergence) then
547  current_state%global_grid%configuration%vertical%w_subs(:) = &
548  -1.0*current_state%global_grid%configuration%vertical%w_subs(:)*zgrid(:)
549  end if
550  else
551  call log_master_log(log_error, "timevarying forcing selected but no forcing file provided - STOP")
552  endif
553  deallocate(z_subs_pl, f_subs_pl)
554  else
555  if (use_time_varying_subsidence) then
556  call check_forcing_status(nf90_open(path = trim(input_file), mode = nf90_nowrite, ncid = ncid))
557  if (log_get_logging_level() .ge. log_debug) then
558  call log_master_log(log_debug, "Reading in subsidence velocity profile from:"//trim(input_file) )
559  end if
560 
561  call read_2d_forcing_dimensions(ncid, time_dim,z_dim)
562  allocate(forcing_input_times(time_dim), z_subs_pl(z_dim), f_subs_2d(z_dim, time_dim))
563  call read_2d_forcing_variables(trim(input_file), ncid, time_dim, forcing_input_times, &
564  z_dim, z_subs_pl, wsubs_key, f_subs_2d)
565  call check_forcing_status(nf90_close(ncid))
566 
567  zgrid=current_state%global_grid%configuration%vertical%z(:)
568  ! interpolate forcing levels onto the MONC vertical grid (zgrid), for all forcing times
569  allocate(current_state%global_grid%configuration%vertical%wsubs_time_vary(size(zgrid), time_dim))
570  call piecewise_linear_2d(z_subs_pl(1:z_dim), forcing_input_times(1:time_dim), &
571  f_subs_2d(1:z_dim,1:time_dim), zgrid, &
572  current_state%global_grid%configuration%vertical%wsubs_time_vary)
573  else
574  call log_master_log(log_error, "constant forcing from file selected but not coded - STOP")
575  endif
576  endif
577  end if
578 
579  ! Time independent large-scale forcing (proxy for e.g. advection/radiation)
580  ! This probably isn't the right place to be doing this
581  if (.not. allocated(current_state%l_forceq))then
582  allocate(current_state%l_forceq(current_state%number_q_fields))
583  current_state%l_forceq=.false.
584  end if
585 
586  l_constant_forcing_theta=options_get_logical(current_state%options_database, "l_constant_forcing_theta")
587  l_constant_forcing_q=options_get_logical(current_state%options_database, "l_constant_forcing_q")
588  l_constant_forcing_u=options_get_logical(current_state%options_database, "l_constant_forcing_u")
589  l_constant_forcing_v=options_get_logical(current_state%options_database, "l_constant_forcing_v")
590 
591  if (l_constant_forcing_q) then
592  allocate(names_force_pl_q(options_get_array_size(current_state%options_database, "names_constant_forcing_q")))
593  call options_get_string_array(current_state%options_database, "names_constant_forcing_q", names_force_pl_q)
594  end if
595 
596  if (l_constant_forcing_theta)then
597  constant_forcing_type_theta=options_get_integer(current_state%options_database, "constant_forcing_type_theta")
598  forcing_timescale_theta=options_get_real(current_state%options_database, "forcing_timescale_theta")
599  l_constant_forcing_theta_z2pressure=options_get_logical(current_state%options_database, "l_constant_forcing_theta_z2pressure")
600 
601  allocate(z_force_pl_theta(options_get_array_size(current_state%options_database, "z_force_pl_theta")), &
602  f_force_pl_theta(options_get_array_size(current_state%options_database, "f_force_pl_theta")))
603  call options_get_real_array(current_state%options_database, "z_force_pl_theta", z_force_pl_theta)
604  call options_get_real_array(current_state%options_database, "f_force_pl_theta", f_force_pl_theta)
605  ! Get profiles
606  relax_to_initial_theta_profile=options_get_logical(current_state%options_database, "relax_to_initial_theta_profile")
607  if (relax_to_initial_theta_profile)then
608  current_state%global_grid%configuration%vertical%theta_force(:) = &
609  current_state%global_grid%configuration%vertical%theta_init(:)
610  else
611  if (l_constant_forcing_theta_z2pressure)then
612  zgrid=current_state%global_grid%configuration%vertical%zn(:)
613  else
614  zgrid=current_state%global_grid%configuration%vertical%prefn(:)
615  end if
616  call piecewise_linear_1d(z_force_pl_theta(1:size(z_force_pl_theta)), f_force_pl_theta(1:size(f_force_pl_theta)), zgrid, &
617  current_state%global_grid%configuration%vertical%theta_force)
618  end if
619 
620  ! Unit conversions...
621  convert_input_theta_from_temperature=options_get_logical(current_state%options_database, &
622  "convert_input_theta_from_temperature")
623  if (convert_input_theta_from_temperature)then ! Input is temperature not theta
624  current_state%global_grid%configuration%vertical%theta_force(:) = &
625  current_state%global_grid%configuration%vertical%theta_force(:)* &
626  current_state%global_grid%configuration%vertical%prefrcp(:)
627  end if
628 
629  if (constant_forcing_type_theta==tendency)then
630  units_theta_force=options_get_string(current_state%options_database, "units_theta_force")
631  select case(trim(units_theta_force))
632  case(k_per_day)
633  current_state%global_grid%configuration%vertical%theta_force(:) = &
634  current_state%global_grid%configuration%vertical%theta_force(:)/seconds_in_a_day
635  case default !(k_per_second)
636  end select
637  end if
638  deallocate(z_force_pl_theta, f_force_pl_theta)
639  end if
640 
641 #ifdef U_ACTIVE
642  if (l_constant_forcing_u)then
643  constant_forcing_type_u=options_get_integer(current_state%options_database, "constant_forcing_type_u")
644  forcing_timescale_u=options_get_real(current_state%options_database, "forcing_timescale_u")
645  relax_to_initial_u_profile=options_get_logical(current_state%options_database, "relax_to_initial_u_profile")
646  if (relax_to_initial_u_profile)then
647  current_state%global_grid%configuration%vertical%u_force(:) = &
648  current_state%global_grid%configuration%vertical%u_init(:)
649  else
650  allocate(z_force_pl_u(options_get_array_size(current_state%options_database, "z_force_pl_u")), &
651  f_force_pl_u(options_get_array_size(current_state%options_database, "f_force_pl_u")))
652  call options_get_real_array(current_state%options_database, "z_force_pl_u", z_force_pl_u)
653  call options_get_real_array(current_state%options_database, "f_force_pl_u", f_force_pl_u)
654  ! Get profiles
655  zgrid=current_state%global_grid%configuration%vertical%zn(:)
656  call piecewise_linear_1d(z_force_pl_u(1:size(z_force_pl_u)), f_force_pl_u(1:size(f_force_pl_u)), zgrid, &
657  current_state%global_grid%configuration%vertical%u_force)
658  deallocate(z_force_pl_u, f_force_pl_u)
659  end if
660 
661 
662  if (constant_forcing_type_u==tendency)then
663  ! Unit conversions...
664  units_u_force=options_get_string(current_state%options_database, "units_u_force")
665  select case(trim(units_u_force))
666  case(m_per_second_per_day)
667  current_state%global_grid%configuration%vertical%u_force(:) = &
668  current_state%global_grid%configuration%vertical%u_force(:)/seconds_in_a_day
669  case default !(m_per_second_per_second)
670  end select
671  end if
672  end if
673 #endif
674 
675 #ifdef V_ACTIVE
676  if (l_constant_forcing_v)then
677  constant_forcing_type_v=options_get_integer(current_state%options_database, "constant_forcing_type_v")
678  forcing_timescale_v=options_get_real(current_state%options_database, "forcing_timescale_v")
679  relax_to_initial_v_profile=options_get_logical(current_state%options_database, "relax_to_initial_v_profile")
680  if (relax_to_initial_v_profile)then
681  current_state%global_grid%configuration%vertical%v_force(:) = &
682  current_state%global_grid%configuration%vertical%v_init(:)
683  else
684  allocate(z_force_pl_v(options_get_array_size(current_state%options_database, "z_force_pl_v")), &
685  f_force_pl_v(options_get_array_size(current_state%options_database, "f_force_pl_v")))
686  call options_get_real_array(current_state%options_database, "z_force_pl_v", z_force_pl_v)
687  call options_get_real_array(current_state%options_database, "f_force_pl_v", f_force_pl_v)
688  ! Get profiles
689  zgrid=current_state%global_grid%configuration%vertical%zn(:)
690  call piecewise_linear_1d(z_force_pl_v(1:size(z_force_pl_v)), f_force_pl_v(1:size(f_force_pl_v)), zgrid, &
691  current_state%global_grid%configuration%vertical%v_force)
692  deallocate(z_force_pl_v, f_force_pl_v)
693  end if
694 
695 
696  if (constant_forcing_type_v==tendency)then
697  ! Unit conversions...
698  units_v_force=options_get_string(current_state%options_database, "units_v_force")
699  select case(trim(units_v_force))
700  case(m_per_second_per_day)
701  current_state%global_grid%configuration%vertical%v_force(:) = &
702  current_state%global_grid%configuration%vertical%v_force(:)/seconds_in_a_day
703  case default !(m_per_second_per_second)
704  end select
705  end if
706  end if
707 #endif
708 
709  if (l_constant_forcing_q) then
710  constant_forcing_type_q=options_get_integer(current_state%options_database, "constant_forcing_type_q")
711  forcing_timescale_q=options_get_real(current_state%options_database, "forcing_timescale_q")
712  nq_force=size(names_force_pl_q)
713  allocate(z_force_pl_q(options_get_array_size(current_state%options_database, "z_force_pl_q")))
714  call options_get_real_array(current_state%options_database, "z_force_pl_q", z_force_pl_q)
715  nzq=size(z_force_pl_q)
716  zgrid=current_state%global_grid%configuration%vertical%zn(:)
717  allocate(f_force_pl_q_tmp(nq_force*nzq))
718  call options_get_real_array(current_state%options_database, "f_force_pl_q", f_force_pl_q_tmp)
719  allocate(f_force_pl_q(nzq, nq_force))
720  f_force_pl_q(1:nzq, 1:nq_force)=reshape(f_force_pl_q_tmp, (/nzq, nq_force/))
721 
722  allocate(units_q_force(options_get_array_size(current_state%options_database, "units_q_force")))
723  call options_get_string_array(current_state%options_database, "units_q_force", units_q_force)
724  do n=1, nq_force
725  iq=get_q_index(trim(names_force_pl_q(n)), 'forcing:time-independent')
726  call piecewise_linear_1d(z_force_pl_q(1:nzq), f_force_pl_q(1:nzq,n), zgrid, &
727  current_state%global_grid%configuration%vertical%q_force(:,iq))
728 
729  current_state%l_forceq(iq)=.true.
730 
731  ! Unit conversions...
732  if (constant_forcing_type_q==tendency)then
733  select case(trim(units_q_force(n)))
734  case(kg_per_kg_per_day)
735  current_state%global_grid%configuration%vertical%q_force(:,iq) = &
736  current_state%global_grid%configuration%vertical%q_force(:,iq)/seconds_in_a_day
737  case(g_per_kg_per_day)
738  current_state%global_grid%configuration%vertical%q_force(:,iq) = &
739  0.001*current_state%global_grid%configuration%vertical%q_force(:,iq)/seconds_in_a_day
740  case(g_per_kg_per_second)
741  current_state%global_grid%configuration%vertical%q_force(:,iq) = &
742  0.001*current_state%global_grid%configuration%vertical%q_force(:,iq)
743  case default !(kg_per_kg_per_second)
744  end select
745  end if
746  end do
747  deallocate(f_force_pl_q_tmp, units_q_force, f_force_pl_q, z_force_pl_q)
748  end if
749 
750  deallocate(zgrid)
751 
752  ! Set tendency diagnostic logicals based on availability
753  ! Need to use 3d tendencies to compute the profiles, so they will be allocated
754  ! in the case where profiles are available
755  l_qdiag = (.not. current_state%passive_q .and. current_state%number_q_fields .gt. 0)
756 
757  l_tend_pr_tot_u = current_state%u%active
758  l_tend_pr_tot_v = current_state%v%active
759  l_tend_pr_tot_th = current_state%th%active
760  l_tend_pr_tot_qv = l_qdiag .and. current_state%number_q_fields .ge. 1
761  l_tend_pr_tot_ql = l_qdiag .and. current_state%number_q_fields .ge. 2
762  l_tend_pr_tot_qi = l_qdiag .and. current_state%number_q_fields .ge. 11
763  l_tend_pr_tot_qr = l_qdiag .and. current_state%number_q_fields .ge. 11
764  l_tend_pr_tot_qs = l_qdiag .and. current_state%number_q_fields .ge. 11
765  l_tend_pr_tot_qg = l_qdiag .and. current_state%number_q_fields .ge. 11
766  l_tend_pr_tot_tabs = l_tend_pr_tot_th
767 
768  l_tend_3d_u = current_state%u%active .or. l_tend_pr_tot_u
769  l_tend_3d_v = current_state%v%active .or. l_tend_pr_tot_v
770  l_tend_3d_th = current_state%th%active .or. l_tend_pr_tot_th
771  l_tend_3d_qv = (l_qdiag .and. current_state%number_q_fields .ge. 1) .or. l_tend_pr_tot_qv
772  l_tend_3d_ql = (l_qdiag .and. current_state%number_q_fields .ge. 2) .or. l_tend_pr_tot_ql
773  l_tend_3d_qi = (l_qdiag .and. current_state%number_q_fields .ge. 11) .or. l_tend_pr_tot_qi
774  l_tend_3d_qr = (l_qdiag .and. current_state%number_q_fields .ge. 11) .or. l_tend_pr_tot_qr
775  l_tend_3d_qs = (l_qdiag .and. current_state%number_q_fields .ge. 11) .or. l_tend_pr_tot_qs
776  l_tend_3d_qg = (l_qdiag .and. current_state%number_q_fields .ge. 11) .or. l_tend_pr_tot_qg
777  l_tend_3d_tabs = l_tend_3d_th
778 
779  ! Allocate 3d tendency fields upon availability
780  if (l_tend_3d_u) then
781  allocate( tend_3d_u(current_state%local_grid%size(z_index), &
782  current_state%local_grid%size(y_index), &
783  current_state%local_grid%size(x_index) ) )
784  endif
785  if (l_tend_3d_v) then
786  allocate( tend_3d_v(current_state%local_grid%size(z_index), &
787  current_state%local_grid%size(y_index), &
788  current_state%local_grid%size(x_index) ) )
789  endif
790  if (l_tend_3d_th) then
791  allocate( tend_3d_th(current_state%local_grid%size(z_index), &
792  current_state%local_grid%size(y_index), &
793  current_state%local_grid%size(x_index) ) )
794  endif
795  if (l_tend_3d_qv) then
796  allocate( tend_3d_qv(current_state%local_grid%size(z_index), &
797  current_state%local_grid%size(y_index), &
798  current_state%local_grid%size(x_index) ) )
799  endif
800  if (l_tend_3d_ql) then
801  allocate( tend_3d_ql(current_state%local_grid%size(z_index), &
802  current_state%local_grid%size(y_index), &
803  current_state%local_grid%size(x_index) ) )
804  endif
805  if (l_tend_3d_qi) then
806  allocate( tend_3d_qi(current_state%local_grid%size(z_index), &
807  current_state%local_grid%size(y_index), &
808  current_state%local_grid%size(x_index) ) )
809  endif
810  if (l_tend_3d_qr) then
811  allocate( tend_3d_qr(current_state%local_grid%size(z_index), &
812  current_state%local_grid%size(y_index), &
813  current_state%local_grid%size(x_index) ) )
814  endif
815  if (l_tend_3d_qs) then
816  allocate( tend_3d_qs(current_state%local_grid%size(z_index), &
817  current_state%local_grid%size(y_index), &
818  current_state%local_grid%size(x_index) ) )
819  endif
820  if (l_tend_3d_qg) then
821  allocate( tend_3d_qg(current_state%local_grid%size(z_index), &
822  current_state%local_grid%size(y_index), &
823  current_state%local_grid%size(x_index) ) )
824  endif
825  if (l_tend_3d_tabs) then
826  allocate( tend_3d_tabs(current_state%local_grid%size(z_index), &
827  current_state%local_grid%size(y_index), &
828  current_state%local_grid%size(x_index) ) )
829  endif
830 
831  ! Allocate profile tendency fields upon availability
832  if (l_tend_pr_tot_u) then
833  allocate( tend_pr_tot_u(current_state%local_grid%size(z_index)) )
834  endif
835  if (l_tend_pr_tot_v) then
836  allocate( tend_pr_tot_v(current_state%local_grid%size(z_index)) )
837  endif
838  if (l_tend_pr_tot_th) then
839  allocate( tend_pr_tot_th(current_state%local_grid%size(z_index)) )
840  endif
841  if (l_tend_pr_tot_qv) then
842  allocate( tend_pr_tot_qv(current_state%local_grid%size(z_index)) )
843  endif
844  if (l_tend_pr_tot_ql) then
845  allocate( tend_pr_tot_ql(current_state%local_grid%size(z_index)) )
846  endif
847  if (l_tend_pr_tot_qi) then
848  allocate( tend_pr_tot_qi(current_state%local_grid%size(z_index)) )
849  endif
850  if (l_tend_pr_tot_qr) then
851  allocate( tend_pr_tot_qr(current_state%local_grid%size(z_index)) )
852  endif
853  if (l_tend_pr_tot_qs) then
854  allocate( tend_pr_tot_qs(current_state%local_grid%size(z_index)) )
855  endif
856  if (l_tend_pr_tot_qg) then
857  allocate( tend_pr_tot_qg(current_state%local_grid%size(z_index)) )
858  endif
859  if (l_tend_pr_tot_tabs) then
860  allocate( tend_pr_tot_tabs(current_state%local_grid%size(z_index)) )
861  endif
862 
863  ! Save the sampling_frequency to force diagnostic calculation on select time steps
864  diagnostic_generation_frequency=options_get_integer(current_state%options_database, "sampling_frequency")
865 
866 
Here is the call graph for this function:
Here is the caller graph for this function:

◆ read_2d_forcing_dimensions()

subroutine forcing_mod::read_2d_forcing_dimensions ( integer, intent(in)  ncid,
integer, intent(out)  time_dim,
integer, intent(out)  z_dim 
)
private

Reads the dimensions for forcing from the NetCDF file. This routine assumes the forcing uses only time and height.

Parameters
ncidThe NetCDF file id
time_dimNumber of elements in the time dimension

Definition at line 1406 of file forcing.F90.

1407  integer, intent(in) :: ncid
1408  integer, intent(out) :: time_dim
1409  integer, intent(out) :: z_dim
1410  integer :: time_dimid, z_dimid
1411 
1412  call check_forcing_status(nf90_inq_dimid(ncid, time_key, time_dimid))
1413  call check_forcing_status(nf90_inquire_dimension(ncid, time_dimid, len=time_dim))
1414 
1415  call check_forcing_status(nf90_inq_dimid(ncid, z_key, z_dimid))
1416  call check_forcing_status(nf90_inquire_dimension(ncid, z_dimid, len=z_dim))
1417 
Here is the call graph for this function:
Here is the caller graph for this function:

◆ read_2d_forcing_variables()

subroutine forcing_mod::read_2d_forcing_variables ( character(*), intent(in)  filename,
integer, intent(in)  ncid,
integer, intent(in)  time_dim,
real(kind=default_precision), dimension(:), intent(inout)  time,
integer, intent(in)  z_dim,
real(kind=default_precision), dimension(:), intent(inout)  z_profile,
character(len=*), intent(in)  force_2d_key,
real(kind=default_precision), dimension(:,:), intent(inout), allocatable  force_2d_var 
)
private

Reads the variables from the NetCDF forcing file. The 2d variables are assumed to be time and height.

Parameters
ncidThe id of the NetCDF file
time_dimThe number of elements in the time dimension
timeThe time data field that is to be read
z_dimis the number of elements in the height dimension
force_2d_keyis the string that defines the forcing variable in teh NetCDF file
force_2d_varis the forcing data field that is read with dimension (t_dim, z_dim)

Definition at line 1427 of file forcing.F90.

1429 
1430  character(*), intent(in) :: filename
1431  character(len=*), intent(in) :: force_2d_key
1432  integer, intent(in) :: ncid, time_dim, z_dim
1433  real(kind=default_precision), intent(inout) :: time(:), z_profile(:)
1434  real(kind=default_precision), dimension(:,:), allocatable, intent(inout) :: force_2d_var
1435 
1436  integer :: status, variable_id
1437 
1438  ! Do some checking on the variable contents so that we can deal with different
1439  ! variable names or missing variables
1440 
1441  ! time and height
1442  status=nf90_inq_varid(ncid, time_key, variable_id)
1443  if (status==nf90_noerr)then
1444  call read_single_forcing_variable(ncid, time_key, data1d=time)
1445  else
1446  call log_log(log_error, "No recognized time variable found in"//trim(filename))
1447  end if
1448 
1449  status=nf90_inq_varid(ncid, z_key, variable_id)
1450  if (status==nf90_noerr)then
1451  call read_single_forcing_variable(ncid, z_key, data1d=z_profile)
1452  else
1453  call log_log(log_error, "No recognized height variable found in"//trim(filename))
1454  end if
1455 
1456  status=nf90_inq_varid(ncid, force_2d_key, variable_id)
1457  if (status==nf90_noerr)then
1458  call read_single_forcing_variable(ncid, force_2d_key, data2d=force_2d_var)
1459  else
1460  call log_log(log_error, "No recognized forcing variable found in"//trim(filename))
1461  end if
1462 
Here is the call graph for this function:
Here is the caller graph for this function:

◆ read_single_forcing_variable()

subroutine forcing_mod::read_single_forcing_variable ( integer, intent(in)  ncid,
character(len=*), intent(in)  key,
real(kind=default_precision), dimension(:), intent(inout), optional  data1d,
real(kind=default_precision), dimension(:,:), intent(inout), optional  data2d 
)
private

Reads a single variable out of a NetCDF file.

Parameters
ncidThe NetCDF file id
keyThe variable key (name) to access
data1dOptional one dimensional data to read into
data3dOptional three dimensional data to read into

Definition at line 1470 of file forcing.F90.

1471  integer, intent(in) :: ncid
1472  character(len=*), intent(in) :: key
1473  real(kind=default_precision), dimension(:), intent(inout), optional :: data1d
1474  real(kind=default_precision), dimension(:,:), intent(inout), optional :: data2d
1475 
1476  integer :: variable_id
1477  real(kind=default_precision), dimension(:,:), allocatable :: sdata
1478 
1479  call check_forcing_status(nf90_inq_varid(ncid, key, variable_id))
1480 
1481  if (.not. present(data1d) .and. .not. present(data2d)) return
1482 
1483  if (present(data1d)) then
1484  call check_forcing_status(nf90_get_var(ncid, variable_id, data1d))
1485  else
1486  ! 2D
1487  allocate(sdata(size(data2d,1),size(data2d,2)))
1488  call check_forcing_status(nf90_get_var(ncid, variable_id, sdata))
1489  data2d(:,:)=reshape(sdata(:,:),(/size(data2d,1),size(data2d,2)/))
1490  !deallocate(sdata)
1491  end if
Here is the call graph for this function:
Here is the caller graph for this function:

◆ save_precomponent_tendencies()

subroutine forcing_mod::save_precomponent_tendencies ( type(model_state_type), intent(in), target  current_state,
integer, intent(in)  cxn,
integer, intent(in)  cyn,
integer, intent(in)  txn,
integer, intent(in)  tyn 
)
private

Save the 3d tendencies coming into the component.

Parameters
current_stateCurrent model state
cxnThe current slice, x, index
cynThe current column, y, index.
txntarget_x_index
tyntarget_y_index

Definition at line 1255 of file forcing.F90.

1256  type(model_state_type), target, intent(in) :: current_state
1257  integer, intent(in) :: cxn, cyn, txn, tyn
1258 
1259  ! Save 3d tendency fields upon request (of 3d or profiles) and availability
1260  if (l_tend_3d_u) then
1261  tend_3d_u(:,tyn,txn)=current_state%su%data(:,cyn,cxn)
1262  endif
1263  if (l_tend_3d_v) then
1264  tend_3d_v(:,tyn,txn)=current_state%sv%data(:,cyn,cxn)
1265  endif
1266  if (l_tend_3d_th) then
1267  tend_3d_th(:,tyn,txn)=current_state%sth%data(:,cyn,cxn)
1268  endif
1269  if (l_tend_3d_qv) then
1270  tend_3d_qv(:,tyn,txn)=current_state%sq(iqv)%data(:,cyn,cxn)
1271  endif
1272  if (l_tend_3d_ql) then
1273  tend_3d_ql(:,tyn,txn)=current_state%sq(iql)%data(:,cyn,cxn)
1274  endif
1275  if (l_tend_3d_qi) then
1276  tend_3d_qi(:,tyn,txn)=current_state%sq(iqi)%data(:,cyn,cxn)
1277  endif
1278  if (l_tend_3d_qr) then
1279  tend_3d_qr(:,tyn,txn)=current_state%sq(iqr)%data(:,cyn,cxn)
1280  endif
1281  if (l_tend_3d_qs) then
1282  tend_3d_qs(:,tyn,txn)=current_state%sq(iqs)%data(:,cyn,cxn)
1283  endif
1284  if (l_tend_3d_qg) then
1285  tend_3d_qg(:,tyn,txn)=current_state%sq(iqg)%data(:,cyn,cxn)
1286  endif
1287  if (l_tend_3d_tabs) then
1288  tend_3d_tabs(:,tyn,txn)=current_state%sth%data(:,cyn,cxn) * current_state%global_grid%configuration%vertical%rprefrcp(:)
1289  endif
1290 
Here is the caller graph for this function:

◆ set_published_field_value()

subroutine forcing_mod::set_published_field_value ( type(component_field_value_type), intent(inout)  field_value,
real(kind=default_precision), dimension(:), optional  real_1d_field,
real(kind=default_precision), dimension(:,:), optional  real_2d_field,
real(kind=default_precision), dimension(:,:,:), optional  real_3d_field 
)
private

Sets the published field value from the temporary diagnostic values held by this component.

Parameters
field_valuePopulated with the value of the field
real_1d_fieldOptional one dimensional real of values to publish
real_2d_fieldOptional two dimensional real of values to publish

Definition at line 1376 of file forcing.F90.

1377  type(component_field_value_type), intent(inout) :: field_value
1378  real(kind=default_precision), dimension(:), optional :: real_1d_field
1379  real(kind=default_precision), dimension(:,:), optional :: real_2d_field
1380  real(kind=default_precision), dimension(:,:,:), optional :: real_3d_field
1381 
1382  if (present(real_1d_field)) then
1383  allocate(field_value%real_1d_array(size(real_1d_field)), source=real_1d_field)
1384  else if (present(real_2d_field)) then
1385  allocate(field_value%real_2d_array(size(real_2d_field, 1), size(real_2d_field, 2)), source=real_2d_field)
1386  else if (present(real_3d_field)) then
1387  allocate(field_value%real_3d_array(size(real_3d_field, 1), size(real_3d_field, 2), size(real_3d_field, 3)), &
1388  source=real_3d_field)
1389  end if
Here is the caller graph for this function:

◆ timestep_callback()

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

Called for each data column and will determine the forcing values in x and y which are then applied to the field source terms.

Parameters
current_stateThe current model state_mod

Definition at line 872 of file forcing.F90.

873  type(model_state_type), target, intent(inout) :: current_state
874  integer :: current_x_index, current_y_index, target_x_index, target_y_index
875 
876  current_x_index=current_state%column_local_x
877  current_y_index=current_state%column_local_y
878  target_y_index=current_y_index-current_state%local_grid%halo_size(y_index)
879  target_x_index=current_x_index-current_state%local_grid%halo_size(x_index)
880 
881  if (current_state%first_timestep_column) then
882  du_profile_diag=0.0_default_precision
883  dv_profile_diag=0.0_default_precision
884  dtheta_profile_diag=0.0_default_precision
885  dq_profile_diag=0.0_default_precision
886  du_subs_profile_diag=0.0_default_precision
887  dv_subs_profile_diag=0.0_default_precision
888  dtheta_subs_profile_diag=0.0_default_precision
889  dq_subs_profile_diag=0.0_default_precision
890  if (l_tend_pr_tot_u) then
891  tend_pr_tot_u(:)= 0.0_default_precision
892  endif
893  if (l_tend_pr_tot_v) then
894  tend_pr_tot_v(:)= 0.0_default_precision
895  endif
896  if (l_tend_pr_tot_th) then
897  tend_pr_tot_th(:)=0.0_default_precision
898  endif
899  if (l_tend_pr_tot_qv) then
900  tend_pr_tot_qv(:)=0.0_default_precision
901  endif
902  if (l_tend_pr_tot_ql) then
903  tend_pr_tot_ql(:)=0.0_default_precision
904  endif
905  if (l_tend_pr_tot_qi) then
906  tend_pr_tot_qi(:)=0.0_default_precision
907  endif
908  if (l_tend_pr_tot_qr) then
909  tend_pr_tot_qr(:)=0.0_default_precision
910  endif
911  if (l_tend_pr_tot_qs) then
912  tend_pr_tot_qs(:)=0.0_default_precision
913  endif
914  if (l_tend_pr_tot_qg) then
915  tend_pr_tot_qg(:)=0.0_default_precision
916  end if
917  if (l_tend_pr_tot_tabs) then
918  tend_pr_tot_tabs(:)=0.0_default_precision
919  endif
920  end if
921 
922  if (current_state%halo_column .or. current_state%timestep<3) return
923 
924  if (mod(current_state%timestep, diagnostic_generation_frequency) == 0) then
925  call save_precomponent_tendencies(current_state, current_x_index, current_y_index, target_x_index, target_y_index)
926  end if
927 
928  ! AH: perform subsidence calculation but first determine if time varying or constant
929  ! If timevarying then work out the profile of subsidence for the given time and
930  ! assign to w_subs, which is used in apply_subsidence_to...
931  !
932  if (use_time_varying_subsidence) then
933  call interpolate_point_linear_2d(forcing_input_times, &
934  current_state%global_grid%configuration%vertical%wsubs_time_vary, &
935  current_state%time, current_state%global_grid%configuration%vertical%w_subs)
936  endif ! if not w_subs is constant and set in the init_callback
937 
938  if (l_subs_pl_theta) then
939  call apply_subsidence_to_flow_fields(current_state)
940  call apply_subsidence_to_theta(current_state)
941  end if
942  if (l_subs_pl_q) call apply_subsidence_to_q_fields(current_state)
943 
944  if (l_constant_forcing_theta)call apply_time_independent_forcing_to_theta(current_state)
945 #ifdef U_ACTIVE
946  if (l_constant_forcing_u)call apply_time_independent_forcing_to_u(current_state)
947 #endif
948 #ifdef V_ACTIVE
949  if (l_constant_forcing_v)call apply_time_independent_forcing_to_v(current_state)
950 #endif
951  if (l_constant_forcing_q)call apply_time_independent_forcing_to_q(current_state)
952 
953  if (mod(current_state%timestep, diagnostic_generation_frequency) == 0) then
954  call compute_component_tendencies(current_state, current_x_index, current_y_index, target_x_index, target_y_index)
955  end if
956 
Here is the call graph for this function:
Here is the caller graph for this function:

Variable Documentation

◆ constant_forcing_type_q

integer forcing_mod::constant_forcing_type_q =TENDENCY
private

Definition at line 75 of file forcing.F90.

75  integer :: constant_forcing_type_q=tendency ! Method for large-scale forcing of q

◆ constant_forcing_type_theta

integer forcing_mod::constant_forcing_type_theta =TENDENCY
private

Definition at line 74 of file forcing.F90.

74  integer :: constant_forcing_type_theta=tendency ! Method for large-scale forcing of theta

◆ constant_forcing_type_u

integer forcing_mod::constant_forcing_type_u =RELAXATION
private

Definition at line 76 of file forcing.F90.

76  integer :: constant_forcing_type_u=relaxation ! Method for large-scale forcing of u

◆ constant_forcing_type_v

integer forcing_mod::constant_forcing_type_v =RELAXATION
private

Definition at line 77 of file forcing.F90.

77  integer :: constant_forcing_type_v=relaxation ! Method for large-scale forcing of v

◆ diagnostic_generation_frequency

integer forcing_mod::diagnostic_generation_frequency
private

Definition at line 115 of file forcing.F90.

115  integer :: diagnostic_generation_frequency

◆ divergence

integer, parameter forcing_mod::divergence =0
private

Definition at line 37 of file forcing.F90.

37  integer, parameter :: DIVERGENCE=0 ! Input for subsidence forcing is a divergence profile

◆ dq_profile

real(kind=default_precision), dimension(:), allocatable forcing_mod::dq_profile
private

Definition at line 50 of file forcing.F90.

50  real(kind=default_precision), allocatable :: dq_profile(:) ! Local profile to be used in the time-indpendent forcing

◆ dq_profile_diag

real(kind=default_precision), dimension(:,:), allocatable forcing_mod::dq_profile_diag
private

Definition at line 54 of file forcing.F90.

◆ dq_subs_profile_diag

real(kind=default_precision), dimension(:,:), allocatable forcing_mod::dq_subs_profile_diag
private

Definition at line 57 of file forcing.F90.

◆ dtheta_profile

real(kind=default_precision), dimension(:), allocatable forcing_mod::dtheta_profile
private

Definition at line 49 of file forcing.F90.

49  real(kind=default_precision), allocatable :: dtheta_profile(:) ! Local profile to be used in time-indpendent forcing

◆ dtheta_profile_diag

real(kind=default_precision), dimension(:), allocatable forcing_mod::dtheta_profile_diag
private

Definition at line 54 of file forcing.F90.

◆ dtheta_subs_profile_diag

real(kind=default_precision), dimension(:), allocatable forcing_mod::dtheta_subs_profile_diag
private

Definition at line 57 of file forcing.F90.

◆ du_profile

real(kind=default_precision), dimension(:), allocatable forcing_mod::du_profile
private

Definition at line 51 of file forcing.F90.

51  real(kind=default_precision), allocatable :: du_profile(:) ! Local profile to be used in the time-indpendent forcing

◆ du_profile_diag

real(kind=default_precision), dimension(:), allocatable forcing_mod::du_profile_diag
private

Definition at line 54 of file forcing.F90.

54  real(kind=default_precision), allocatable :: du_profile_diag(:), dv_profile_diag(:), dtheta_profile_diag(:), &
55  dq_profile_diag(:,:)

◆ du_subs_profile_diag

real(kind=default_precision), dimension(:), allocatable forcing_mod::du_subs_profile_diag
private

Definition at line 57 of file forcing.F90.

57  real(kind=default_precision), allocatable :: du_subs_profile_diag(:), dv_subs_profile_diag(:), &
58  dtheta_subs_profile_diag(:), dq_subs_profile_diag(:,:)

◆ dv_profile

real(kind=default_precision), dimension(:), allocatable forcing_mod::dv_profile
private

Definition at line 52 of file forcing.F90.

52  real(kind=default_precision), allocatable :: dv_profile(:) ! Local profile to be used in the time-indpendent forcing

◆ dv_profile_diag

real(kind=default_precision), dimension(:), allocatable forcing_mod::dv_profile_diag
private

Definition at line 54 of file forcing.F90.

◆ dv_subs_profile_diag

real(kind=default_precision), dimension(:), allocatable forcing_mod::dv_subs_profile_diag
private

Definition at line 57 of file forcing.F90.

◆ forcing_input_times

real(kind=default_precision), dimension(:), allocatable forcing_mod::forcing_input_times
private

Definition at line 62 of file forcing.F90.

62  real(kind=default_precision), allocatable :: forcing_input_times(:)

◆ forcing_timescale_q

real(kind=default_precision) forcing_mod::forcing_timescale_q
private

Definition at line 65 of file forcing.F90.

65  real(kind=default_precision) :: forcing_timescale_q ! Timescale for forcing of q

◆ forcing_timescale_theta

real(kind=default_precision) forcing_mod::forcing_timescale_theta
private

Definition at line 64 of file forcing.F90.

64  real(kind=default_precision) :: forcing_timescale_theta ! Timescale for forcing of theta

◆ forcing_timescale_u

real(kind=default_precision) forcing_mod::forcing_timescale_u
private

Definition at line 66 of file forcing.F90.

66  real(kind=default_precision) :: forcing_timescale_u ! Timescale for forcing of u

◆ forcing_timescale_v

real(kind=default_precision) forcing_mod::forcing_timescale_v
private

Definition at line 67 of file forcing.F90.

67  real(kind=default_precision) :: forcing_timescale_v ! Timescale for forcing of v

◆ increment

integer, parameter forcing_mod::increment =2
private

Definition at line 42 of file forcing.F90.

42  integer, parameter :: INCREMENT=2 ! Input for large-scale forcing: values are increments (deltas) over timescale

◆ input_file

character(max_file_len) forcing_mod::input_file
private

Definition at line 35 of file forcing.F90.

35  character(MAX_FILE_LEN) :: input_file

◆ iqg

integer forcing_mod::iqg =0
private

Definition at line 113 of file forcing.F90.

◆ iqi

integer forcing_mod::iqi =0
private

Definition at line 113 of file forcing.F90.

◆ iql

integer forcing_mod::iql =0
private

Definition at line 113 of file forcing.F90.

◆ iqr

integer forcing_mod::iqr =0
private

Definition at line 113 of file forcing.F90.

◆ iqs

integer forcing_mod::iqs =0
private

Definition at line 113 of file forcing.F90.

◆ iqv

integer forcing_mod::iqv =0
private

Definition at line 113 of file forcing.F90.

113  integer :: iqv=0, iql=0, iqr=0, iqi=0, iqs=0, iqg=0

◆ l_constant_forcing_q

logical forcing_mod::l_constant_forcing_q
private

Definition at line 70 of file forcing.F90.

70  logical :: l_constant_forcing_q ! Use a time-independent forcing for q

◆ l_constant_forcing_theta

logical forcing_mod::l_constant_forcing_theta
private

Definition at line 69 of file forcing.F90.

69  logical :: l_constant_forcing_theta ! Use a time-independent forcing for theta

◆ l_constant_forcing_theta_z2pressure

logical forcing_mod::l_constant_forcing_theta_z2pressure
private

Definition at line 79 of file forcing.F90.

79  logical :: l_constant_forcing_theta_z2pressure ! profile is a function of pressure not height

◆ l_constant_forcing_u

logical forcing_mod::l_constant_forcing_u
private

Definition at line 71 of file forcing.F90.

71  logical :: l_constant_forcing_u ! Use a time-independent forcing for u

◆ l_constant_forcing_v

logical forcing_mod::l_constant_forcing_v
private

Definition at line 72 of file forcing.F90.

72  logical :: l_constant_forcing_v ! Use a time-independent forcing for v

◆ l_subs_local_q

logical forcing_mod::l_subs_local_q
private

Definition at line 91 of file forcing.F90.

91  logical :: l_subs_local_q ! if .true. then subsidence applied locally (i.e. not with mean fields) to q fields

◆ l_subs_local_theta

logical forcing_mod::l_subs_local_theta
private

Definition at line 90 of file forcing.F90.

90  logical :: l_subs_local_theta ! if .true. then subsidence applied locally (i.e. not with mean fields) to theta field

◆ l_subs_pl_q

logical forcing_mod::l_subs_pl_q
private

Definition at line 88 of file forcing.F90.

88  logical :: l_subs_pl_q ! if .true. then subsidence applied to q fields

◆ l_subs_pl_theta

logical forcing_mod::l_subs_pl_theta
private

Definition at line 87 of file forcing.F90.

87  logical :: l_subs_pl_theta ! if .true. then subsidence applied to theta field

◆ l_tend_3d_qg

logical forcing_mod::l_tend_3d_qg
private

Definition at line 101 of file forcing.F90.

◆ l_tend_3d_qi

logical forcing_mod::l_tend_3d_qi
private

Definition at line 101 of file forcing.F90.

◆ l_tend_3d_ql

logical forcing_mod::l_tend_3d_ql
private

Definition at line 101 of file forcing.F90.

◆ l_tend_3d_qr

logical forcing_mod::l_tend_3d_qr
private

Definition at line 101 of file forcing.F90.

◆ l_tend_3d_qs

logical forcing_mod::l_tend_3d_qs
private

Definition at line 101 of file forcing.F90.

◆ l_tend_3d_qv

logical forcing_mod::l_tend_3d_qv
private

Definition at line 101 of file forcing.F90.

◆ l_tend_3d_tabs

logical forcing_mod::l_tend_3d_tabs
private

Definition at line 101 of file forcing.F90.

◆ l_tend_3d_th

logical forcing_mod::l_tend_3d_th
private

Definition at line 101 of file forcing.F90.

◆ l_tend_3d_u

logical forcing_mod::l_tend_3d_u
private

Definition at line 101 of file forcing.F90.

101  logical :: l_tend_3d_u, l_tend_3d_v, l_tend_3d_th,l_tend_3d_qv, &
102  l_tend_3d_ql,l_tend_3d_qi,l_tend_3d_qr,l_tend_3d_qs,l_tend_3d_qg, &
103  l_tend_3d_tabs

◆ l_tend_3d_v

logical forcing_mod::l_tend_3d_v
private

Definition at line 101 of file forcing.F90.

◆ l_tend_pr_tot_qg

logical forcing_mod::l_tend_pr_tot_qg
private

Definition at line 109 of file forcing.F90.

◆ l_tend_pr_tot_qi

logical forcing_mod::l_tend_pr_tot_qi
private

Definition at line 109 of file forcing.F90.

◆ l_tend_pr_tot_ql

logical forcing_mod::l_tend_pr_tot_ql
private

Definition at line 109 of file forcing.F90.

◆ l_tend_pr_tot_qr

logical forcing_mod::l_tend_pr_tot_qr
private

Definition at line 109 of file forcing.F90.

◆ l_tend_pr_tot_qs

logical forcing_mod::l_tend_pr_tot_qs
private

Definition at line 109 of file forcing.F90.

◆ l_tend_pr_tot_qv

logical forcing_mod::l_tend_pr_tot_qv
private

Definition at line 109 of file forcing.F90.

◆ l_tend_pr_tot_tabs

logical forcing_mod::l_tend_pr_tot_tabs
private

Definition at line 109 of file forcing.F90.

◆ l_tend_pr_tot_th

logical forcing_mod::l_tend_pr_tot_th
private

Definition at line 109 of file forcing.F90.

◆ l_tend_pr_tot_u

logical forcing_mod::l_tend_pr_tot_u
private

Definition at line 109 of file forcing.F90.

109  logical :: l_tend_pr_tot_u, l_tend_pr_tot_v,l_tend_pr_tot_th,l_tend_pr_tot_qv, &
110  l_tend_pr_tot_ql,l_tend_pr_tot_qi,l_tend_pr_tot_qr,l_tend_pr_tot_qs,l_tend_pr_tot_qg, &
111  l_tend_pr_tot_tabs

◆ l_tend_pr_tot_v

logical forcing_mod::l_tend_pr_tot_v
private

Definition at line 109 of file forcing.F90.

◆ max_file_len

integer, parameter forcing_mod::max_file_len =200
private

Maximum length of surface condition input filename.

Definition at line 34 of file forcing.F90.

34  integer, parameter :: MAX_FILE_LEN=200

◆ names_force_pl_q

character(len=string_length), dimension(:), allocatable forcing_mod::names_force_pl_q
private

Definition at line 93 of file forcing.F90.

93  character(len=STRING_LENGTH), dimension(:), allocatable :: names_force_pl_q ! names of q variables to force

◆ q_profile

real(kind=default_precision), dimension(:), allocatable forcing_mod::q_profile
private

Definition at line 45 of file forcing.F90.

45  real(kind=default_precision), allocatable :: q_profile(:) ! Local profile to be used in the subsidence calculation

◆ relax_to_initial_theta_profile

logical forcing_mod::relax_to_initial_theta_profile
private

Definition at line 83 of file forcing.F90.

83  logical :: relax_to_initial_theta_profile ! For relaxation, use initial profile as the target

◆ relax_to_initial_u_profile

logical forcing_mod::relax_to_initial_u_profile
private

Definition at line 81 of file forcing.F90.

81  logical :: relax_to_initial_u_profile ! For relaxation, use initial profile as the target

◆ relax_to_initial_v_profile

logical forcing_mod::relax_to_initial_v_profile
private

Definition at line 82 of file forcing.F90.

82  logical :: relax_to_initial_v_profile ! For relaxation, use initial profile as the target

◆ relaxation

integer, parameter forcing_mod::relaxation =1
private

Definition at line 41 of file forcing.F90.

41  integer, parameter :: RELAXATION=1 ! Input for large-scale forcing: values are target values to relax to over timescale

◆ subsidence

integer, parameter forcing_mod::subsidence =1
private

Definition at line 38 of file forcing.F90.

38  integer, parameter :: SUBSIDENCE=1 ! Input for subsidence forcing is the subsidence velocity profile

◆ tend_3d_qg

real(kind=default_precision), dimension(:,:,:), allocatable forcing_mod::tend_3d_qg
private

Definition at line 97 of file forcing.F90.

◆ tend_3d_qi

real(kind=default_precision), dimension(:,:,:), allocatable forcing_mod::tend_3d_qi
private

Definition at line 97 of file forcing.F90.

◆ tend_3d_ql

real(kind=default_precision), dimension(:,:,:), allocatable forcing_mod::tend_3d_ql
private

Definition at line 97 of file forcing.F90.

◆ tend_3d_qr

real(kind=default_precision), dimension(:,:,:), allocatable forcing_mod::tend_3d_qr
private

Definition at line 97 of file forcing.F90.

◆ tend_3d_qs

real(kind=default_precision), dimension(:,:,:), allocatable forcing_mod::tend_3d_qs
private

Definition at line 97 of file forcing.F90.

◆ tend_3d_qv

real(kind=default_precision), dimension(:,:,:), allocatable forcing_mod::tend_3d_qv
private

Definition at line 97 of file forcing.F90.

◆ tend_3d_tabs

real(kind=default_precision), dimension(:,:,:), allocatable forcing_mod::tend_3d_tabs
private

Definition at line 97 of file forcing.F90.

◆ tend_3d_th

real(kind=default_precision), dimension(:,:,:), allocatable forcing_mod::tend_3d_th
private

Definition at line 97 of file forcing.F90.

◆ tend_3d_u

real(kind=default_precision), dimension(:,:,:), allocatable forcing_mod::tend_3d_u
private

Definition at line 97 of file forcing.F90.

97  real(kind=default_precision), dimension(:,:,:), allocatable :: &
98  tend_3d_u, tend_3d_v, tend_3d_th,tend_3d_qv, &
99  tend_3d_ql,tend_3d_qi,tend_3d_qr,tend_3d_qs,tend_3d_qg, &
100  tend_3d_tabs

◆ tend_3d_v

real(kind=default_precision), dimension(:,:,:), allocatable forcing_mod::tend_3d_v
private

Definition at line 97 of file forcing.F90.

◆ tend_pr_tot_qg

real(kind=default_precision), dimension(:), allocatable forcing_mod::tend_pr_tot_qg
private

Definition at line 105 of file forcing.F90.

◆ tend_pr_tot_qi

real(kind=default_precision), dimension(:), allocatable forcing_mod::tend_pr_tot_qi
private

Definition at line 105 of file forcing.F90.

◆ tend_pr_tot_ql

real(kind=default_precision), dimension(:), allocatable forcing_mod::tend_pr_tot_ql
private

Definition at line 105 of file forcing.F90.

◆ tend_pr_tot_qr

real(kind=default_precision), dimension(:), allocatable forcing_mod::tend_pr_tot_qr
private

Definition at line 105 of file forcing.F90.

◆ tend_pr_tot_qs

real(kind=default_precision), dimension(:), allocatable forcing_mod::tend_pr_tot_qs
private

Definition at line 105 of file forcing.F90.

◆ tend_pr_tot_qv

real(kind=default_precision), dimension(:), allocatable forcing_mod::tend_pr_tot_qv
private

Definition at line 105 of file forcing.F90.

◆ tend_pr_tot_tabs

real(kind=default_precision), dimension(:), allocatable forcing_mod::tend_pr_tot_tabs
private

Definition at line 105 of file forcing.F90.

◆ tend_pr_tot_th

real(kind=default_precision), dimension(:), allocatable forcing_mod::tend_pr_tot_th
private

Definition at line 105 of file forcing.F90.

◆ tend_pr_tot_u

real(kind=default_precision), dimension(:), allocatable forcing_mod::tend_pr_tot_u
private

Definition at line 105 of file forcing.F90.

105  real(kind=default_precision), dimension(:), allocatable :: &
106  tend_pr_tot_u, tend_pr_tot_v, tend_pr_tot_th,tend_pr_tot_qv, &
107  tend_pr_tot_ql,tend_pr_tot_qi,tend_pr_tot_qr,tend_pr_tot_qs,tend_pr_tot_qg, &
108  tend_pr_tot_tabs

◆ tend_pr_tot_v

real(kind=default_precision), dimension(:), allocatable forcing_mod::tend_pr_tot_v
private

Definition at line 105 of file forcing.F90.

◆ tendency

integer, parameter forcing_mod::tendency =0
private

Definition at line 40 of file forcing.F90.

40  integer, parameter :: TENDENCY=0 ! Input for large-scale forcing: values are tendencies (time derivatives)

◆ theta_profile

real(kind=default_precision), dimension(:), allocatable forcing_mod::theta_profile
private

Definition at line 44 of file forcing.F90.

44  real(kind=default_precision), allocatable :: theta_profile(:) ! Local profile to be used in the subsidence calculation

◆ time_key

character(len=*), parameter forcing_mod::time_key = "time"
private

NetCDF data time key.

Definition at line 29 of file forcing.F90.

29  character(len=*), parameter :: &
30  TIME_KEY = "time", & !< NetCDF data time key
31  z_key = "z", &
32  wsubs_key = "wsubs"

◆ u_profile

real(kind=default_precision), dimension(:), allocatable forcing_mod::u_profile
private

Definition at line 46 of file forcing.F90.

46  real(kind=default_precision), allocatable :: u_profile(:) ! Local profile to be used in the subsidence calculation

◆ use_time_varying_subsidence

logical forcing_mod::use_time_varying_subsidence
private

Definition at line 85 of file forcing.F90.

85  logical :: use_time_varying_subsidence ! Use time dependent subsidence veocity (read from file)

◆ v_profile

real(kind=default_precision), dimension(:), allocatable forcing_mod::v_profile
private

Definition at line 47 of file forcing.F90.

47  real(kind=default_precision), allocatable :: v_profile(:) ! Local profile to be used in the subsidence calculation

◆ w_subs_varies_with_time

real(kind=default_precision), dimension(:,:), allocatable forcing_mod::w_subs_varies_with_time
private

Definition at line 61 of file forcing.F90.

61  real(kind=default_precision), allocatable :: w_subs_varies_with_time(:,:)

◆ wsubs_key

character(len=*), parameter forcing_mod::wsubs_key = "wsubs"
private

NetCDF data subsidence velocity.

Definition at line 29 of file forcing.F90.

◆ z_key

character(len=*), parameter forcing_mod::z_key = "z"
private

NetCDF data height(z) key.

Definition at line 29 of file forcing.F90.

logging_mod::log_error
integer, parameter, public log_error
Only log ERROR messages.
Definition: logging.F90:11
logging_mod::log_log
subroutine, public log_log(level, message, str)
Logs a message at the specified level. If the level is above the current level then the message is ig...
Definition: logging.F90:75
logging_mod::log_get_logging_level
integer function, public log_get_logging_level()
Retrieves the current logging level.
Definition: logging.F90:122
logging_mod::log_debug
integer, parameter, public log_debug
Log DEBUG, INFO, WARNING and ERROR messages.
Definition: logging.F90:14
logging_mod::log_master_log
subroutine, public log_master_log(level, message)
Will log just from the master process.
Definition: logging.F90:47
datadefn_mod::default_precision
integer, parameter, public default_precision
MPI communication type which we use for the prognostic and calculation data.
Definition: datadefn.F90:17