regrid_from_vertex_to_cell Function

private pure function regrid_from_vertex_to_cell(cell_index, cell_level, nverticesoncell, verticesoncell, kiteforcell, kiteareasonvertex, areacell, vertex_value) result(cell_value)

This function computes the area weighted average (i.e., cell_value) at the specified cell point (i.e., cell_index and cell_level) from the values at its surrounding vertex points (i.e., vertex_value). The formulation used here is adapted and generalized from the atm_compute_solve_diagnostics subroutine in MPAS.

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: cell_index
integer, intent(in) :: cell_level
integer, intent(in) :: nverticesoncell(:)
integer, intent(in) :: verticesoncell(:,:)
integer, intent(in) :: kiteforcell(:,:)
real(kind=rkind), intent(in) :: kiteareasonvertex(:,:)
real(kind=rkind), intent(in) :: areacell(:)
real(kind=rkind), intent(in) :: vertex_value(:,:)

Return Value real(kind=rkind)


Called by

proc~~regrid_from_vertex_to_cell~~CalledByGraph proc~regrid_from_vertex_to_cell regrid_from_vertex_to_cell proc~dyn_mpas_compute_cell_relative_vorticity mpas_dynamical_core_type%dyn_mpas_compute_cell_relative_vorticity proc~dyn_mpas_compute_cell_relative_vorticity->proc~regrid_from_vertex_to_cell

Variables

Type Visibility Attributes Name Initial
integer, private :: i
integer, private :: j
integer, private :: vertex_index

Source Code

    pure function regrid_from_vertex_to_cell(cell_index, cell_level, &
            nverticesoncell, verticesoncell, kiteforcell, kiteareasonvertex, areacell, &
            vertex_value) result(cell_value)
        integer, intent(in) :: cell_index, cell_level
        integer, intent(in) :: nverticesoncell(:), verticesoncell(:, :), kiteforcell(:, :)
        real(rkind), intent(in) :: kiteareasonvertex(:, :), areacell(:)
        real(rkind), intent(in) :: vertex_value(:, :)
        real(rkind) :: cell_value

        integer :: i, j, vertex_index

        cell_value = 0.0_rkind

        do i = 1, nverticesoncell(cell_index)
            j = kiteforcell(i, cell_index)
            vertex_index = verticesoncell(i, cell_index)

            cell_value = cell_value + &
                kiteareasonvertex(j, vertex_index) * vertex_value(cell_level, vertex_index)
        end do

        cell_value = cell_value / areacell(cell_index)
    end function regrid_from_vertex_to_cell