SAMSIM
|
Contains subroutines and functions related to multi-phase thermodynamics. More...
Functions/Subroutines | |
subroutine, public | gett (H, S_bu, T_in, T, phi, k) |
Determines equilibrium Temperature of a layer for given S_bu and H as well as solid fraction. More... | |
subroutine, public | expulsion (phi, thick, m, psi_s, psi_l, psi_g, V_ex) |
Determines Brine flux expelled from out of a layer due to freezing. More... | |
subroutine, public | sub_fl_q (psi_s_1, psi_l_1, psi_g_1, thick_1, T_1, psi_s_2, psi_l_2, psi_g_2, thick_2, T_2, fl_Q) |
Determines conductive heat flux between two layers. More... | |
subroutine, public | sub_fl_q_0 (psi_s, psi_l, psi_g, thick, T, T_bound, direct_flag, fl_Q) |
Determines conductive Heat flux between layer and boundary temperatures. More... | |
subroutine, public | sub_fl_q_styropor (k_styropor, fl_Q) |
Niels, 2017 add: Determines conductive Heat flux below styropor cover. More... | |
real(wp) function, public | func_s_br (T, S_bu) |
Computes salinity of brine pockets for given temperature in Celsius of mushy layer. More... | |
real(wp) function, public | func_ddt_s_br (T) |
Computes temperature derivative of brine pocket salinity for given temperature in Celsius of mushy layer. More... | |
Contains subroutines and functions related to multi-phase thermodynamics.
See the subroutine and function descriptions for details.
COPYRIGHT
This file is part of SAMSIM.
SAMSIM is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
SAMSIM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
You should have received a copy of the GNU General Public License along with SAMSIM. If not, see http://www.gnu.org/licenses/.
subroutine, public mo_thermo_functions::expulsion | ( | real(wp), intent(in) | phi, |
real(wp), intent(in) | thick, | ||
real(wp), intent(inout) | m, | ||
real(wp), intent(out) | psi_s, | ||
real(wp), intent(out) | psi_l, | ||
real(wp), intent(out) | psi_g, | ||
real(wp), intent(out) | V_ex | ||
) |
Determines Brine flux expelled from out of a layer due to freezing.
If the volume of ice and brine exceed the Volume of the layer brine is expelled. The volume of the ejected brine is calculated and exported. The volume fractions are also calculated.
real(wp) function, public mo_thermo_functions::func_ddt_s_br | ( | real(wp), intent(in) | T | ) |
Computes temperature derivative of brine pocket salinity for given temperature in Celsius of mushy layer.
Subroutine computes one ddT_S_br for one given T NaCl solutions and seawater produce slight variations. Which solution is used is specified by salt_flag. Based on Notz 2005 p. 36 ddT_S_br = c2+2*c3*T+3*c4*T^2
For T below -20 we simply use a linear extension based on "Composition of sea ice and its tensile strength". The actual salinity function is much more complicated and depends on the salt composition, but the linear fit is far better than using the polynomial fit.
The NaCL precipitates at -22, leading to ice/salt kristall mix below -22. Ideally the whole code would be modified to take the non-continous transition at -22 but given that there is currently little interest I can't be bothered to put in the effor.
[in] | t | Temperature in Celsius |
real(wp) function, public mo_thermo_functions::func_s_br | ( | real(wp), intent(in) | T, |
real(wp), intent(in), optional | S_bu | ||
) |
Computes salinity of brine pockets for given temperature in Celsius of mushy layer.
Subroutine computes one S_br for one given T in Celsius by third-order polynomial. NaCl solutions and seawater produce slight variations. Which solution is used is determined by salt_flag. S_br = c1+c2*T+c3*T^2+c4*T^3 Based on Notz 2005 p. 36
For T below -20 we simply use a linear extension based on "Composition of sea ice and its tensile strength". The actual salinity function is much more complicated and depends on the salt composition, but the linear fit is far better than using the polynomial fit.
The NaCL precipitates at -22, leading to ice/salt kristall mix below -22. Ideally the whole code would be modified to take the non-continous transition at -22 but given that there is currently little interest I can't be bothered to put in the effor.
[in] | t | Temperature in Celsius |
subroutine, public mo_thermo_functions::gett | ( | real(wp), intent(in) | H, |
real(wp), intent(in) | S_bu, | ||
real(wp), intent(in) | T_in, | ||
real(wp), intent(out) | T, | ||
real(wp), intent(out) | phi, | ||
integer, intent(in), optional | k | ||
) |
Determines equilibrium Temperature of a layer for given S_bu and H as well as solid fraction.
The temperature of a fully liquid layer is used to see if the resulting brine salinity is lower than the bulk salinity. After checking if the layer is a fluid or a mushy layer the temperature is calculated by solving f(T) = 0 using the Newton method. f(T) = -latent_heat-H+latent_heat*S_bu/S_br(T) + c_s*T+c_s_beta*T^2/2 f'(T) = c_s+c_s_beta*T-latent_heat*S_bu*S_br'(T)/S_br^2 Described in Notz2005, subsubsection 5.6.1. See func_S_br(T) and func_ddT_S_br(T). First guess T_0 must be given, low first guess lead to overshooting which would lead to very high Temperatures. To avoid this, an if loop sets T to freezing T when T>0. Freezing T is also calculated at the beginning using the Newton-Method. If S_bu<0.001 then it is treated as pure ice.
[in] | h | Enthalpy [J/kg] |
[in] | s_bu | Bulk Salinity [g/kg] |
[in] | t_in | input Temperature for T_0 [C] |
[out] | t | Temperature [C] |
[out] | phi | solid fraction |
subroutine, public mo_thermo_functions::sub_fl_q | ( | real(wp), intent(in) | psi_s_1, |
real(wp), intent(in) | psi_l_1, | ||
real(wp), intent(in) | psi_g_1, | ||
real(wp), intent(in) | thick_1, | ||
real(wp), intent(in) | T_1, | ||
real(wp), intent(in) | psi_s_2, | ||
real(wp), intent(in) | psi_l_2, | ||
real(wp), intent(in) | psi_g_2, | ||
real(wp), intent(in) | thick_2, | ||
real(wp), intent(in) | T_2, | ||
real(wp), intent(out) | fl_Q | ||
) |
Determines conductive heat flux between two layers.
Details can be found in Notz 2005, especially equation 5.7. The gas volume is assumed to have no thermal properties at all. First the thermal resistance R is calculated using the approximated thermal conductivity of the mushy layer (see Notz 2005 eq. 3.41.). Then the heat flux Q is simply (T_1-T_2)/R "_1" denotes the upper layer and "_2" the lower layer. A positive heat flux is from lower to upper layer.
subroutine, public mo_thermo_functions::sub_fl_q_0 | ( | real(wp), intent(in) | psi_s, |
real(wp), intent(in) | psi_l, | ||
real(wp), intent(in) | psi_g, | ||
real(wp), intent(in) | thick, | ||
real(wp), intent(in) | T, | ||
real(wp), intent(in) | T_bound, | ||
integer | direct_flag, | ||
real(wp), intent(out) | fl_Q | ||
) |
Determines conductive Heat flux between layer and boundary temperatures.
Details can be found in Notz 2005, especially equation 5.10 and 5.11. The gas volume is assumed to have no thermal properties. direct_flag denotes if the boundary layer is above or below the layer. 1 : = layer above boundary -1: = layer below boundary
[in] | t_bound | T_bound temperature of boundary layer |
subroutine, public mo_thermo_functions::sub_fl_q_styropor | ( | real(wp), intent(in) | k_styropor, |
real(wp), intent(inout) | fl_Q | ||
) |
Niels, 2017 add: Determines conductive Heat flux below styropor cover.
Standard approach.