m_parameters.py

This module contains the most important computation parameters. It is divided into several categories:

Running parameters

1. Output files settings

2. Reloading options

3. Termination criteria

Geometry and material composition

4. Mesh geometry

5. Material composition

Physics

6. Stokes problem settings

7. Heat transport and melting settings

8. Rheology

m_parameters_docs.BC_Stokes_problem

Var:

Boundary conditions for the Stokes problem m_equations.Equations.equation_Stokes() in the following order

• top boundary

  • ["free_slip/no_slip/free surface"],

  • ["velocity_x/velocity_y", value] for specifying one component of the velocity,

  • or ["velocity", value_x, value_y] for specifying both components of the velocity

• bottom boundary

  • ["free_slip/no_slip/free surface"],

  • ["velocity_x/velocity_y", value] for specifying one component of the velocity,

  • or ["velocity", value_x, value_y] for specifying both components of the velocity

• left boundary

  • ["free_slip/no_slip"],

  • ["velocity_x/velocity_y", value] for specifying one component of the velocity,

  • or ["velocity", value_x, value_y] for specifying both components of the velocity

• right boundary

  • ["free_slip/no_slip"],

  • ["velocity_x/velocity_y", value] for specifying one component of the velocity,

  • or ["velocity", value_x, value_y] for specifying both components of the velocity

Vartype:

• list of four lists (one for each boundary)

• unit m s\(^{-1}\)

Hint

For example, for a tectonic simulation with free surface at the top boundary, ouflow through the side boundaries by 10 km/Myr and strictly vertical inflow through the bottom boundary at the same rate, the boundary condition will look as follows:

BC_Stokes_problem = [["free_surface"],          # top boundary       (1)
               ["velocity", 0.0, +10e3/Myr],    # bottom boundary    (2)
               ["velocity_x", -10e3/Myr],       # left boundary      (3)
               ["velocity_x", +10e3/Myr]]       # right boundary     (4)

m_parameters_docs.BC_heat_transfer

Var:

Boundary conditions for heat transfer equation in the following order

• top boundary ["temp/heat_flux/radiation", value]

Note

If "radiation" boundary condition is selected, the value serves as an initial estimate of the top boudnary temperature for the nonlinear solver.

• bottom boundary ["temp/heat_flux", value]

• left boundary ["temp/heat_flux", value]

• right boundary ["temp/heat_flux", value]

Vartype:

list of four lists (one for each boundary)

• unit K for temperature

• unit W m\(^{-2}\) for the heat flux

Hint

For example, in order to prescribe temperature 90 K at the top boundary, 265 K at the bottom boundary and insulating side boundaries, the boundary condition will look as follows:

BC_heat_transfer = [["temp", 90.0],    # top boundary    (1)
                  ["temp", 265.0],     # bottom boundary (2)
                  ["heat_flux", 0.0],  # left boundary   (3)
                  ["heat_flux", 0.0]]  # right boundary  (4)

m_parameters_docs.DAL_factor

Var:

Strength of the phase transition at the ice-water boundary due to the DAL effect

Vartype:

float, unit W m\(^{-3}\)

m_parameters_docs.G_ice

Var:

Shear modulus of ice.

Vartype:

float, unit Pa

m_parameters_docs.Lt

Var:

Latent heat of the material.

Vartype:

float, unit J kg\(^{-1}\)

m_parameters_docs.Paraview_Output

List of strings indicating which properties will be saved, e.g.:

Paraview_Output = ["viscosity", "temperature", "plastic_strain"]

The following are defined:

String	Description
temperature	\(T\) (K)
velocity	\(v\) (m s \(^{-1}\))
pressure	\(p\) (Pa)
viscosity	\(\eta\) (Pa s)
viscosity_log	log\(_{10}(\eta)\) (-)
eta_v	\(\eta_{ductile}\) (Pa s)
composition_0	1st material in materials[]
composition_1	2nd material in materials[]
composition_2	3rd material in materials[]
tracers	number of tracers in the cell
tidal_heating	\(H_{tidal}\) (W m \(^{-3}\))
melt_fraction	\(\chi_m\) (%)
mesh_displacement	\(\boldsymbol{u}_{mesh}\) (m)
topography_bottom	\(h_{bot}\) (-)
topography_top	\(h_{top}\) (-)
iteration_error
plastic_strain	\(\varepsilon_p\) (-)
yield_stress	\(\sigma_Y\) (Pa)
stress_dev_inv	\(\sigma_{II}\) (Pa)
strain_rate_inv	\(\dot{\varepsilon}_{II}\) (s \(^{-1}\))
yield_function	1 if yielding, 0 otherwise
density	\(\rho\) (kg m \(^{-3}\))
z_function	\(Z\), visco-elastic parameter
shear_modulus	\(G\) (Pa)
cohesion	\(C\) (Pa)

Tip

A function can be added by passing it to __init__ in m_filenames.py. Follow the procedure for already existing functions.

m_parameters_docs.Paraview_Output_Ini

List of strings indicating which properties will be saved at initial condition.

The following are defined:

Input	Note
temperature	\(T\) (K)
topography_top	\(h_{top}\) (-)
topography_bottom	\(h_{bot}\) (-)
tracers	number of tracers in the cell
density	\(\rho\) (kg m \(^{-3}\))
melt_fraction	\(\chi_m\) (%)

m_parameters_docs.Picard_iter_error

Var:

Minimum relative error in velocity field sufficient to exit the Stokes solver Picard iterations, see m_equations.Equations.solve_Stokes_problem().

Vartype:

float

m_parameters_docs.Picard_iter_max

Var:

Maximum number of Stokes solver Picard iterations in case of nonlinear rheology (stress-dependent viscosity, plasticity), see m_equations.Equations.solve_Stokes_problem().

Vartype:

integer

m_parameters_docs.Q_activ

Var:

Activation energy for temp-dep viscosity type. See m_rheology.eta_ductile().

Vartype:

float, unit J mol\(^{-1}\)

m_parameters_docs.T_melt

Var:

Melting temperature of the solid.

Vartype:

float, unit K

m_parameters_docs.Tracers_Output

List of strings indicating which properties carried by the tracers will be saved, e.g.:

Tracers_Output = ["rank", "dev_stress_xx", "melt_fraction"]

Properties to save on tracers

String	Description
rank	the process the tracer belongs to
dev_stress_xx	xx component of the deviatoric stress
dev_stress_xz	xz component of the deviatoric stress
plastic_strain	amount of plastic strain
original_y	original z-coordinate of the tracer
composition_n	concentration of material n (n being a number)
melt_fraction	amount of partial melt on the tracer
origin	0 if the tracer is original, 1 if added later
ID	unique ID of the tracer

m_parameters_docs.cfl

Var:

The CFL parameter, determining the length of the adaptive time step, see m_timestep.time_step().

Vartype:

float

m_parameters_docs.cohesion_strong

Var:

Cohesion of an undamaged material. See m_rheology.cohesion().

Vartype:

float, unit Pa

m_parameters_docs.cohesion_weak

Var:

Cohesion of a fully damaged material. See m_rheology.cohesion().

Vartype:

float, unit Pa

m_parameters_docs.cos_perturbation

Var:

Whether to perturb the temperature initial condition by a function. See m_equations.Equations.thermal_perturbation().

Vartype:

boolean

m_parameters_docs.dT_melt_tresh

Var:

Threshold for adding or removing tracers when melting.

Vartype:

float, unit K

m_parameters_docs.d_grain

Var:

Constant grain size. See m_rheology.eta_ductile().

Vartype:

float, unit m

m_parameters_docs.elasticity

Var:

Whether the material behaves as a Maxwell visco-elastic medium.

Vartype:

boolean

m_parameters_docs.empty_cells_allowed

Var:

Whether to allow elements without tracers. If True, tracers will not be added in m_tracers.Tracers.introduce_tracers() nor m_tracers.Tracers.add_tracers().

Vartype:

boolean

Warning

Use True only when material composition is the only tracer-requiring feature.

m_parameters_docs.empty_cells_composition

Var:

If empty_cells_allowed = True, this determines the composition of empty cells, see m_interpolation.composition_interpolation().

Vartype:

integer

m_parameters_docs.empty_cells_region

Var:

Specifies geometry of the region where the empty cells will be given empty_cells_composition even without tracers inside, see m_tracers.Tracers.introduce_tracers().

Vartype:

list of lists consisting of the region shape and its dimensions. For example

empty_cells_region = [["rectangle",  0, length, 20e3, height]]

Note

Only "rectangle" type is implemented at this point.

m_parameters_docs.eps_strong

Var:

Value of the plastic strain at which the material starts to accumulate damage. See m_rheology.cohesion().

Vartype:

float, unit -

m_parameters_docs.eps_weak

Var:

Critical value of the plastic strain beyond which the material is considered as fully damaged. See m_rheology.cohesion().

Vartype:

float, unit -

m_parameters_docs.error_type

Determines the way of calculating the iteration error.

Var:

• "integrated" computes the error as \(\, \int_\Omega |\boldsymbol{v} - \boldsymbol{v_k}|/|\boldsymbol{v}| \textrm{d}\Omega\)

• "maximum" computes the error as \(\, \textrm{max}(\boldsymbol{v} - \boldsymbol{v_k}|/| \boldsymbol{v}|)\)

Vartype:

string

m_parameters_docs.eta_0

Var:

• Constant viscosity for constant viscosity type.

• Viscosity prefactor for temp-dep viscosity type. See m_rheology.eta_ductile().

Vartype:

float, unit Pa s

m_parameters_docs.eta_max

Var:

Upper cut-off value for viscoplastic viscosity \(\eta_{vp}\). See m_rheology.eta_ductile().

Vartype:

float, unit Pa s

m_parameters_docs.eta_min

Var:

Lower cut-off value for plastic viscosity \(\eta_{p}\). See m_rheology.eta_eff().

Vartype:

float, unit Pa s

m_parameters_docs.g

Var:

Gravity acceleration.

Vartype:

float, unit m s\(^{-2}\)

m_parameters_docs.h_bot_ini

Var:

Expression defining initial bottom topography. Will be take into account only if initial_topography = True.

Vartype:

Expression

m_parameters_docs.h_top_ini

Var:

Expression defining initial top topography. Will be take into account only if initial_topography = True.

Vartype:

Expression

m_parameters_docs.healing

Var:

Whether the accumulated plastic strain decays in time due to microscopic healing processes. See m_rheology.plastic_strain_integration().

Vartype:

boolean

m_parameters_docs.healing_timescale

Var:

Characteristic time scale for the microscopic healing processes. See m_rheology.plastic_strain_integration().

Vartype:

float

m_parameters_docs.height

Height of the rectangular mesh.

Vartype:

float

m_parameters_docs.init_cond_profile

Var:

Whether to find an initial conductive profile as an initial condition for the heat transfer problem.

Vartype:

boolean

Note

If False, then a linear temperature profile between the top and bottom temperature will be prescribed.

m_parameters_docs.initial_topography

Var:

Whether there will be an initial topography or not.`.

Vartype:

boolean

m_parameters_docs.int_friction_angle2

Var:

Angle of internal friction in degrees.

Vartype:

float, unit degrees

m_parameters_docs.integration_method

Method for integration of the trajectories of Lagrangian tracers.

Var:

• "RK2" for the 2nd order Runge-Kutta scheme

• "RK4" for the 4nd order Runge-Kutta scheme

Vartype:

string

m_parameters_docs.interface(x)

A function that describes the shape of an interface.

:meta hide-value:.

m_parameters_docs.internal_melting

Var:

Whether generate partial melt if the temperature of the solid reaches \(T_{melt}\).

Vartype:

boolean

m_parameters_docs.length

Length of the rectangular mesh.

Vartype:

float

m_parameters_docs.loading_mesh

Whether to read an external mesh.

Var:

• True reads a mesh from file mesh_name.

• False creates a new mesh based on x_div, z_div and triangle_types

Vartype:

boolean

m_parameters_docs.materials

Var:

Simple material assignment using “interface”, “circle” or “rectangle” geometries.

Vartype:

List of lists consisting of the region shape and its dimensions:

• Circle: ["circle", center_x, center_y, radius]

• Rectangle: ["rectangle", left, right, bottom, top]

• Interface: ["interface", "below/above"]

For example, in order to prescribe an ice shell with 1-km-thick distinct layer at the surface, enter the following:

materials = [["rectangle",  0, length, 0, height],  ["rectangle",  0, length, height-1e3, height]]

Attention

The n-th material overrides the (n-1)th.

Tip

For more complicated geometries, modify introduce_tracers function in the m_tracers.py

m_parameters_docs.mesh_displacement_laplace

Determines whether full Laplace equation or only its z-part will be solved to obtain mesh displacement on the mesh between the top and bottom boundary.

Var:

• "full" distributes the displacement by solving \(\nabla^2 h_2 = 0\). Useful for simulation involving tectonics, mesh is more stable.

• "z_only" distributes the displacement linearly by solving \(\nabla^2_z h_2 = 0\)

Vartype:

string

m_parameters_docs.mesh_name

Mesh to be loaded.

Var:

path to a mesh in xml format, e.g. mesh_25x100km.xml

Vartype:

string

m_parameters_docs.monitor_cache

Var:

If True, returns the number of cache in each time step strored in /home/username/.cache/dijitso/lib/.

Vartype:

boolean

Warning

Number of caches should stop increasing after a few time steps. If not, FEnICS is making a new cache in each time step, which may lead to crash of the code after some time.

m_parameters_docs.name

Var:

Name of the directory with the results. The directory with the results will be named data_name.

Vartype:

string

m_parameters_docs.nonlinear_heat_equation

Var:

Whether the heat transfer m_equations.Equations.equation_heat() should be solved with nonlinear solver.

Vartype:

boolean

m_parameters_docs.output_frequency

Var:

Frequency for the output to ParaView, HDF5 data and mesh, and tracer files. The first time step is always saved, the last one if the simulation succesfully ends through the termination criteria.

Vartype:

list consisting of a string "steps" or "time", and an integer/float For example:

output_frequency = ["steps", 10]

outputs the data every 10 steps, while:

output_frequency = ["time", 100*kyr]

outputs the data every 100 kyr.

m_parameters_docs.perturb_ampl

Var:

Thermal perturbation amplitude \(A\) in m_equations.Equations.thermal_perturbation().

Vartype:

float, unit K

m_parameters_docs.perturb_freq

Var:

Number \(N\)of half-cosine waves in the thermal perturbation formula m_equations.Equations.thermal_perturbation().

Vartype:

float, unit -

m_parameters_docs.phase_transition

Var:

Whether there is a phase transition at the ice-water boundary.

Vartype:

boolean

m_parameters_docs.plasticity

Var:

Whether the material yields when the yield stress is reached.

Vartype:

boolean

m_parameters_docs.protect_directory

• True: Repeated run protects the original data by creating saving into a directory data_name_new instead

• False: Repeated run with the same name overwrites the data

m_parameters_docs.refinement

Var:

Minimum and maximum x- and y-coordinates of a rectangle area of the mesh to be refined, see m_mesh.MeshModule.refine_mesh().

Vartype:

list of 4x n floats, where n is the number of refined areas

Hint

In order to refine the upper half of the mesh:

refinement = [0.0, length, height/2.0, height]

By repeating the sequence, the areas of refinement can be superposed. For example

refinement = [0.0, length, 0, height/2.0, 0.0, length, 0, height/4.0]

gives first-level refinement in the bottom half domain and second-level one in the bottom quarter of the domain.

Warning

For tectonic simulations, the top boundary nodes need to be equidistant because of detecting and smoothing node oscillations, see m_mesh.detect_oscillations().

m_parameters_docs.reload

Var:

Whether or not to load HDF5 data from reload_name/HDF5/data.h5.

Vartype:

boolean

m_parameters_docs.reload_name

Var:

The directory from which the data will be reloaded when reloading_HDF5 = True.

Vartype:

string

m_parameters_docs.reload_step

Var:

Time stamp of the HDF5 file which will be loaded when reload = True.

Vartype:

integer or “last” (reloads the last time step)

m_parameters_docs.restart_time

Whether to reset time when reloading_HDF5 = True.

Var:

• True sets time to 0.

• False continues from the given checkpoint.

Vartype:

boolean

m_parameters_docs.rho_l

Var:

Density of the liquid below the domain (subsurface ocean).

Vartype:

float, unit kg m\(^{-3}\)

m_parameters_docs.rho_m

Var:

Density of the melt produced by the solid.

Vartype:

float, unit kg m\(^{-3}\)

m_parameters_docs.rho_s

Var:

Density of the domain solid.

Vartype:

float, unit kg m\(^{-3}\)

m_parameters_docs.solve_energy_problem

Var:

Whether to solve heat transfer m_equations.Equations.equation_heat().

Vartype:

boolean

m_parameters_docs.stat_output

List of strings indicating which properties will be saved in the text file at each time step.

The following are defined:

Input	Note
nusselt	\(Nu\) (-)
q_top	\(q_{top}\) (W m \(^{-2}\))
vrms	\(v_{rms}\) (m s \(^{-1}\))
avg_h_bot	average bottom topography (m)
h_top_max	maximum top topography (m)
time	wallclock time of the simulation (hours)
timestep	one time step duration (seconds)

m_parameters_docs.stokes_null

Var:

Method for correcting velocity field in case of free surface conditions on both top and bottom boundary.

Vartype:

• "boundary" for subtracting the average \(v_z\)integrated along the top boundary

• "volume" for subtracting the average \(v_z\)integrated over the full domain

m_parameters_docs.termination_condition

Time or step criterion for ending the simulation naturally.

Var:

• ["step", 1000] exits the time loop after the 100th step

• ["time", 100*kyr] exits the time loop as soon as 100 kyr is reached

• ["initial_condition", 100*kyr] exits the time loop after writing the initial condition

m_parameters_docs.tidal_dissipation

Var:

Whether there is a volumetric source term in the heat transfer equation.

Vartype:

boolean

m_parameters_docs.time_step_position

Determines where to compute new time step.

Var:

"stokes" after the Stokes solver or "end" at the end of the time loop

Vartype:

string

m_parameters_docs.time_step_strategy

Specifies the method for computing a new time step, see m_timestep.time_step().

Var:

• "domain" computes a new time step based on specific criteria in every step of the time loop (see m_timestep module)

• "convective" prescribes a convective time step given by \(v_{max}/\Delta x_{min}\)

• "constant" prescribes a constant time step given by the dt_const parameter

Vartype:

string

m_parameters_docs.time_units

The units of time at which the output will be written.

Var:

Use 1.0 for seconds or in case of non-dimensional problems, or yr, kyr or Myr

Vartype:

float

m_parameters_docs.time_units_string

Var:

String representation of time_units, up to the user.

Vartype:

string, e.g. "s", "yr", "kyr" or "Myr"

m_parameters_docs.triangle_types

Method of dividing basic squares into mesh triangle elements.

Var:

"crossed", "left", "right", "right/left" or "left/right"

Vartype:

string

m_parameters_docs.viscosity_type

Var:

Defines which formula for ductile viscosity will be used, see m_rheology.eta_ductile().

Vartype:

string

• "constant" for constant viscosity of value \(\eta_0\).

• "temp-dep" for temperature-dependent viscosity

• "GK_2001" for temperature and stress-dependent viscosity of ice-I

• "composition" for composition-dependent viscosity (used for selected benchmarks)

m_parameters_docs.x_div

Number of nodes in horizontal direction.

Var:

default int(z_div*(length/height)) which keeps acpect ratio of the elements equal to 1

Vartype:

integer

m_parameters_docs.yield_stress_max

Var:

Upper and lower cut-off values for the yield stress, respectively.

Vartype:

float, float

m_parameters_docs.yield_stress_min

Var:

Upper and lower cut-off values for the yield stress, respectively.

Vartype:

float, float

m_parameters_docs.z_div

Number of nodes in vertical direction.

Vartype:

integer