plume_mod Module

Provides a concrete implementation of the basal_surface data type, representing a buoyant plume beneath an ice shelf.


Uses

  • module~~plume_mod~~UsesGraph module~plume_mod plume_mod logger_mod logger_mod module~plume_mod->logger_mod h5lt h5lt module~plume_mod->h5lt module~melt_relationship_mod melt_relationship_mod module~plume_mod->module~melt_relationship_mod module~equation_of_state_mod equation_of_state_mod module~plume_mod->module~equation_of_state_mod module~ode_solvers_mod ode_solvers_mod module~plume_mod->module~ode_solvers_mod factual_mod factual_mod module~plume_mod->factual_mod module~upstream_plume_mod upstream_plume_mod module~plume_mod->module~upstream_plume_mod module~jenkins1991_entrainment_mod jenkins1991_entrainment_mod module~plume_mod->module~jenkins1991_entrainment_mod hdf5 hdf5 module~plume_mod->hdf5 module~basal_surface_mod basal_surface_mod module~plume_mod->module~basal_surface_mod module~simple_plume_boundary_mod simple_plume_boundary_mod module~plume_mod->module~simple_plume_boundary_mod module~plume_boundary_mod plume_boundary_mod module~plume_mod->module~plume_boundary_mod iso_fortran_env iso_fortran_env module~plume_mod->iso_fortran_env module~pseudospectral_block_mod pseudospectral_block_mod module~plume_mod->module~pseudospectral_block_mod module~ambient_mod ambient_mod module~plume_mod->module~ambient_mod module~linear_eos_mod linear_eos_mod module~plume_mod->module~linear_eos_mod penf penf module~plume_mod->penf module~uniform_ambient_mod uniform_ambient_mod module~plume_mod->module~uniform_ambient_mod module~entrainment_mod entrainment_mod module~plume_mod->module~entrainment_mod module~boundary_types_mod boundary_types_mod module~plume_mod->module~boundary_types_mod module~dallaston2015_melt_mod dallaston2015_melt_mod module~plume_mod->module~dallaston2015_melt_mod module~melt_relationship_mod->factual_mod module~melt_relationship_mod->iso_fortran_env module~equation_of_state_mod->factual_mod module~equation_of_state_mod->iso_fortran_env module~ode_solvers_mod->logger_mod module~ode_solvers_mod->iso_fortran_env module~ode_solvers_mod->penf module~nitsol_mod nitsol_mod module~ode_solvers_mod->module~nitsol_mod module~upstream_plume_mod->logger_mod module~upstream_plume_mod->factual_mod module~upstream_plume_mod->module~plume_boundary_mod module~upstream_plume_mod->iso_fortran_env module~upstream_plume_mod->penf module~upstream_plume_mod->module~boundary_types_mod module~rksuite_90 rksuite_90 module~upstream_plume_mod->module~rksuite_90 module~uniform_gradient_field_mod uniform_gradient_field_mod module~upstream_plume_mod->module~uniform_gradient_field_mod module~jenkins1991_entrainment_mod->factual_mod module~jenkins1991_entrainment_mod->iso_fortran_env module~jenkins1991_entrainment_mod->module~entrainment_mod module~basal_surface_mod->factual_mod module~basal_surface_mod->hdf5 module~basal_surface_mod->iso_fortran_env module~basal_surface_mod->module~nitsol_mod module~simple_plume_boundary_mod->factual_mod module~simple_plume_boundary_mod->module~plume_boundary_mod module~simple_plume_boundary_mod->iso_fortran_env module~simple_plume_boundary_mod->module~boundary_types_mod module~plume_boundary_mod->factual_mod module~plume_boundary_mod->iso_fortran_env module~plume_boundary_mod->module~boundary_types_mod module~pseudospectral_block_mod->logger_mod module~pseudospectral_block_mod->factual_mod module~pseudospectral_block_mod->iso_fortran_env module~pseudospectral_block_mod->penf chebyshev_mod chebyshev_mod module~pseudospectral_block_mod->chebyshev_mod module~ambient_mod->factual_mod module~ambient_mod->iso_fortran_env module~linear_eos_mod->module~equation_of_state_mod module~linear_eos_mod->factual_mod module~linear_eos_mod->iso_fortran_env module~uniform_ambient_mod->factual_mod module~uniform_ambient_mod->iso_fortran_env module~uniform_ambient_mod->module~ambient_mod module~entrainment_mod->factual_mod module~entrainment_mod->iso_fortran_env module~dallaston2015_melt_mod->module~melt_relationship_mod module~dallaston2015_melt_mod->factual_mod module~dallaston2015_melt_mod->iso_fortran_env module~nitsol_mod->iso_fortran_env module~rksuite_90_prec rksuite_90_prec module~rksuite_90->module~rksuite_90_prec module~uniform_gradient_field_mod->factual_mod module~uniform_gradient_field_mod->iso_fortran_env utils_mod utils_mod module~uniform_gradient_field_mod->utils_mod module~rksuite_90_prec->iso_fortran_env

Contents


Variables

TypeVisibility AttributesNameInitial
character(len=9), public, parameter:: hdf_type_name ='plume'
character(len=9), public, parameter:: hdf_thickness ='thickness'
character(len=8), public, parameter:: hdf_velocity ='velocity'
character(len=11), public, parameter:: hdf_temperature ='temperature'
character(len=8), public, parameter:: hdf_salinity ='salinity'
character(len=5), public, parameter:: hdf_delta ='delta'
character(len=2), public, parameter:: hdf_nu ='nu'
character(len=2), public, parameter:: hdf_mu ='mu'
character(len=5), public, parameter:: hdf_r ='r_val'
character(len=3), public, parameter:: hdf_phi ='phi'

Abstract Interfaces

abstract interface

  • private pure function scalar_func(location) result(scalar)

    Abstract interface for function providing the initial values for the scalar properties of a plume object when it is being instantiated.

    Arguments

    Type IntentOptional AttributesName
    real(kind=r8), intent(in), dimension(:):: location

    The position $\vec{x}$ at which to compute the property

    Return Value real(kind=r8)

    The value of the scalar quantity at location

abstract interface

  • private pure function velocity_func(location) result(vector)

    Abstract interface for function providing the plume velocity when an object is being instantiated.

    Arguments

    Type IntentOptional AttributesName
    real(kind=r8), intent(in), dimension(:):: location

    The position $\vec{x}$ at which to compute the thickness

    Return Value real(kind=r8), dimension(:), allocatable

    The velocity vector of the water in the plume at location


Derived Types

type, public, extends(basal_surface) :: plume

A concrete implementation of the basal_surface abstract data type, representing the buoyant plume beneath an ice shelf.

Components

TypeVisibility AttributesNameInitial
type(cheb1d_scalar_field), private :: thickness

The thickness of the plume

type(cheb1d_vector_field), private :: velocity

The velocity of the plume

type(cheb1d_vector_field), private :: velocity_dx

The derivative of the velocity field

type(cheb1d_scalar_field), private :: temperature

The temperature of the plume

type(cheb1d_scalar_field), private :: temperature_dx

The derivative of the temperature of the plume

type(cheb1d_scalar_field), private :: salinity

The salinity of the plume

type(cheb1d_scalar_field), private :: salinity_dx

The derivative of the salinity of the plume

class(abstract_entrainment), private, allocatable:: entrainment_formulation

An object which provides the parameterisation for entrainment of water into the plume.

class(abstract_melt_relationship), private, allocatable:: melt_formulation

An object which provides the parameterisation for melting, salt, and heat fluxes from the plume to the ice.

class(ambient_conditions), private, allocatable:: ambient_conds

An object specifying the temperature and salinity of the ambient ocean at its interface with the plume.

class(equation_of_state), public, allocatable:: eos

An object specifying the equation of state relating the plume water's density to its temperature and salinity.

class(plume_boundary), private, allocatable:: boundaries

An object specifying the boundary conditions for the plume.

real(kind=r8), private :: delta

The dimensionless ratio

real(kind=r8), public :: nu

The dimensionless ratio

real(kind=r8), private :: mu

The dimensionless ratio

real(kind=r8), private :: r_val

The dimensionless ratio of the ocean water density to the density of the overlying ice shelf.

real(kind=r8), public :: phi

The inverse Rossby number,

real(kind=r8), private :: time

The time at which the ice shelf is in this state

integer, private :: thickness_size

The number of data values in the thickness field

integer, private :: velocity_size

The number of data values in the velocity field

integer, private :: temperature_size

The number of data values in the temperature field

integer, private :: salinity_size

the number of data values in the salinity field

logical, private, dimension(7):: lower_bounds =.false.

Which variables have boundary conditions at the grounding line.

logical, private, dimension(7):: upper_bounds =.false.

Which variables have boundary conditions at the calving front.

type(pseudospec_block), private :: precond

A pseudospectral differentiation block which can be used for preconditioning.

Type-Bound Procedures

procedure, public :: initialise => plume_initialise
procedure, public :: basal_melt => plume_melt
procedure, public :: basal_drag_parameter => plume_drag_parameter
procedure, public :: water_density => plume_water_density
procedure, public :: update => plume_update
procedure, public :: data_size => plume_data_size
procedure, public :: state_vector => plume_state_vector
procedure, public :: read_data => plume_read_data
procedure, public :: write_data => plume_write_data
procedure, public :: solve => plume_solve

Functions

private function plume_melt(this) result(melt)

Author
Christopher MacMackin
Date
April 2016

Computes and returns the melt rate at the bottom of the ice shelf due to interaction with the plume.

Arguments

Type IntentOptional AttributesName
class(plume), intent(in) :: this

Return Value class(scalar_field), pointer

The melt rate at the base of the ice shelf.

private function plume_drag_parameter(this) result(drag)

Author
Christopher MacMackin
Date
April 2016

Computes and returns a quantity which may be necessary to determine the frictional drag the plume exerts on the bottom of the ice shelf. The plume would actually tend to exert no drag on the bottom of the ice shelf, but this method is present so that there is a consistent interface with the ground data type.

Arguments

Type IntentOptional AttributesName
class(plume), intent(in) :: this

Return Value class(scalar_field), pointer

The melt rate at the base of the ice sheet.

private function plume_water_density(this) result(density)

Author
Christopher MacMackin
Date
April 2016

Computes and returns the density of the plume water beneath the ice shelf. The density of this water would vary depending on how much saline ambient water has been entrained into the plume versus how much fresh water has been released due to melting. However, the Boussinesq approximation is used here and only a single reference density is returned.

Read more…

Arguments

Type IntentOptional AttributesName
class(plume), intent(in) :: this

Return Value real(kind=r8)

The density of the water at the base of the ice sheet.

private function plume_data_size(this)

Author
Christopher MacMackin
Date
August 2016

Returns the number of elements in the plume's state vector. This is the size of the vector returned by state_vector and taken as an argument by update.

Arguments

Type IntentOptional AttributesName
class(plume), intent(in) :: this

Return Value integer

The number of elements in the plume's state vector.

private function plume_state_vector(this) result(state_vector)

Author
Christopher MacMackin
Date
April 2016

Returns the state vector for the current state of the plume. This takes the form of a 1D array.

Arguments

Type IntentOptional AttributesName
class(plume), intent(in) :: this

Return Value real(kind=r8), dimension(:), allocatable

The state vector describing the plume.


Subroutines

private subroutine plume_initialise(this, domain, resolution, thickness, velocity, temperature, salinity, entrainment_formulation, melt_formulation, ambient_conds, eos, boundaries, delta, nu, mu, r_val, phi)

Author
Christopher MacMackin
Date
April 2016

Instantiates a plume object with initial coniditions provided by the arguments.At present only a 1D model is supported. If information is provided for higher dimensions then it will be ignored.

Arguments

Type IntentOptional AttributesName
class(plume), intent(out) :: this

A plume object with its domain and initial conditions set according to the arguments of the constructor function.

real(kind=r8), intent(in), dimension(:,:):: domain

An array containing the upper and lower limits of the domain for the plume. The first index represents the dimension for which the boundaries apply. If the second index is 1 then it corresponds to the lower bound. If the second index is 2 then it corresponds to the upper bound.

integer, intent(in), dimension(:):: resolution

The number of data points in each dimension

procedure(scalar_func) :: thickness

A function which calculates the initial value of the thickness of the plume at a given location.

procedure(velocity_func) :: velocity

A function which calculates the initial value of the velocity (vector) of the water at a given location in a plume.

procedure(scalar_func) :: temperature

A function which calculates the initial value of the temperature of the plume at a given location.

procedure(scalar_func) :: salinity

A function which calculates the initial value of the salinity of the plume at a given location.

class(abstract_entrainment), intent(inout), optional allocatable:: entrainment_formulation

An object which calculates entrainment into the plume. Will be unallocated on exit. Defaults to that used by Jenkins (1991) with the coefficient $E_0 = 1$.

class(abstract_melt_relationship), intent(inout), optional allocatable:: melt_formulation

An object which calculates melting and the resulting thermal transfer into/out of the plume. Will be unallocated on exit. Defaults to that used by Dallaston et al. (2015), scaled to be consistent with the nondimensionalisation used here.

class(ambient_conditions), intent(inout), optional allocatable:: ambient_conds

An object specifying the salinity and temperature of the ambient ocean. Will be unallocated on exit. Defaults to uniform ambient salinity and temperature, both of which are set to 0 (as temperature and salinity are measured relative to some reference value).

class(equation_of_state), intent(inout), optional allocatable:: eos

An object specifying the equation of state for the water in the plume. Will be unallocated on exit. Defaults to linearised equation of state with no temperature dependence and a haline contraction coefficient of 1. The reference density is set to be 1 in the dimensionless units when salinity and temeprature are 0.

class(plume_boundary), intent(inout), optional allocatable:: boundaries

An object providing the boundary conditions for the plume. Will be unallocated on exit. Defaults to those used by Dallaston et al. (2015).

real(kind=r8), intent(in), optional :: delta

The dimensionless ratio . Defaults to 0.036.

real(kind=r8), intent(in), optional :: nu

The dimensionless ratio . Defaults to 0.

real(kind=r8), intent(in), optional :: mu

The dimensionless ratio . Defaults to 0.

real(kind=r8), intent(in), optional :: r_val

The dimensionless ratio of the water density to the ice shelf density, Defaults to 1.12.

real(kind=r8), intent(in), optional :: phi

The inverse Rossby number, . Defaults to 0.

private subroutine plume_update(this, state_vector, ice_thickness)

Author
Christopher MacMackin
Date
April 2016

Updates the state of the plume from its state vector. The state vector is a real array containing the value of each of the plume's properties at each of the locations on the grid used in discretization.

Arguments

Type IntentOptional AttributesName
class(plume), intent(inout) :: this
real(kind=r8), intent(in), dimension(:):: state_vector

A real array containing the data describing the state of the plume.

class(scalar_field), intent(in), optional :: ice_thickness

The ice thickness which, if present, will be used to update the calculation of the melt rate.

private subroutine plume_read_data(this, file_id, group_name, error)

Author
Chris MacMackin
Date
April 2017

Reads the state of the plume object from an HDF file in the specified group. This sets the thickness, velocity, temperature, salinity dataset, and parameter values.

Arguments

Type IntentOptional AttributesName
class(plume), intent(inout) :: this
integer(kind=hid_t), intent(in) :: file_id

The identifier for the HDF5 file/group in which this data is meant to be written.

character(len=*), intent(in) :: group_name

The name to give the group in the HDF5 file storing the ice shelf's data.

integer, intent(out) :: error

Flag indicating whether routine ran without error. If no error occurs then has value 0.

private subroutine plume_write_data(this, file_id, group_name, error)

Author
Chris MacMackin
Date
November 2016

Writes the state of the plume object to an HDF file in the specified group. This will consist of a thickness, a velocity, a temperature, and a salinity dataset.

Arguments

Type IntentOptional AttributesName
class(plume), intent(in) :: this
integer(kind=hid_t), intent(in) :: file_id

The identifier for the HDF5 file/group in which this data is meant to be written.

character(len=*), intent(in) :: group_name

The name to give the group in the HDF5 file storing the ice shelf's data.

integer, intent(out) :: error

Flag indicating whether routine ran without error. If no error occurs then has value 0.

private subroutine plume_solve(this, ice_thickness, ice_density, ice_temperature, time, success)

Author
Chris MacMackin
Date
March 2017

Solves the state of the plume for the specified ice properties, at the specified time. This is done using the a quasilinearisation method and a GMRES iterative linear solver.

Arguments

Type IntentOptional AttributesName
class(plume), intent(inout) :: this
class(scalar_field), intent(in) :: ice_thickness

Thickness of the ice above the basal surface

real(kind=r8), intent(in) :: ice_density

The density of the ice above the basal surface, assumed uniform

real(kind=r8), intent(in) :: ice_temperature

The temperature of the ice above the basal surface, assumed uniform

real(kind=r8), intent(in) :: time

The time to which the basal surface should be solved

logical, intent(out) :: success

True if the solver is successful, false otherwise