30 type(model_state_type),
target,
intent(inout) :: current_state
32 if (.not.
allocated(current_state%p%data))
then 34 allocate(current_state%p%data(current_state%local_grid%size(z_index) + current_state%local_grid%halo_size(z_index) * 2, &
35 current_state%local_grid%size(y_index) + current_state%local_grid%halo_size(y_index) * 2, &
36 current_state%local_grid%size(x_index) + current_state%local_grid%halo_size(x_index) * 2))
37 current_state%p%data=0.0_default_precision
38 current_state%p%active=.true.
45 type(model_state_type),
target,
intent(inout) :: current_state
47 current_state%p%data(:, current_state%column_local_y, current_state%column_local_x) = 0.0_default_precision
49 if (current_state%halo_column)
return 51 if (current_state%field_stepping == forward_stepping)
then 62 type(model_state_type),
target,
intent(inout) :: current_state
64 current_state%p%active=.false.
65 deallocate(current_state%p%data)
71 type(model_state_type),
target,
intent(inout) :: current_state
75 if (current_state%first_timestep_column) current_state%local_divmax=0.0_default_precision
77 do k=2,current_state%local_grid%size(z_index)
78 if (abs(current_state%p%data(k, current_state%column_local_y, current_state%column_local_x)) .gt. &
79 current_state%local_divmax)
then 80 current_state%local_divmax=abs(current_state%p%data(k,current_state%column_local_y, current_state%column_local_x))
88 type(model_state_type),
target,
intent(inout) :: current_state
90 real(kind=DEFAULT_PRECISION) :: timec
92 timec=1.0_default_precision/current_state%dtm
94 call calculate_p(current_state, current_state%p, current_state%u, current_state%v, current_state%w, timec, &
95 current_state%column_local_y, current_state%column_local_x)
101 type(model_state_type),
target,
intent(inout) :: current_state
103 real(kind=DEFAULT_PRECISION) :: timec
105 timec=1.0_default_precision/(2.0_default_precision*current_state%dtm)
107 call calculate_p(current_state, current_state%p, current_state%zu, current_state%zv, current_state%zw, timec, &
108 current_state%column_local_y, current_state%column_local_x)
120 subroutine calculate_p(current_state, p, u, v, w, timec, y_local, x_local)
121 type(model_state_type),
target,
intent(inout) :: current_state
122 type(prognostic_field_type),
intent(inout) :: u, v, w, p
123 real(kind=DEFAULT_PRECISION),
intent(in) :: timec
124 integer,
intent(in) :: y_local, x_local
128 p%data(:, y_local, x_local)=0.0_default_precision
129 do k=2,current_state%local_grid%size(z_index)
131 p%data(k, y_local, x_local)= p%data(k, y_local, x_local)+ &
132 current_state%global_grid%configuration%horizontal%cx*(u%data(k, y_local, x_local)-u%data(k, y_local, x_local-1))
135 p%data(k, y_local, x_local)= p%data(k, y_local, x_local)+ &
136 current_state%global_grid%configuration%horizontal%cy*(v%data(k, y_local, x_local)-v%data(k, y_local-1, x_local))
139 p%data(k, y_local, x_local)= p%data(k, y_local, x_local)+ &
140 4.0_default_precision*(current_state%global_grid%configuration%vertical%tzc2(k)*w%data(k, y_local, x_local)-&
141 current_state%global_grid%configuration%vertical%tzc1(k)* w%data(k-1, y_local, x_local))
143 p%data(k, y_local, x_local)=p%data(k, y_local, x_local) * timec
subroutine init_callback(current_state)
The initialisation callback will allocate memory for the P field and initialise it.
integer, parameter, public forward_stepping
subroutine handle_forward_stepping(current_state)
Handles the calculation of P when the model is forward stepping.
Contains prognostic field definitions and functions.
A prognostic field which is assumed to be 3D.
integer, parameter, public default_precision
MPI communication type which we use for the prognostic and calculation data.
integer, parameter, public z_index
Grid index parameters.
Contains common definitions for the data and datatypes used by MONC.
The ModelState which represents the current state of a run.
Description of a component.
Calculates the local divergence error.
Interfaces and types that MONC components must specify.
subroutine find_max_divergence_error(current_state)
Finds the maximum divergence error locally and writes this into the local variable.
Functionality to support the different types of grid and abstraction between global grids and local o...
type(component_descriptor_type) function, public diverr_get_descriptor()
Descriptor of this component for registration.
subroutine handle_centred_stepping(current_state)
Handles the calculation of P when the model is centred stepping.
subroutine finalisation_callback(current_state)
Called at finalisation this will deallocate the P field as the model shuts down.
subroutine calculate_p(current_state, p, u, v, w, timec, y_local, x_local)
Calculates P based upon flow fields and the stepping.
The model state which represents the current state of a run.
integer, parameter, public y_index
subroutine timestep_callback(current_state)
Called each timestep this will initialise P and update the local divergence error for the current col...
integer, parameter, public x_index