MONC
lwrad_exponential.F90
Go to the documentation of this file.
1 
12  use grids_mod, only : z_index, y_index, x_index
13  use science_constants_mod, only : cp
15 
16  implicit none
17 
18  ! Indices for cloud
19  integer :: iql
20  ! Index for the top of the domain
21  integer :: k_top, x_local, y_local
22 
23  ! local column longwave flux variables, calculated using prescribed values from
24  ! the global/user_config file
25  real(DEFAULT_PRECISION), allocatable :: lwrad_flux_top(:), lwrad_flux_base(:), lwrad_flux(:)
26  ! local column liquid water content
27  real(DEFAULT_PRECISION), allocatable :: qc_col(:)
28  ! local density factor and raidiation factor used in conversions
29  real(DEFAULT_PRECISION), allocatable :: density_factor(:), radiation_factor(:)
30  ! declare radiative heating rate for longwave. This is declared as a 3D array so
31  ! that the heating rates can be stored for diagnostics or use when the radiation
32  ! timestep is longer than the model timestep
33  real(DEFAULT_PRECISION), allocatable :: sth_lw(:,:,:)
34 
35  ! declare radiative flux variables that are read for global or user_config
37 
38 
40 contains
41 
45  lwrad_exponential_get_descriptor%name="lwrad_exponential"
51 
55  subroutine initialisation_callback(current_state)
56  type(model_state_type), target, intent(inout) :: current_state
57 
58  integer :: kkp, k ! look counter
59  k_top = current_state%local_grid%size(z_index) + current_state%local_grid%halo_size(z_index) * 2
60  x_local = current_state%local_grid%size(x_index) + current_state%local_grid%halo_size(x_index) * 2
61  y_local = current_state%local_grid%size(y_index) + current_state%local_grid%halo_size(x_index) * 2
62 
63  allocate(lwrad_flux_top(k_top))
64  allocate(lwrad_flux_base(k_top))
65  allocate(lwrad_flux(k_top))
66  allocate(qc_col(k_top))
67  allocate(density_factor(k_top))
68  allocate(radiation_factor(k_top))
69  ! NOTE: this may be the wrong declaration as sth_lw may need to declared on the whole domain
70  allocate(sth_lw(k_top,y_local,x_local))
71 
72  iql=get_q_index(standard_q_names%CLOUD_LIQUID_MASS, 'lwrad_exponential')
73 
74  longwave_exp_decay=options_get_real(current_state%options_database, "longwave_exp_decay")
75  cltop_longwave_flux=options_get_real(current_state%options_database, "cltop_longwave_flux")
76  clbase_longwave_flux=options_get_real(current_state%options_database, "clbase_longwave_flux")
77 
78  density_factor(:) = current_state%global_grid%configuration%vertical%rhon(:)* &
79  current_state%global_grid%configuration%vertical%dz(:)
80 
83 
84  end subroutine initialisation_callback
85 
89  subroutine timestep_callback(current_state)
90  type(model_state_type), target, intent(inout) :: current_state
91 
92  integer :: k ! Loop counter
93  integer :: icol, jcol ! Shorthand column indices
94 
95  real(DEFAULT_PRECISION) :: dtm ! Local timestep variable
96 
97  if (current_state%halo_column) return
98 
99  dtm = current_state%dtm*2.0
100  if (current_state%field_stepping == forward_stepping) dtm=current_state%dtm! Should this be revised to scalar_stepping
101 
102  icol=current_state%column_local_x
103  jcol=current_state%column_local_y
104 
105  ! set the column liquid water content
106  if (current_state%field_stepping == forward_stepping) then ! Should this be revised to scalar_stepping
107  qc_col(:) = current_state%q(iql)%data(:, jcol, icol) + current_state%sq(iql)%data(:, jcol, icol)*dtm
108  else
109  qc_col(:)= current_state%zq(iql)%data(:, jcol, icol) + current_state%sq(iql)%data(:, jcol, icol)*dtm
110 
111  end if
112 
113  ! initialise the flux top and base to 0.0
114  lwrad_flux_top(:) = 0.0
115  lwrad_flux_base(:) = 0.0
116 
117  ! First workout cloud top cooling
119 
120  do k = k_top-1, 1, -1
121  if (qc_col(k+1) > 1.e-10) then
122  lwrad_flux_top(k) = lwrad_flux_top(k+1)* &
124  else
126  endif
127  enddo
128 
129  ! Second workout the cloud base warming
130  do k = 2, k_top
131  if (qc_col(k) > 1.e-10) then
132  lwrad_flux_base(k) = lwrad_flux_base(k-1)* &
134  else
136  endif
137  enddo
138 
139  ! Next, sum the fluxes
140  do k = 1, k_top
142  enddo
143 
144  ! workout radiative heating rate, and stored on the processors local grid - is this correct??
145  ! or should the declaration be on the global grid?
146  do k = 2, k_top
147  sth_lw(k, jcol, icol) = -(lwrad_flux(k) - lwrad_flux(k-1))*radiation_factor(k)
148  enddo
149 
150  ! update the current_state sth
151  current_state%sth%data(:, jcol, icol) = current_state%sth%data(:, jcol, icol) + sth_lw(:, jcol, icol)
152 
153  end subroutine timestep_callback
154 
155  subroutine finalisation_callback(current_state)
157  type(model_state_type), target, intent(inout) :: current_state
159 
160  end subroutine finalisation_callback
161 
162 end module lwrad_exponential_mod
163 
164 
165 
166 
167 
168 
real(default_precision), dimension(:), allocatable radiation_factor
integer, parameter, public forward_stepping
Definition: state.F90:15
type(standard_q_names_type), public standard_q_names
Definition: q_indices.F90:59
subroutine initialisation_callback(current_state)
The initialisation callback sets up the prescribed longwave fluxes and the exponential decay factor...
real(default_precision), dimension(:,:,:), allocatable sth_lw
real(default_precision) cltop_longwave_flux
real(kind=default_precision), public cp
subroutine finalisation_callback(current_state)
real(default_precision), dimension(:), allocatable lwrad_flux
integer, parameter, public default_precision
MPI communication type which we use for the prognostic and calculation data.
Definition: datadefn.F90:17
integer, parameter, public z_index
Grid index parameters.
Definition: grids.F90:14
Contains common definitions for the data and datatypes used by MONC.
Definition: datadefn.F90:2
The ModelState which represents the current state of a run.
Definition: state.F90:39
Scientific constant values used throughout simulations. Each has a default value and this can be over...
This manages the Q variables and specifically the mapping between names and the index that they are s...
Definition: q_indices.F90:2
Simple exponential scheme to calculate the longwave radiation associated with cloud. The scheme is based on the methods used in GASS intercomparison cases, e.g. DYCOMS, ISDAC.
Interfaces and types that MONC components must specify.
real(default_precision) longwave_exp_decay
real(kind=default_precision) function, public options_get_real(options_database, key, index)
Retrieves a real value from the database that matches the provided key.
Functionality to support the different types of grid and abstraction between global grids and local o...
Definition: grids.F90:5
Manages the options database. Contains administration functions and deduce runtime options from the c...
type(component_descriptor_type) function, public lwrad_exponential_get_descriptor()
Provides the descriptor back to the caller and is used in component registration. ...
real(default_precision), dimension(:), allocatable lwrad_flux_top
The model state which represents the current state of a run.
Definition: state.F90:2
real(default_precision), dimension(:), allocatable qc_col
integer, parameter, public y_index
Definition: grids.F90:14
subroutine timestep_callback(current_state)
Called for each column per timestep this will apply a forcing term to the aerosol fields...
real(default_precision), dimension(:), allocatable lwrad_flux_base
real(default_precision) clbase_longwave_flux
integer, parameter, public x_index
Definition: grids.F90:14
integer function, public get_q_index(name, assigning_component)
Add in a new entry into the register if the name does not already exist or return the index of the pr...
Definition: q_indices.F90:112
real(default_precision), dimension(:), allocatable density_factor