namespace mom_eos

Overview

Provides subroutines for quantities specific to the equation of state. More…

namespace mom_eos {

// interfaces

interface calculate_density;
interface calculate_density_derivs;
interface calculate_density_second_derivs;
interface calculate_spec_vol;
interface calculate_tfreeze;

// global variables

integer, parameter, public eos_linear = 1;
integer, parameter, public eos_unesco = 2;
integer, parameter, public eos_wright = 3;
integer, parameter, public eos_teos10 = 4;
integer, parameter, public eos_nemo = 5;

// global functions

subroutine, public calculate_specific_vol_derivs(
    T T,
    S S,
    pressure pressure,
    dSV_dT dSV_dT,
    dSV_dS dSV_dS,
    start start,
    npts npts,
    EOS EOS
    );

subroutine, public calculate_compress(
    T T,
    S S,
    pressure pressure,
    rho rho,
    drho_dp drho_dp,
    start start,
    npts npts,
    EOS EOS
    );

subroutine, public int_specific_vol_dp(
    T T,
    S S,
    p_t p_t,
    p_b p_b,
    alpha_ref alpha_ref,
    HI HI,
    EOS EOS,
    dza dza,
    intp_dza intp_dza,
    intx_dza intx_dza,
    inty_dza inty_dza,
    halo_size halo_size,
    bathyP bathyP,
    dP_tiny dP_tiny,
    useMassWghtInterp useMassWghtInterp
    );

subroutine, public int_density_dz(
    T T,
    S S,
    z_t z_t,
    z_b z_b,
    rho_ref rho_ref,
    rho_0 rho_0,
    G_e G_e,
    HII HII,
    HIO HIO,
    EOS EOS,
    dpa dpa,
    intz_dpa intz_dpa,
    intx_dpa intx_dpa,
    inty_dpa inty_dpa,
    bathyT bathyT,
    dz_neglect dz_neglect,
    useMassWghtInterp useMassWghtInterp
    );

logical function, public query_compressible(EOS EOS);
subroutine, public eos_init(param_file param_file, EOS EOS);

subroutine, public eos_manual_init(
    EOS EOS,
    form_of_EOS form_of_EOS,
    form_of_TFreeze form_of_TFreeze,
    EOS_quadrature EOS_quadrature,
    Compressible Compressible,
    Rho_T0_S0 Rho_T0_S0,
    drho_dT drho_dT,
    dRho_dS dRho_dS,
    TFr_S0_P0 TFr_S0_P0,
    dTFr_dS dTFr_dS,
    dTFr_dp dTFr_dp
    );

subroutine, public eos_allocate(EOS EOS);
subroutine, public eos_end(EOS EOS);

subroutine, public eos_use_linear(
    Rho_T0_S0 Rho_T0_S0,
    dRho_dT dRho_dT,
    dRho_dS dRho_dS,
    EOS EOS,
    use_quadrature use_quadrature
    );

subroutine, public int_density_dz_generic(
    T T,
    S S,
    z_t z_t,
    z_b z_b,
    rho_ref rho_ref,
    rho_0 rho_0,
    G_e G_e,
    HII HII,
    HIO HIO,
    EOS EOS,
    dpa dpa,
    intz_dpa intz_dpa,
    intx_dpa intx_dpa,
    inty_dpa inty_dpa,
    bathyT bathyT,
    dz_neglect dz_neglect,
    useMassWghtInterp useMassWghtInterp
    );

subroutine, public int_density_dz_generic_plm(
    T_t T_t,
    T_b T_b,
    S_t S_t,
    S_b S_b,
    z_t z_t,
    z_b z_b,
    rho_ref rho_ref,
    rho_0 rho_0,
    G_e G_e,
    dz_subroundoff dz_subroundoff,
    bathyT bathyT,
    HII HII,
    HIO HIO,
    EOS EOS,
    dpa dpa,
    intz_dpa intz_dpa,
    intx_dpa intx_dpa,
    inty_dpa inty_dpa,
    useMassWghtInterp useMassWghtInterp
    );

subroutine, public find_depth_of_pressure_in_cell(
    T_t T_t,
    T_b T_b,
    S_t S_t,
    S_b S_b,
    z_t z_t,
    z_b z_b,
    P_t P_t,
    P_tgt P_tgt,
    rho_ref rho_ref,
    G_e G_e,
    EOS EOS,
    P_b P_b,
    z_out z_out,
    z_tol z_tol
    );

subroutine, public int_density_dz_generic_ppm(
    T T,
    T_t T_t,
    T_b T_b,
    S S,
    S_t S_t,
    S_b S_b,
    z_t z_t,
    z_b z_b,
    rho_ref rho_ref,
    rho_0 rho_0,
    G_e G_e,
    HII HII,
    HIO HIO,
    EOS EOS,
    dpa dpa,
    intz_dpa intz_dpa,
    intx_dpa intx_dpa,
    inty_dpa inty_dpa
    );

subroutine, public int_spec_vol_dp_generic(
    T T,
    S S,
    p_t p_t,
    p_b p_b,
    alpha_ref alpha_ref,
    HI HI,
    EOS EOS,
    dza dza,
    intp_dza intp_dza,
    intx_dza intx_dza,
    inty_dza inty_dza,
    halo_size halo_size,
    bathyP bathyP,
    dP_neglect dP_neglect,
    useMassWghtInterp useMassWghtInterp
    );

subroutine, public int_spec_vol_dp_generic_plm(
    T_t T_t,
    T_b T_b,
    S_t S_t,
    S_b S_b,
    p_t p_t,
    p_b p_b,
    alpha_ref alpha_ref,
    dP_neglect dP_neglect,
    bathyP bathyP,
    HI HI,
    EOS EOS,
    dza dza,
    intp_dza intp_dza,
    intx_dza intx_dza,
    inty_dza inty_dza,
    useMassWghtInterp useMassWghtInterp
    );

subroutine, public convert_temp_salt_for_teos10(T T, S S, press press, G G, kd kd, mask_z mask_z, EOS EOS);

subroutine, public extract_member_eos(
    EOS EOS,
    form_of_EOS form_of_EOS,
    form_of_TFreeze form_of_TFreeze,
    EOS_quadrature EOS_quadrature,
    Compressible Compressible,
    Rho_T0_S0 Rho_T0_S0,
    drho_dT drho_dT,
    dRho_dS dRho_dS,
    TFr_S0_P0 TFr_S0_P0,
    dTFr_dS dTFr_dS,
    dTFr_dp dTFr_dp
    );

} // namespace mom_eos

Detailed Documentation

Provides subroutines for quantities specific to the equation of state.

The MOM_EOS module is a wrapper for various equations of state (e.g. Linear, Wright, UNESCO) and provides a uniform interface to the rest of the model independent of which equation of state is being used.

Global Variables

integer, parameter, public eos_linear = 1

A named integer specifying an equation of state.

integer, parameter, public eos_unesco = 2

A named integer specifying an equation of state.

integer, parameter, public eos_wright = 3

A named integer specifying an equation of state.

integer, parameter, public eos_teos10 = 4

A named integer specifying an equation of state.

integer, parameter, public eos_nemo = 5

A named integer specifying an equation of state.

Global Functions

subroutine, public calculate_specific_vol_derivs(
    T T,
    S S,
    pressure pressure,
    dSV_dT dSV_dT,
    dSV_dS dSV_dS,
    start start,
    npts npts,
    EOS EOS
    )

Calls the appropriate subroutine to calculate specific volume derivatives for an array.

Parameters:

t

Potential temperature referenced to the surface [degC]

s

Salinity [ppt]

pressure

Pressure [Pa]

dsv_dt

The partial derivative of specific volume with potential temperature [m3 kg-1 degC-1].

dsv_ds

The partial derivative of specific volume with salinity [m3 kg-1 ppt-1].

start

Starting index within the array

npts

The number of values to calculate

eos

Equation of state structure

subroutine, public calculate_compress(
    T T,
    S S,
    pressure pressure,
    rho rho,
    drho_dp drho_dp,
    start start,
    npts npts,
    EOS EOS
    )

Calls the appropriate subroutine to calculate the density and compressibility for 1-D array inputs.

Parameters:

t

Potential temperature referenced to the surface [degC]

s

Salinity [ppt]

pressure

Pressure [Pa]

rho

In situ density [kg m-3].

drho_dp

The partial derivative of density with pressure (also the inverse of the square of sound speed) in s2 m-2.

start

Starting index within the array

npts

The number of values to calculate

eos

Equation of state structure

subroutine, public int_specific_vol_dp(
    T T,
    S S,
    p_t p_t,
    p_b p_b,
    alpha_ref alpha_ref,
    HI HI,
    EOS EOS,
    dza dza,
    intp_dza intp_dza,
    intx_dza intx_dza,
    inty_dza inty_dza,
    halo_size halo_size,
    bathyP bathyP,
    dP_tiny dP_tiny,
    useMassWghtInterp useMassWghtInterp
    )

Calls the appropriate subroutine to alculate analytical and nearly-analytical integrals in pressure across layers of geopotential anomalies, which are required for calculating the finite-volume form pressure accelerations in a non-Boussinesq model. There are essentially no free assumptions, apart from the use of Bode’s rule to do the horizontal integrals, and from a truncation in the series for log(1-eps/1+eps) that assumes that |eps| < .

Parameters:

hi

The horizontal index structure

t

Potential temperature referenced to the surface [degC]

s

Salinity [ppt]

p_t

Pressure at the top of the layer [Pa].

p_b

Pressure at the bottom of the layer [Pa].

alpha_ref

A mean specific volume that is subtracted out to reduce the magnitude of each of the integrals, m3 kg-1. The calculation is mathematically identical with different values of alpha_ref, but this reduces the effects of roundoff.

eos

Equation of state structure

dza

The change in the geopotential anomaly across

intp_dza

The integral in pressure through the layer of the

intx_dza

The integral in x of the difference between the

inty_dza

The integral in y of the difference between the

halo_size

The width of halo points on which to calculate dza.

bathyp

The pressure at the bathymetry [Pa]

dp_tiny

A miniscule pressure change with the same units as p_t (Pa?)

usemasswghtinterp

If true, uses mass weighting to interpolate T/S for top and bottom integrals.

subroutine, public int_density_dz(
    T T,
    S S,
    z_t z_t,
    z_b z_b,
    rho_ref rho_ref,
    rho_0 rho_0,
    G_e G_e,
    HII HII,
    HIO HIO,
    EOS EOS,
    dpa dpa,
    intz_dpa intz_dpa,
    intx_dpa intx_dpa,
    inty_dpa inty_dpa,
    bathyT bathyT,
    dz_neglect dz_neglect,
    useMassWghtInterp useMassWghtInterp
    )

This subroutine calculates analytical and nearly-analytical integrals of pressure anomalies across layers, which are required for calculating the finite-volume form pressure accelerations in a Boussinesq model.

Parameters:

hii

Ocean horizontal index structures for the input arrays

hio

Ocean horizontal index structures for the output arrays

t

Potential temperature referenced to the surface [degC]

s

Salinity [ppt]

z_t

Height at the top of the layer in depth units [Z ~> m].

z_b

Height at the bottom of the layer [Z ~> m].

rho_ref

A mean density [kg m-3], that is subtracted out to reduce the magnitude of each of the integrals.

rho_0

A density [kg m-3], that is used to calculate the pressure (as p~=-z*rho_0*G_e) used in the equation of state.

g_e

The Earth’s gravitational acceleration [m2 Z-1 s-2 ~> m s-2].

eos

Equation of state structure

dpa

The change in the pressure anomaly across the layer [Pa].

intz_dpa

The integral through the thickness of the layer of

intx_dpa

The integral in x of the difference between the

inty_dpa

The integral in y of the difference between the

bathyt

The depth of the bathymetry [Z ~> m].

dz_neglect

A miniscule thickness change [Z ~> m].

usemasswghtinterp

If true, uses mass weighting to interpolate T/S for top and bottom integrals.

logical function, public query_compressible(EOS EOS)

Returns true if the equation of state is compressible (i.e. has pressure dependence)

Parameters:

eos

Equation of state structure

subroutine, public eos_init(param_file param_file, EOS EOS)

Initializes EOS_type by allocating and reading parameters.

Parameters:

param_file

Parameter file structure

eos

Equation of state structure

subroutine, public eos_manual_init(
    EOS EOS,
    form_of_EOS form_of_EOS,
    form_of_TFreeze form_of_TFreeze,
    EOS_quadrature EOS_quadrature,
    Compressible Compressible,
    Rho_T0_S0 Rho_T0_S0,
    drho_dT drho_dT,
    dRho_dS dRho_dS,
    TFr_S0_P0 TFr_S0_P0,
    dTFr_dS dTFr_dS,
    dTFr_dp dTFr_dp
    )

Manually initialized an EOS type (intended for unit testing of routines which need a specific EOS)

Parameters:

eos

Equation of state structure

form_of_eos

A coded integer indicating the equation of state to use.

form_of_tfreeze

A coded integer indicating the expression for the potential temperature of the freezing point.

eos_quadrature

If true, always use the generic (quadrature) code for the integrals of density.

compressible

If true, in situ density is a function of pressure.

rho_t0_s0

Density at T=0 degC and S=0 ppt [kg m-3]

drho_dt

Partial derivative of density with temperature in [kg m-3 degC-1]

drho_ds

Partial derivative of density with salinity in [kg m-3 ppt-1]

tfr_s0_p0

The freezing potential temperature at S=0, P=0 [degC].

dtfr_ds

The derivative of freezing point with salinity in [degC ppt-1].

dtfr_dp

The derivative of freezing point with pressure in [degC Pa-1].

subroutine, public eos_allocate(EOS EOS)

Allocates EOS_type.

Parameters:

eos

Equation of state structure

subroutine, public eos_end(EOS EOS)

Deallocates EOS_type.

Parameters:

eos

Equation of state structure

subroutine, public eos_use_linear(
    Rho_T0_S0 Rho_T0_S0,
    dRho_dT dRho_dT,
    dRho_dS dRho_dS,
    EOS EOS,
    use_quadrature use_quadrature
    )

Set equation of state structure (EOS) to linear with given coefficients.

This routine is primarily for testing and allows a local copy of the EOS_type (EOS argument) to be set to use the linear equation of state independent from the rest of the model.

Parameters:

rho_t0_s0

Density at T=0 degC and S=0 ppt [kg m-3]

drho_dt

Partial derivative of density with temperature [kg m-3 degC-1]

drho_ds

Partial derivative of density with salinity [kg m-3 ppt-1]

use_quadrature

If true, always use the generic (quadrature) code for the integrals of density.

eos

Equation of state structure

subroutine, public int_density_dz_generic(
    T T,
    S S,
    z_t z_t,
    z_b z_b,
    rho_ref rho_ref,
    rho_0 rho_0,
    G_e G_e,
    HII HII,
    HIO HIO,
    EOS EOS,
    dpa dpa,
    intz_dpa intz_dpa,
    intx_dpa intx_dpa,
    inty_dpa inty_dpa,
    bathyT bathyT,
    dz_neglect dz_neglect,
    useMassWghtInterp useMassWghtInterp
    )

This subroutine calculates (by numerical quadrature) integrals of pressure anomalies across layers, which are required for calculating the finite-volume form pressure accelerations in a Boussinesq model.

Parameters:

hii

Horizontal index type for input variables.

hio

Horizontal index type for output variables.

t

Potential temperature of the layer [degC].

s

Salinity of the layer [ppt].

z_t

Height at the top of the layer in depth units [Z ~> m].

z_b

Height at the bottom of the layer [Z ~> m].

rho_ref

A mean density [kg m-3], that is subtracted out to reduce the magnitude of each of the integrals.

rho_0

A density [kg m-3], that is used to calculate the pressure (as p~=-z*rho_0*G_e) used in the equation of state.

g_e

The Earth’s gravitational acceleration [m2 Z-1 s-2 ~> m s-2].

eos

Equation of state structure

dpa

The change in the pressure anomaly

intz_dpa

The integral through the thickness of the

intx_dpa

The integral in x of the difference between

inty_dpa

The integral in y of the difference between

bathyt

The depth of the bathymetry [Z ~> m].

dz_neglect

A miniscule thickness change [Z ~> m].

usemasswghtinterp

If true, uses mass weighting to interpolate T/S for top and bottom integrals.

subroutine, public int_density_dz_generic_plm(
    T_t T_t,
    T_b T_b,
    S_t S_t,
    S_b S_b,
    z_t z_t,
    z_b z_b,
    rho_ref rho_ref,
    rho_0 rho_0,
    G_e G_e,
    dz_subroundoff dz_subroundoff,
    bathyT bathyT,
    HII HII,
    HIO HIO,
    EOS EOS,
    dpa dpa,
    intz_dpa intz_dpa,
    intx_dpa intx_dpa,
    inty_dpa inty_dpa,
    useMassWghtInterp useMassWghtInterp
    )

Compute pressure gradient force integrals by quadrature for the case where T and S are linear profiles.

Parameters:

hii

Ocean horizontal index structures for the input arrays

hio

Ocean horizontal index structures for the output arrays

t_t

Potential temperatue at the cell top [degC]

t_b

Potential temperatue at the cell bottom [degC]

s_t

Salinity at the cell top [ppt]

s_b

Salinity at the cell bottom [ppt]

z_t

The geometric height at the top of the layer,

z_b

The geometric height at the bottom of the layer [Z ~> m].

rho_ref

A mean density [kg m-3], that is subtracted out to reduce the magnitude of each of the integrals.

rho_0

A density [kg m-3], that is used to calculate the pressure (as p~=-z*rho_0*G_e) used in the equation of state.

g_e

The Earth’s gravitational acceleration [m2 Z-1 s-2 ~> m s-2].

dz_subroundoff

A miniscule thickness change [Z ~> m].

bathyt

The depth of the bathymetry [Z ~> m].

eos

Equation of state structure

dpa

The change in the pressure anomaly across the layer [Pa].

intz_dpa

The integral through the thickness of the layer of

intx_dpa

The integral in x of the difference between the

inty_dpa

The integral in y of the difference between the

usemasswghtinterp

If true, uses mass weighting to interpolate T/S for top and bottom integrals.

subroutine, public find_depth_of_pressure_in_cell(
    T_t T_t,
    T_b T_b,
    S_t S_t,
    S_b S_b,
    z_t z_t,
    z_b z_b,
    P_t P_t,
    P_tgt P_tgt,
    rho_ref rho_ref,
    G_e G_e,
    EOS EOS,
    P_b P_b,
    z_out z_out,
    z_tol z_tol
    )

Find the depth at which the reconstructed pressure matches P_tgt.

Parameters:

t_t

Potential temperatue at the cell top [degC]

t_b

Potential temperatue at the cell bottom [degC]

s_t

Salinity at the cell top [ppt]

s_b

Salinity at the cell bottom [ppt]

z_t

Absolute height of top of cell [Z ~> m]. (Boussinesq ????)

z_b

Absolute height of bottom of cell [Z ~> m].

p_t

Anomalous pressure of top of cell, relative to g*rho_ref*z_t [Pa]

p_tgt

Target pressure at height z_out, relative to g*rho_ref*z_out [Pa]

rho_ref

Reference density with which calculation are anomalous to

g_e

Gravitational acceleration [m2 Z-1 s-2 ~> m s-2]

eos

Equation of state structure

p_b

Pressure at the bottom of the cell [Pa]

z_out

Absolute depth at which anomalous pressure = p_tgt [Z ~> m].

z_tol

The tolerance in finding z_out [Z ~> m].

subroutine, public int_density_dz_generic_ppm(
    T T,
    T_t T_t,
    T_b T_b,
    S S,
    S_t S_t,
    S_b S_b,
    z_t z_t,
    z_b z_b,
    rho_ref rho_ref,
    rho_0 rho_0,
    G_e G_e,
    HII HII,
    HIO HIO,
    EOS EOS,
    dpa dpa,
    intz_dpa intz_dpa,
    intx_dpa intx_dpa,
    inty_dpa inty_dpa
    )

Compute pressure gradient force integrals for the case where T and S are parabolic profiles.

Parameters:

hii

Ocean horizontal index structures for the input arrays

hio

Ocean horizontal index structures for the output arrays

t

Potential temperature referenced to the surface [degC]

t_t

Potential temperatue at the cell top [degC]

t_b

Potential temperatue at the cell bottom [degC]

s

Salinity [ppt]

s_t

Salinity at the cell top [ppt]

s_b

Salinity at the cell bottom [ppt]

z_t

Height at the top of the layer [Z ~> m].

z_b

Height at the bottom of the layer [Z ~> m].

rho_ref

A mean density [kg m-3], that is subtracted out to reduce the magnitude of each of the integrals.

rho_0

A density [kg m-3], that is used to calculate the pressure (as p~=-z*rho_0*G_e) used in the equation of state.

g_e

The Earth’s gravitational acceleration [m s-2]

eos

Equation of state structure

dpa

The change in the pressure anomaly across the layer [Pa].

intz_dpa

The integral through the thickness of the layer of

intx_dpa

The integral in x of the difference between the

inty_dpa

The integral in y of the difference between the

subroutine, public int_spec_vol_dp_generic(
    T T,
    S S,
    p_t p_t,
    p_b p_b,
    alpha_ref alpha_ref,
    HI HI,
    EOS EOS,
    dza dza,
    intp_dza intp_dza,
    intx_dza intx_dza,
    inty_dza inty_dza,
    halo_size halo_size,
    bathyP bathyP,
    dP_neglect dP_neglect,
    useMassWghtInterp useMassWghtInterp
    )

This subroutine calculates integrals of specific volume anomalies in pressure across layers, which are required for calculating the finite-volume form pressure accelerations in a non-Boussinesq model. There are essentially no free assumptions, apart from the use of Bode’s rule quadrature to do the integrals.

Parameters:

hi

A horizontal index type structure.

t

Potential temperature of the layer [degC].

s

Salinity of the layer [ppt].

p_t

Pressure atop the layer [Pa].

p_b

Pressure below the layer [Pa].

alpha_ref

A mean specific volume that is subtracted out to reduce the magnitude of each of the integrals [m3 kg-1]. The calculation is mathematically identical with different values of alpha_ref, but alpha_ref alters the effects of roundoff, and answers do change.

eos

Equation of state structure

dza

The change in the geopotential anomaly

intp_dza

The integral in pressure through the

intx_dza

The integral in x of the difference

inty_dza

The integral in y of the difference

halo_size

The width of halo points on which to calculate dza.

bathyp

The pressure at the bathymetry [Pa]

dp_neglect

A miniscule pressure change with the same units as p_t (Pa?)

usemasswghtinterp

If true, uses mass weighting to interpolate T/S for top and bottom integrals.

subroutine, public int_spec_vol_dp_generic_plm(
    T_t T_t,
    T_b T_b,
    S_t S_t,
    S_b S_b,
    p_t p_t,
    p_b p_b,
    alpha_ref alpha_ref,
    dP_neglect dP_neglect,
    bathyP bathyP,
    HI HI,
    EOS EOS,
    dza dza,
    intp_dza intp_dza,
    intx_dza intx_dza,
    inty_dza inty_dza,
    useMassWghtInterp useMassWghtInterp
    )

This subroutine calculates integrals of specific volume anomalies in pressure across layers, which are required for calculating the finite-volume form pressure accelerations in a non-Boussinesq model. There are essentially no free assumptions, apart from the use of Bode’s rule quadrature to do the integrals.

Parameters:

hi

A horizontal index type structure.

t_t

Potential temperature at the top of the layer [degC].

t_b

Potential temperature at the bottom of the layer [degC].

s_t

Salinity at the top the layer [ppt].

s_b

Salinity at the bottom the layer [ppt].

p_t

Pressure atop the layer [Pa].

p_b

Pressure below the layer [Pa].

alpha_ref

A mean specific volume that is subtracted out to reduce the magnitude of each of the integrals [m3 kg-1]. The calculation is mathematically identical with different values of alpha_ref, but alpha_ref alters the effects of roundoff, and answers do change.

dp_neglect

A miniscule pressure change with the same units as p_t (Pa?)

bathyp

The pressure at the bathymetry [Pa]

eos

Equation of state structure

dza

The change in the geopotential anomaly

intp_dza

The integral in pressure through the

intx_dza

The integral in x of the difference

inty_dza

The integral in y of the difference

usemasswghtinterp

If true, uses mass weighting to interpolate T/S for top and bottom integrals.

subroutine, public convert_temp_salt_for_teos10(
    T T,
    S S,
    press press,
    G G,
    kd kd,
    mask_z mask_z,
    EOS EOS
    )

Convert T&S to Absolute Salinity and Conservative Temperature if using TEOS10.

Parameters:

g

The ocean’s grid structure

t

Potential temperature referenced to the surface [degC]

s

Salinity [ppt]

press

Pressure at the top of the layer [Pa].

eos

Equation of state structure

mask_z

3d mask regulating which points to convert.

kd

The number of layers to work on

subroutine, public extract_member_eos(
    EOS EOS,
    form_of_EOS form_of_EOS,
    form_of_TFreeze form_of_TFreeze,
    EOS_quadrature EOS_quadrature,
    Compressible Compressible,
    Rho_T0_S0 Rho_T0_S0,
    drho_dT drho_dT,
    dRho_dS dRho_dS,
    TFr_S0_P0 TFr_S0_P0,
    dTFr_dS dTFr_dS,
    dTFr_dp dTFr_dp
    )

Extractor routine for the EOS type if the members need to be accessed outside this module.

Parameters:

eos

Equation of state structure

form_of_eos

A coded integer indicating the equation of state to use.

form_of_tfreeze

A coded integer indicating the expression for the potential temperature of the freezing point.

eos_quadrature

If true, always use the generic (quadrature) code for the integrals of density.

compressible

If true, in situ density is a function of pressure.

rho_t0_s0

Density at T=0 degC and S=0 ppt [kg m-3]

drho_dt

Partial derivative of density with temperature in [kg m-3 degC-1]

drho_ds

Partial derivative of density with salinity in [kg m-3 ppt-1]

tfr_s0_p0

The freezing potential temperature at S=0, P=0 [degC].

dtfr_ds

The derivative of freezing point with salinity [degC PSU-1].

dtfr_dp

The derivative of freezing point with pressure [degC Pa-1].