Chamber.jl
This is the API documentation for Chamber.jl. See the package's README for instructions about general usage.
Functions
Chamber.IC_Finder
— MethodIC_Finder(composition::Mafic, M_h2o::Float64, M_co2::Float64, M_tot::Float64, P::Float64, T::Float64, V::Float64, rho_m::Float64, param_IC::ParamICFinder{Float64})::NamedTuple{(:eps_g0, :X_co20, :mco2_diss, :phase),NTuple{4,Float64}}
An iterative function that uses thermodynamic modeling to determine the gas phase composition of a mafic magma chamber.
Arguments
composition
:Mafic()
M_h2o
: total mass of water in magma (kg)M_co2
: total mass of CO2 in magma (kg)M_tot
: total mass of magma (kg)P
: Pressure (Pa)T
: Temperature (K)V
: chamber volume (m³)rho_m
: density of the melt (kg/m³)param_IC
: An instance of theParamICFinder
composite type containing the parameters for IC_Finder
Returns
A NamedTuple with the following fields:
eps_g0
: Initial gas fraction of the magma chamberX_co20
: Initial mole fraction of CO2 in the gas phasemco2_diss
: Total mass of CO2 dissolved in the melt (kg)phase
: Integer representing the state of the magma
Details
The function uses a thermodynamic model to find the initial gas fraction (eps_g0
) and initial mole fraction of CO2 in the gas phase (X_co20
) for a mafic magma chamber. The model assumes that the chamber contains two phases: a liquid melt and a gas phase. The model calculates the saturation state of the magma with respect to CO2, and determines whether the gas phase contains only CO2 or a mixture of CO2 and H2O. If the chamber is in two-phase state, the model returns eps_g0 = 0.0
, X_co20 = 0.0
, mco2_diss = Total mass of CO2 dissolved in the melt
, and phase = 2
. If the chamber is in gas or liquid state, the function iteratively solves for eps_g0
and X_co20
using the solve_X_co2
and get_eps_g
helper functions until the relative difference between eps_g0
and X_co20
in two consecutive iterations is less than the tolerance Tol
, or until the maximum number of iterations max_count
is reached.
Chamber.IC_Finder
— MethodIC_Finder(composition::Silicic, M_h2o::Float64, M_co2::Float64, M_tot::Float64, P::Float64, T::Float64, V::Float64, rho_m::Float64, param_IC::ParamICFinder{Float64})::NamedTuple{(:eps_g0, :X_co20, :mco2_diss, :phase),NTuple{4,Float64}}
An iterative function that uses thermodynamic modeling to determine the gas phase composition of a silicic magma chamber.
Arguments
composition
:Silicic()
M_h2o
: total mass of water in magma (kg)M_co2
: total mass of CO2 in magma (kg)M_tot
: total mass of magma (kg)P
: Pressure (Pa)T
: Temperature (K)V
: chamber volume (m³)rho_m
: density of the melt (kg/m³)param_IC
: An instance of theParamICFinder
composite type containing the parameters for IC_Finder
Returns
A NamedTuple with the following fields:
eps_g0
: Initial gas fraction of the magma chamberX_co20
: Initial mole fraction of CO2 in the gas phasemco2_diss
: Total mass of CO2 dissolved in the melt (kg)phase
: Integer representing the state of the magma
Details
The function uses a thermodynamic model to find the initial gas fraction (eps_g0
) and initial mole fraction of CO2 in the gas phase (X_co20
) for a silicic magma chamber. The model assumes that the chamber contains two phases: a liquid melt and a gas phase. The model calculates the saturation state of the magma with respect to CO2, and determines whether the gas phase contains only CO2 or a mixture of CO2 and H2O. If the chamber is in two-phase state, the model returns eps_g0 = 0.0
, X_co20 = 0.0
, mco2_diss = Total mass of CO2 dissolved in the melt
, and phase = 2
. If the chamber is in gas or liquid state, the function iteratively solves for eps_g0
and X_co20
using the solve_X_co2
and get_eps_g
helper functions until the relative difference between eps_g0
and X_co20
in two consecutive iterations is less than the tolerance Tol
, or until the maximum number of iterations max_count
is reached.
Chamber.a1x_f
— MethodThis is for a11 and a12 a11: rho, drhodP, V, dVdP a12: rho, drhodT, V, dVdT
Chamber.affect!
— Methodaffect!(int, idx, sw::SW{Int8}, param::Param{Float64}, param_saved_var::ParamSaved{Float64}, param_IC_Finder::ParamICFinder{Float64})
Re-initialize the condition when the event happens. This function modifies the current state of the integrator (int.u
) when a particular event occurs during the simulation. The function adjusts various parameters based on the current state of the integrator and the custom parameters that were passed in.
Arguments
int
: The current state of the integrator. It's format is from the DifferentialEquations.jl packageidx
: The index of the event that caused the function to be called.sw
: A custom parameter used to control simulation behavior.param
: A custom parameter containing physical constants and other model parameters.param_saved_var
: A custom parameter used to store values from the previous time step.param_IC_Finder
: A custom parameter used to control the behavior of the IC_Finder function.
The arguments int
and idx
are from the DifferentialEquations.jl package. These argument formats are specific to the DifferentialEquations.jl package.
Chamber.boundary_conditions_new
— Methodboundary_conditions_new(P::Float64, T::Float64, V::Float64, rho_m::Float64, rho_x::Float64, c::Float64, sw::SW{Int8}, T_in::Float64, M_h2o::Float64, M_co2::Float64, total_Mass::Float64, param::Param{Float64}, param_saved_var::ParamSaved{Float64})
Arguments
P
: Pressure (Pa)T
: Temperature (K)V
: chamber volume (m³)rho_m
: density of meltrho_x
: density of crystal of magmac
: heat of magmasw
: eruption/coolingmodule/viscousrelaxation controlT_in
: TemperatureM_h2o
: total mess of H2O in the magmaM_co2
: total mess of CO2 in the magmatotal_Mass
: total mess of magma chamber
Chamber.build_co2
— Methodbuild_co2(Pw::T, Pc::T, Temp::T, dPwdP::T, dPcdP::T, dPwdXco2::T, dPcdXco2::T)::NamedTuple{(:C_co2, :dC_co2dT, :dC_co2dP, :dC_co2dXco2),NTuple{4,T}} where {T<:Float64}
- This function is used within the
exsolve
function.
Returns
A NamedTuple with the following fields:
C_co2
dC_co2dT
dC_co2dP
dC_co2dXco2
Chamber.build_mdot_in
— Methodbuild_mdot_in(fluxing::Bool, rho_m0::Float64, log_vfr::Float64, P_0::Float64, T_in::Float64)::Float64
- This function is used within the
chamber
function.
Returns
mdot_in
: mass inflow rate
Chamber.build_meq
— Methodbuild_meq(composition::Mafic, P::T, Temp::T, X_co2::T)::NamedTuple{(:meq, :dmeqdT, :dmeqdP, :dmeqdXco2),NTuple{4,T}} where {T<:Float64}
- This function is used within the
exsolve
function.
Returns
A NamedTuple with the following fields:
meq
dmeqdT
dmeqdP
dmeqdXco2
Chamber.build_meq
— Methodbuild_meq(composition::Silicic, Pw::T, Pc::T, Temp::T, dPwdP::T, dPcdP::T, dPwdXco2::T, dPcdXco2::T)::NamedTuple{(:meq, :dmeqdT, :dmeqdP, :dmeqdXco2),NTuple{4,T}} where {T<:Float64}
- This function is used within the
exsolve
function.
Returns
A NamedTuple with the following fields:
meq
dmeqdT
dmeqdP
dmeqdXco2
Chamber.build_rho_rc
— Methodbuild_rho_rc(eps_m::T, eps_g::T, eps_x::T, rho_m::T, rho_g::T, rho_x::T, drho_m_dP::T, drho_g_dP::T, drho_x_dP::T, drho_m_dT::T, drho_g_dT::T, drho_x_dT::T, c_x::T, c_m::T, c_g::T, deps_x_dP::T, deps_x_dT::T)::Vector{T} where {T<:Float64}
- This function is used within the
odeChamber
function.
Returns
[rho
, drho_dP
, drho_dT
, drho_deps_g
, rc
, drc_dP
, drc_dT
]
Chamber.chamber
— Functionchamber(composition::Union{Silicic,Mafic}, end_time::Float64, log_volume_km3::Float64, InitialConc_H2O::Float64, InitialConc_CO2::Float64, log_vfr::Float64, depth::Float64, output_dirname::String; kwargs...)
Simulate the eruption of a volcano using a model for the frequency of eruptions of upper crustal magma chambers based on Degruyter and Huber (2014).
Arguments
composition
: The magma composition. UseSilicic()
for rhyolite composition (arc magma) orMafic()
for basalt composition (rift).end_time
: Maximum magma chamber evolution duration in seconds.log_volume_km3
: The initial volume of the chamber in logarithmic scale. The actual initial chamber volume is calculated as 10^(log_volume_km3
) in km³.InitialConc_H2O
: The initial weight fraction of water in the magma (exsolved + dissolved).InitialConc_CO2
: The initial weight fraction of CO₂ in the magma (exsolved + dissolved).log_vfr
: Magma recharge rate in km³/yr calculated as 10^(log_vfr
).depth
: Depth of the magma chamber in meters.output_dirname
(optional): Name of the output directory. Defaults to current timestamp.
Keyword Arguments
plotfig
(optional): (default:true
). Generate and plot figures for each result if true.
Returns
A DataFrame
containing the solution with columns:
time
: Simulation timestamps in seconds.P+dP
: Pressure in Pa.T
: Temperature in K.eps_g
: Gas volume fraction.V
: Volume of the magma chamber in m³.rho_m
: Density of the melt in kg/m³.rho_x
: Density of magma crystal in kg/m³.X_CO2
: Mole fraction of CO2 in the gas.total_mass
: Total mass of magma chamber in kg.total_mass_H2O
: Total mass of water in the magma in kg.total_mass_CO2
: Total mass of CO₂ in the magma in kg.eps_x
: Crystal volume fraction.
Outputs
A directory named after output_dirname
or the default value, containing the following files:
out.csv
: a CSV file containing the solution columns listed above.eruptions.csv
, A CSV file containing the datas of eruptions with the following columns: timeoferuption (sec), durationoferuption (sec), masserupted (kg) and volumeerupted (km³).- Figures for P+dP(t), T(t), epsg(t), V(t), XCO2(t), total_mass(t).
References
- W. Degruyter and C. Huber (2014). A model for eruption frequency of upper crustal silicic magma chambers. Earth Planet. Sci. Lett. https://doi.org/10.1016/j.epsl.2014.06.047
Examples
# Run a simulation with silicic magma chamber
julia> composition = Silicic()
julia> end_time = 3e9
julia> log_volume_km3 = 0.2
julia> InitialConc_H2O = 0.04
julia> InitialConc_CO2 = 0.001
julia> log_vfr = -3.3
julia> depth = 8e3
julia> chamber(composition, end_time, log_volume_km3, InitialConc_H2O, InitialConc_CO2, log_vfr, depth)
# Run a simulation with mafic magma chamber, with custom directory name "MyDirname"
julia> composition = Mafic()
julia> end_time = 3e9
julia> log_volume_km3 = 0.2
julia> InitialConc_H2O = 0.01
julia> InitialConc_CO2 = 0.001
julia> log_vfr = -3.3
julia> depth = 8e3
julia> output_dirname = "MyDirname"
julia> chamber(composition, end_time, log_volume_km3, InitialConc_H2O, InitialConc_CO2, log_vfr, depth, output_dirname)
Chamber.chamber
— Functionchamber(composition::Union{Silicic,Mafic}, end_time::Float64, log_volume_km3_vector::Union{Float64,Vector{Float64}}, InitialConc_H2O_vector::Union{Float64,Vector{Float64}}, InitialConc_CO2_vector::Union{Float64,Vector{Float64}}, log_vfr_vector::Union{Float64,Vector{Float64}}, depth_vector::Union{Float64,Vector{Float64}}, output_dirname::String; kwargs...)
Simulate the eruption of a volcano using a model for the frequency of eruptions of upper crustal magma chambers based on Degruyter and Huber (2014).
Arguments
composition
: The magma composition. UseSilicic()
for rhyolite composition (arc magma) orMafic()
for basalt composition (rift).end_time
: Maximum magma chamber evolution duration in seconds.log_volume_km3
: The initial volume of the chamber in logarithmic scale. The actual initial chamber volume is calculated as 10^(log_volume_km3
) in km³.InitialConc_H2O
: The initial weight fraction of water in the magma (exsolved + dissolved).InitialConc_CO2
: The initial weight fraction of CO₂ in the magma (exsolved + dissolved).log_vfr
: Magma recharge rate in km³/yr calculated as 10^(log_vfr
).depth
: Depth of the magma chamber in meters.output_dirname
(optional): Name of the output directory. Defaults to current timestamp.
Keyword Arguments
plotfig
(optional): (default:true
). Generate and plot figures for each result if true.
Returns
A DataFrame
containing the solution with columns:
time
: Simulation timestamps in seconds.P+dP
: Pressure in Pa.T
: Temperature in K.eps_g
: Gas volume fraction.V
: Volume of the magma chamber in m³.rho_m
: Density of the melt in kg/m³.rho_x
: Density of magma crystal in kg/m³.X_CO2
: Mole fraction of CO2 in the gas.total_mass
: Total mass of magma chamber in kg.total_mass_H2O
: Total mass of water in the magma in kg.total_mass_CO2
: Total mass of CO₂ in the magma in kg.eps_x
: Crystal volume fraction.
Outputs
A directory named after output_dirname
or the default value, containing the following files:
out.csv
: a CSV file containing the solution columns listed above.eruptions.csv
, A CSV file containing the datas of eruptions with the following columns: timeoferuption (sec), durationoferuption (sec), masserupted (kg) and volumeerupted (km³).- Figures for P+dP(t), T(t), epsg(t), V(t), XCO2(t), total_mass(t).
References
- W. Degruyter and C. Huber (2014). A model for eruption frequency of upper crustal silicic magma chambers. Earth Planet. Sci. Lett. https://doi.org/10.1016/j.epsl.2014.06.047
Examples
# Run a simulation with silicic magma chamber
julia> composition = Silicic()
julia> end_time = 3e9
julia> log_volume_km3 = 0.2
julia> InitialConc_H2O = [0.04, 0.045]
julia> InitialConc_CO2 = [0.001, 0.0012]
julia> log_vfr = -3.3
julia> depth = 8e3
julia> chamber(composition, end_time, log_volume_km3, InitialConc_H2O, InitialConc_CO2, log_vfr, depth)
# Run a simulation with mafic magma chamber, with custom directory name "MyDirname"
julia> composition = Mafic()
julia> end_time = 3e9
julia> log_volume_km3 = 0.2
julia> InitialConc_H2O = [0.01, 0.012]
julia> InitialConc_CO2 = 0.001
julia> log_vfr = -3.3
julia> depth = [7e3, 8e3]
julia> output_dirname = "MyDirname"
julia> chamber(composition, end_time, log_volume_km3, InitialConc_H2O, InitialConc_CO2, log_vfr, depth, output_dirname)
Chamber.check_for_duplicates
— Methodcheck_for_duplicates(log_volume_km3_vector::Union{Float64,Vector{Float64}}, InitialConc_H2O_vector::Union{Float64,Vector{Float64}}, InitialConc_CO2_vector::Union{Float64,Vector{Float64}}, log_vfr_vector::Union{Float64,Vector{Float64}}, depth_vector::Union{Float64,Vector{Float64}})::Nothing
This function checks if any of the input arrays contain duplicate elements. If duplicates are found, it raises an error indicating which arguments have duplicates.
- This function is used within the
chamber
function
Chamber.compute_dXdP_dXdT
— Methodcompute_dXdP_dXdT(u::Float64, param::Param, var::String)::Tuple{Float64,Float64,Float64}
- This function is used within the
odeChamber
function.
Chamber.crystal_fraction
— Methodcrystal_fraction(composition::Mafic, T::Float64, P::Float64, mH2O::Float64, mCO2::Float64)::NamedTuple{(:eps_x, :deps_x_dP, :deps_x_dT, :deps_x_deps_g, :deps_x_dmco2_t, :deps_x_dmh2o_t), NTuple{6, Float64}}
Calculate crystal volume fraction and its derivatives in magma at given temperature, pressure, and water and CO2 contents.
Arguments
composition
:Mafic()
T
: Temperature (K)P
: Pressure (Pa)mH2O
: Weight fration of the H2O in magma.mCO2
: Weight fration of the CO2 in magma.
Returns
A NamedTuple with the following fields:
eps_x
: Crystal volume fraction.deps_x_dP
: Derivative of eps_x with respect to pressure.deps_x_dT
: Derivative of eps_x with respect to temperature.deps_x_deps_g
: Derivative of eps_x with respect to gas volume fraction.deps_x_dmco2_t
: Derivative of eps_x with respect to CO2 content.deps_x_dmh2o_t
: Derivative of eps_x with respect to H2O content.
Chamber.crystal_fraction
— Methodcrystal_fraction(composition::Silicic, T::Float64, P::Float64, mH2O::Float64, mCO2::Float64)::NamedTuple{(:eps_x, :deps_x_dP, :deps_x_dT, :deps_x_deps_g, :deps_x_dmco2_t, :deps_x_dmh2o_t), NTuple{6, Float64}}
Calculate crystal volume fraction and its derivatives in magma at given temperature, pressure, and water and CO2 contents.
Arguments
composition
:Silicic()
T
: Temperature (K)P
: Pressure (Pa)mH2O
: Weight fration of the H2O in magma.mCO2
: Weight fration of the CO2 in magma.
Returns
A NamedTuple with the following fields:
eps_x
: Crystal volume fraction.deps_x_dP
: Derivative of eps_x with respect to pressure.deps_x_dT
: Derivative of eps_x with respect to temperature.deps_x_deps_g
: Derivative of eps_x with respect to gas volume fraction.deps_x_dmco2_t
: Derivative of eps_x with respect to CO2 content.deps_x_dmh2o_t
: Derivative of eps_x with respect to H2O content.
Chamber.crystal_fraction_eps_x
— Methodcrystal_fraction_eps_x(composition::Mafic, T::Float64, P::Float64, mH2O::Float64, mCO2::Float64)::Float64
- Calculates crystal volume fraction at given temperature, pressure, and water and CO2 contents.
- Spetialized version of
crystal_fraction
that computes
eps_x` only.
Returns
eps_x
: Crystal volume fraction.
Chamber.crystal_fraction_eps_x
— Methodcrystal_fraction_eps_x(composition::Silicic, T::Float64, P::Float64, mH2O::Float64, mCO2::Float64)::Float64
- Calculates crystal volume fraction at given temperature, pressure, and water and CO2 contents.
- Spetialized version of
crystal_fraction
that computeseps_x
only.
Returns
eps_x
: Crystal volume fraction.
Chamber.dM_X_t_dt_f
— MethoddM_X_t_dt_f(rho, V, Mdot_v_in, Mdot_v_out, m_h2o, Mdot_in, Mdot_out)
For dM_h2o_t_dt
& dM_co2_t_dt
- MdotXin:
Mdot_v_in
orMdot_c_in
- MdotXout:
Mdot_v_out
orMdot_c_out
- mX: `mh2o
or
m_co2`
Chamber.dX_dxdydz
— MethoddX_dxdydz(composition::Union{Silicic,Mafic}, s::String, x::T, y::T, z::T)::NamedTuple{(:X, :dXdx, :dXdy, :dXdz),NTuple{4,T}} where {T<:Float64}
- This function is used within the
parameters_melting_curve
function. - Calculate value from different parameter sets (# of parameter sets: Silicic: 3, Mafic: 2)
Arguments
composition
:Silicic()
orMafic()
s
: String that represents the parameter set to use. ForSilicic
,s
can be"a"
,"b"
, or"c"
. ForMafic
,s
can be"a"
or"b"
.x
: Water (H2O)y
: Gas (CO2)z
: Pressure (P)
Returns
A NamedTuple with the following fields:
X
: the calculated value for the given parameter setdXdx
: the derivative ofX
with respect tox
dXdy
: the derivative ofX
with respect toy
dXdz
: the derivative ofX
with respect toz
Chamber.drc_dX_f
— Methoddrc_dX_f(;eps_x::T, c_x::T, drho_x_dX::T, eps_g::T, c_g::T, drho_g_dX::T, eps_m::T, c_m::T, drho_m_dX::T, rho_x::T, rho_m::T, deps_x_dX::T)::T where {T<:Float64}
Computing the derivatives of the product of density and specific heat for the mixture
X
: representP
orT
- This function is used within the
build_rho_rc
function.
Returns
drc_dX
: representdrc_dP
ordrc_dT
Chamber.drho_dX_f
— Methoddrho_dX_f(;eps_m::T, eps_g::T, eps_x::T, drho_m_dX::T, drho_g_dX::T, drho_x_dX::T, rho_x::T, rho_m::T, deps_x_dX::T)::T where {T<:Float64}
X
: representP
orT
- This function is used within the
build_rho_rc
function.
Returns
drho_dX
: representdrho_dP
ordrho_dT
Chamber.dwater_dx
— Methoddwater_dx(composition::Mafic, p::T, t::T, x::T)::T where {T<:Float64}
- Derivative of Water Partitioning Function with respect to X_CO2
- This function is used within the
exsolve3
function.
Arguments
p
: pressure (Pa)t
: temperature (K)x
: The previous mole fraction of CO2 (X_CO2)
Chamber.dwater_dx
— Methoddwater_dx(composition::Silicic, p::T, t::T, x::T)::T where {T<:Float64}
- Function for the derivative of water with reference to X_CO2
- This function is used within the
exsolve3
function.
Arguments
p
: pressure (Pa)t
: temperature (K)x
: The previous mole fraction of CO2 (X_CO2)
Chamber.eos_g
— Methodeos_g(P::Float64, T::Float64)::EosG{Float64}
parametrization of redlich kwong taken from Huber et al. 2010
Arguments
P
: Pressure (Pa)T
: Temperature (K)
Returns
An instance of the EosG
struct, containing the following fields:
rho_g
: density of the gas (kg/m³).drho_g_dP
: partial derivative of the density with respect to pressure (kg/m³/Pa).drho_g_dT
: partial derivative of the density with respect to temperature (kg/m³/K).
Chamber.eos_g_rho_g
— Methodeos_g_rho_g(P::Float64, T::Float64)::Float64
Spetialized version of eos_g
that computes rho_g
only.
Arguments
P
: Pressure (Pa)T
: Temperature (K)
Returns
rho_g
: density of the gas in kg/m³.
Chamber.exsolve
— Methodexsolve(composition::Mafic, P::Float64, T::Float64, X_co2::Float64)::NamedTuple{(:meq, :dmeqdP, :dmeqdT, :dmeqdXco2, :C_co2, :dC_co2dP, :dC_co2dT, :dC_co2dXco2), NTuple{8, Float64}}
This script uses Liu et al. (2006) to calculate the solubility of water
Arguments
composition
:Mafic()
P
: Pressure (Pa)T
: Temperature (K)X_co2
: mole fraction of CO2 in gas.
Chamber.exsolve
— Methodexsolve(composition::Silicic, P::Float64, T::Float64, X_co2::Float64)::NamedTuple{(:meq, :dmeqdP, :dmeqdT, :dmeqdXco2, :C_co2, :dC_co2dP, :dC_co2dT, :dC_co2dXco2), NTuple{8, Float64}}
This script uses Liu et al. (2006) to calculate the solubility of water
Arguments
composition
:Silicic()
P
: Pressure (Pa)T
: Temperature (K)X_co2
: mole fraction of CO2 in gas.
Chamber.exsolve3
— Methodexsolve3(composition::Mafic, P::Float64, T::Float64, m_eq::Float64)::NamedTuple{(:C_co2, :X_co2), NTuple{2, Float64}}
Takes pressure, temperature, and amount of water to solve for the concentration of CO2 and X_CO2 using a Newton Raphson scheme
Arguments
composition
:Mafic()
P
: pressure (Pa)T
: temperature (K)m_eq
: amount of water
Returns
A NamedTuple with the following fields:
C_co2
: The concentration of CO2 in the melt.X_co2
: The mole fraction of CO2 in the gas.
Chamber.exsolve3
— Methodexsolve3(composition::Silicic, P::Float64, T::Float64, m_eq::Float64)::NamedTuple{(:C_co2, :X_co2), NTuple{2, Float64}}
Takes pressure, temperature, and amount of water to solve for the concentration of CO2 and X_CO2 using a Newton Raphson scheme
Arguments
composition
:Silicic()
P
: pressure (Pa)T
: temperature (K)m_eq
: amount of water
Returns
A NamedTuple with the following fields:
C_co2
: The concentration of CO2 in the melt.X_co2
: The mole fraction of CO2 in the gas.
Chamber.exsolve_meq
— Methodexsolve_meq(composition::Mafic, P::Float64, T::Float64, X_co2::Float64)::Float64
Spetialized version of exsolve
that computes meq
only.
Arguments
composition
:Mafic()
P
: Pressure (Pa)T
: Temperature (K)X_co2
: mole fraction of CO2 in gas.
Returns
meq
: amount of water
Chamber.exsolve_meq
— Methodexsolve_meq(composition::Silicic, P::Float64, T::Float64, X_co2::Float64)::Float64
Spetialized version of exsolve
that computes meq
only.
Arguments
composition
:Silicic()
P
: Pressure (Pa)T
: Temperature (K)X_co2
: mole fraction of CO2 in gas.
Returns
meq
: amount of water
Chamber.find_liq
— Methodfind_liq(composition::Mafic, water::Float64, co2::Float64, P::Float64, ini_eps_x::Float64)::Float64
Given the weight fractions of H2O and CO2 in magma, and the pressure, finds the temperature at which the magma will crystallize and reach the desired volume fraction of crystals. The function uses the parameters of the melting curve and an error function to estimate the liquidus temperature.
Arguments
composition
:Mafic()
water
: Weight fration of the H2O in magma.co2
: Weight fration of the CO2 in magma.P
: Pressure (Pa)ini_eps_x
: The starting volumn fraction of crystal.
Returns
Tl
: the estimated liquidus temperature (K).
Chamber.find_liq
— Methodfind_liq(composition::Silicic, water::Float64, co2::Float64, P::Float64, ini_eps_x::Float64)::Float64
Given the weight fractions of H2O and CO2 in magma, and the pressure, finds the temperature at which the magma will crystallize and reach the desired volume fraction of crystals. The function uses the parameters of the melting curve and an error function to estimate the liquidus temperature.
Arguments
composition
:Silicic()
water
: Weight fration of the H2O in magma.co2
: Weight fration of the CO2 in magma.P
: Pressure (Pa)ini_eps_x
: The starting volumn fraction of crystal.
Returns
Tl
: the estimated liquidus temperature (K).
Chamber.fun
— Methodfun(x::Float64, P::Float64, T::Float64, eps_g0::Float64, eps_x0::Float64, V::Float64, rho_m::Float64, rho_g::Float64, mm_co2::Float64, mm_h2o::Float64, M_co2::Float64)::Float64
- This function is used within the
solve_X_co2
function to solve X_co20.
Chamber.gas_heat_capacity
— Methodgas_heat_capacity(X_co2::Float64)
Arguments
X_co2
: mole fraction of CO2 in the gas.
Chamber.get_eps_g
— Methodget_eps_g(composition::Union{Silicic,Mafic}, eps_g_prev::Float64, X_co20::Float64, P::Float64, T::Float64, eps_x0::Float64, V::Float64, rho_m::Float64, rho_g::Float64, M_h2o::Float64, M_co2::Float64)::NamedTuple{(:eps_g0, :mco2_diss),NTuple{2,Float64}}
- This function is used within the
IC_Finder
function to get epsg0 and mco2diss.
Returns
A NamedTuple with the following fields:
eps_g0
mco2_diss
Chamber.get_phase
— Methodget_phase(composition::Union{Silicic,Mafic}, P::Float64, T::Float64, V::Float64, rho_m::Float64, M_h2o::Float64, M_co2::Float64, eps_x0::Float64)::NamedTuple{(:phase, :m_co2_melt),NTuple{2,Float64}}
- This function is used within the
IC_Finder
function to determine the initial phase.
Returns
A NamedTuple with the following fields:
phase
: 2 or 3m_co2_melt
Chamber.heat_conduction_chamberCH
— Methodheat_conduction_chamberCH(maxn::Int64, a::Float64, c::Float64, dr::Float64, kappa::Float64, rho::Float64, cp::Float64, Tb::Float64, param_sv::ParamSaved{Float64})::Float64
The function calculates the heat conduction within a magma chamber by solving the heat equation in cylindrical coordinates.
Arguments
maxn
: number of termsa
: radius of magma chamber (m)c
: radius of outer shell (m)dr
: grid spacing for calculating the heat transferkappa
: thermal diffusivity of host rocks (m^2/s)rho
: density of the magma (kg/m³)cp
: specific heat of magma (J/(kg K))Tb
: boundary temperature of the outer shell (K)param_sv
: a structParamSaved
stores the previous results of the function
Returns
Q
: heat flow rate through the magma chamber (J/s)
Chamber.heat_conduction_chamber_profileCH
— Methodheat_conduction_chamber_profileCH(maxn::Int64, a::Float64, c::Float64, r::Float64, kappa::Float64, Tb::Float64, param_sv::ParamSaved{Float64})::Float64
Calculates the heat transfer in a magma chamber based on the thermal diffusivity of the host rocks and the temperature profile of the outer shell.
Arguments
maxn
: number of termsa
: radius of magma chamber (m)c
: radius of outer shell (m)r
: grid spacing of calculate the heat transferkappa
: thermal diffusivity of host rocksTb
: boundary temperature of the outer shell (K)param_sv
: a structParamSaved
stores the previous results of the function
Returns
Trt
Chamber.ic_phase_conversion
— Methodic_phase_conversion(phase_here::T, composition::Union{Silicic,Mafic}, M_h2o::T, M_co2::T, M_tot::T, P::T, Temp::T, V::T, rho_m::T, param_IC::ParamICFinder{T})::NamedTuple{(:eps_g_temp, :X_co2_temp, :C_co2_temp, :phase),NTuple{4,T}} where {T<:Float64}
This function is used to handle phase conversion by iteratively modifying the max_count
or Tol
parameters of function IC_Finder
until the correct phase is obtained.
- This function is used within the
affect!
function.
Chamber.mco2_dissolved_sat
— Methodmco2_dissolved_sat(X::Float64, P::Float64, T::Float64)::Float64
- This function is used within the
solve_X_co2
function
Arguments
X
: mole fraction of CO2 in gasP
: Pressure (Pa)T
: Temperature (K)
Chamber.meq_water
— Methodmeq_water(composition::Mafic, X::Float64, P::Float64, T::Float64)::Float64
- This function is used within the
get_eps_g
function
Arguments
X
: mole fration of H2O in gasP
: Pressure (Pa)T
: Temperature (K)
Chamber.meq_water
— Methodmeq_water(composition::Silicic, X::Float64, P::Float64, T::Float64)::Float64
- This function is used within the
get_eps_g
function
Arguments
X
: mole fration of H2O in gasP
: Pressure (Pa)T
: Temperature (K)
Chamber.odeChamber
— MethododeChamber(du::Vector{Float64}, u::Vector{Float64}, params::Tuple{Param{Float64}, ParamSaved{Float64}, SW{Int8}}, t::Float64)
- Define the ODE equation.
- Solve the model for eruption frequency of upper crustal magma chambers using an ODE solver.
Arguments
du
: An array that stores the output of the ODE solver, i.e., the values of the derivatives of the solutionu
with respect to timet
.u
: An array that stores the values of the solution at each time stept
.p
: A tuple that stores the model parameters and some saved variables, which are described in more detail below.t
: The time points corresponding to the saved values of the ODE solution.
The arguments du
, u
, p
, and t
are from the DifferentialEquations.jl package. These argument formats are specific to the DifferentialEquations.jl package.
Parameters
param
: A custom parameter containing physical constants and other model parameters.param_saved_var
: A custom parameter used to store values from the previous time step.sw
: A custom parameter used to control simulation behavior.
Returns
The function modifies du
in place to store the values of the derivatives of the solution u
with respect to time t
.
Chamber.parameters_melting_curve
— Methodparameters_melting_curve(composition::Mafic, mH2O::Float64, mCO2::Float64, P::Float64)::NamedTuple{(:a, :dadx, :dady, :dadz, :b, :dbdx, :dbdy, :dbdz), NTuple{8, Float64}}
- This function is used within the
find_liq
,crystal_fraction
,crystal_fraction_eps_x
function.
Arguments
composition
:Mafic()
mH2O
: Weight fration of the H2O in magma.mCO2
: Weight fration of the CO2 in magma.P
: Pressure (Pa)
Chamber.parameters_melting_curve
— Methodparameters_melting_curve(composition::Silicic, mH2O::Float64, mCO2::Float64, P::Float64)::NamedTuple{(:a, :dadx, :dady, :dadz, :b, :dbdx, :dbdy, :dbdz, :c, :dcdx, :dcdy, :dcdz), NTuple{12, Float64}}
- This function is used within the
find_liq
,crystal_fraction
,crystal_fraction_eps_x
function.
Arguments
composition
:Silicic()
mH2O
: Weight fration of the H2O in magma.mCO2
: Weight fration of the CO2 in magma.P
: Pressure (Pa)
Chamber.rc_f
— Methodrc_f(;rho_x::T, eps_x::T, c_x::T, rho_m::T, eps_m::T, c_m::T, rho_g::T, eps_g::T, c_g::T)::T where {T<:Float64}
Computing the product of density and specific heat
- This function is used within the
build_rho_rc
function.
Returns
rc
: the product of density and specific heat
Chamber.record_erupt_end
— Methodrecord_erupt_end(time::T, erupt_saved::EruptSaved{T}, param::Param{T}) where {T<:Float64}
Record duration, mass and volume of eruptions to EruptSaved
- This function is used within the
affect!
function.
Chamber.record_erupt_start
— Methodrecord_erupt_start(time::T, eps_g::T, eps_x::T, rho_m::T, rho_x::T, rho_g::T, erupt_saved::EruptSaved{T}) where {T<:Float64}
Record time, density of eruptions to EruptSaved
- This function is used within the
affect!
function.
Chamber.rho_0_f
— Methodrho_0_f(eps_g0::Float64, eps_x0::Float64, rho_g0::Float64, rho_m0::Float64, rho_x0::Float64)::Float64
- This function is used within the
chamber
function.
Returns
rho_0
: initial bulk density (kg/m³)
Chamber.rho_f
— Methodrho_f(;eps_m::T, eps_g::T, eps_x::T, rho_m::T, rho_g::T, rho_x::T)::T where {T<:Float64}
- This function is used within the
build_rho_rc
function.
Returns
rho
: bulk density (kg/m³)
Chamber.solve_NR
— Methodsolve_NR(f, f_prime, errorTol::T, count_max::T, Xc_initial::T)::T where {T<:Float64}
- Finding mole fraction of CO2 in gas (X_CO2).
- This function is used within the
exsolve3
function.
Arguments
f
: Water Paritioning Functionf_prime
: Function for the derivative of water with reference to X_CO2errorTol
: error tolerance, the default value is 1e-10count_max
: Maximum loop count, the default value is 1e2Xc_initial
: initial guess of X_CO2, the dedault value is 1e-2
Returns
X_CO2
: mole fraction of CO2 in gas
Chamber.solve_X_co2
— Methodsolve_X_co2(composition::Mafic, eps_g0::Float64, X_co2_prev::Float64, P::Float64, T::Float64, eps_x0::Float64, V::Float64, rho_m::Float64, rho_g::Float64, M_co2::Float64, Tol::Float64)::Float64
- This function is used within the
IC_Finder
function to solve X_co20.
Returns
X_co20
Chamber.solve_X_co2
— Methodsolve_X_co2(composition::Silicic, eps_g0::Float64, X_co2_prev::Float64, P::Float64, T::Float64, eps_x0::Float64, V::Float64, rho_m::Float64, rho_g::Float64, M_co2::Float64, Tol::Float64)::Float64
- This function is used within the
IC_Finder
function to solve X_co20.
Returns
X_co20
Chamber.stopChamber_MT
— MethodstopChamber_MT(out, u::Vector{Float64}, t::Float64, int, sw::SW{Int8}, param::Param{Float64})
Define the stopping criteria for an ODE solver that simulates a magma chamber.
Arguments
out
: An array where the function should save the condition value at the right index. The maximum index ofout
should be specified in thelen
property ofcallback
, which allows for a chain oflen
events, triggering theith
event whenout[i] = 0
. The function returns the value ofout[8]
as the last condition. Checking Event Handling and Callback Functions page ofDifferentialEquations.jl
for more details.u
: A vector containing the state of the system at timet
.t
: The current time of the ODE solver.int
: The current state of the integrator. It's format is from the DifferentialEquations.jl packagesw
: A custom parameter used to control simulation behavior.param
: A custom parameter containing physical constants and other model parameters.
The arguments out
, u
, t
, and int
are from the DifferentialEquations.jl package. These argument formats are specific to the DifferentialEquations.jl package.
Returns
The out
array is modified in-place to contain the condition values at the current state of the system. The function computes several quantities based on the current state of the system and the model parameters, such as the crystal fraction, the maximum amount of exsolved fluid, and the amount of CO2 in the melt. It then uses these values to calculate the condition values to be stored in the out
array, as follows:
out[1]
: the crystal fraction.out[2]
: the ratio of crystal fraction to liquid fraction, minus 0.8.out[3]
: if an eruption is not occurring, the pressure difference between the lithostatic pressure and the current pressure, minus the critical pressure difference. Otherwise, negative the critical pressure difference.out[4]
: if an eruption is occurring, the lithostatic pressure minus the current pressure. Otherwise, negative the critical pressure difference.out[5]
: the crystal fraction minus 0.5.out[6]
: the difference between the amount of water in the melt and the maximum amount of exsolved fluid.out[7]
: negative the difference between the initial lithostatic pressure and the current pressure plus the critical pressure difference.out[8]
: the difference between the amount of CO2 in the melt and the saturation concentration.
Note that the out
and u
arguments are in the format expected by the DifferentialEquations.jl package, and the function is intended to be used as a condition for a callback function.
Chamber.water
— Methodwater(composition::Silicic, p::T, t::T, x::T, c::T)::T where {T<:Float64}
- Water Paritioning Function
- This function is used within the
exsolve3
function.
Arguments
p
: pressure (Pa)t
: temperature (K)x
: The previous mole fraction of CO2 (X_CO2)c
: amount of water
Chamber.water
— Methodwater(composition::Silicic, p::T, t::T, x::T, c::T)::T where {T<:Float64}
- Water Paritioning Function
- This function is used within the
exsolve3
function.
Arguments
p
: pressure (Pa)t
: temperature (K)x
: The previous mole fraction of CO2 (X_CO2)c
: amount of water
Chamber.OdeSetting
— TypeSettings for ODE Solver Calculation
Chamber.Param
— TypeParameters
Chamber.ParamICFinder
— TypeSettings for IC_Finder
Chamber.RheolComposition
— TypeRheology of the crust
Chamber.SW
— Typeeruption/coolingmodule/viscousrelaxation control
Index
Chamber.OdeSetting
Chamber.Param
Chamber.ParamICFinder
Chamber.RheolComposition
Chamber.SW
Chamber.IC_Finder
Chamber.IC_Finder
Chamber.a1x_f
Chamber.affect!
Chamber.boundary_conditions_new
Chamber.build_co2
Chamber.build_mdot_in
Chamber.build_meq
Chamber.build_meq
Chamber.build_rho_rc
Chamber.chamber
Chamber.chamber
Chamber.check_for_duplicates
Chamber.compute_dXdP_dXdT
Chamber.crystal_fraction
Chamber.crystal_fraction
Chamber.crystal_fraction_eps_x
Chamber.crystal_fraction_eps_x
Chamber.dM_X_t_dt_f
Chamber.dX_dxdydz
Chamber.drc_dX_f
Chamber.drho_dX_f
Chamber.dwater_dx
Chamber.dwater_dx
Chamber.eos_g
Chamber.eos_g_rho_g
Chamber.exsolve
Chamber.exsolve
Chamber.exsolve3
Chamber.exsolve3
Chamber.exsolve_meq
Chamber.exsolve_meq
Chamber.find_liq
Chamber.find_liq
Chamber.fun
Chamber.gas_heat_capacity
Chamber.get_eps_g
Chamber.get_phase
Chamber.heat_conduction_chamberCH
Chamber.heat_conduction_chamber_profileCH
Chamber.ic_phase_conversion
Chamber.mco2_dissolved_sat
Chamber.meq_water
Chamber.meq_water
Chamber.odeChamber
Chamber.parameters_melting_curve
Chamber.parameters_melting_curve
Chamber.rc_f
Chamber.record_erupt_end
Chamber.record_erupt_start
Chamber.rho_0_f
Chamber.rho_f
Chamber.solve_NR
Chamber.solve_X_co2
Chamber.solve_X_co2
Chamber.stopChamber_MT
Chamber.water
Chamber.water