List functions
These are all the functions that are provided for the LaMEM Julia Setup interface
LaMEM.LaMEM_Model.BoundaryConditions
— TypeStructure that contains the LaMEM boundary conditions information.
noslip::Vector{Int64}
: No-slip boundary flag mask (left right front back bottom top)open_top_bound::Int64
: Stress-free (free surface/infinitely fast erosion) top boundary flagtemp_top::Float64
: Constant temperature on the top boundarytemp_bot::Float64
: Constant temperature on the bottom boundaryexx_num_periods::Int64
: number intervals of constant background strain rate (x-axis)exx_time_delims::Vector{Float64}
: time delimiters (one less than number of intervals, not required for one interval)exx_strain_rates::Vector{Float64}
: strain rates for each intervaleyy_num_periods::Int64
eyy_time_delims::Vector{Float64}
eyy_strain_rates::Vector{Float64}
exy_num_periods::Int64
exy_time_delims::Vector{Float64}
exy_strain_rates::Vector{Float64}
exz_num_periods::Int64
exz_time_delims::Vector{Float64}
exz_strain_rates::Vector{Float64}
eyz_num_periods::Int64
eyz_time_delims::Vector{Float64}
eyz_strain_rates::Vector{Float64}
bg_ref_point::Vector{Float64}
: background strain rate reference point (fixed)
LaMEM.LaMEM_Model.Dike
— TypeDefines the properties related to inserting dikes
ID::Int64
: Material phase IDMf::Float64
: value for dike/magma- accommodated extension, between 0 and 1, in the front of the box, for phase dikeMc::Float64
: [optional] value for dike/magma- accommodate extension, between 0 and 1, for dike phase; M is linearly interpolated between Mf & Mc and Mc & Mb, if not set, Mc default is set to -1 so it is not usedy_Mc::Union{Nothing, Float64}
: [optional], location for Mc, must be between front and back boundaries of dike box, if not set, default value to 0.0, but not usedMb::Union{Nothing, Float64}
: value for dike/magma-accommodated extension, between 0 and 1, in the back of the box, for phase dikePhaseID::Union{Nothing, Int64}
: Phase IDPhaseTransID::Union{Nothing, Int64}
: Phase transition ID
LaMEM.LaMEM_Model.FreeSurface
— TypeStructure that contains the LaMEM free surface information.
surf_use::Int64
: Free surface activation flagsurf_corr_phase::Int64
: air phase ratio correction flag (phases in an element that contains are modified based on the surface position)surf_level::Float64
: initial level of the free surfacesurf_air_phase::Int64
: phase ID of sticky air layersurf_max_angle::Float64
: maximum angle with horizon (smoothed if larger)surf_topo_file::String
: initial topography file (redundant)erosion_model::Int64
: erosion model [0-none (default), 1-infinitely fast, 2-prescribed rate with given level]er_num_phases::Int64
: number of erosion phaseser_time_delims::Vector{Float64}
: erosion time delimiters (one less than number)er_rates::Vector{Float64}
: constant erosion rates in different time periodser_levels::Vector{Int64}
: levels above which we apply constant erosion rates in different time periodssediment_model::Int64
: sedimentation model [0-none (dafault), 1-prescribed rate with given level, 2-cont. margin]sed_num_layers::Int64
: number of sediment layerssed_time_delims::Vector{Float64}
: sediment layers time delimiters (one less than number)sed_rates::Vector{Float64}
: sediment rates in different time periodssed_levels::Vector{Float64}
: levels below which we apply constant sediment rates in different time periodssed_phases::Vector{Int64}
: sediment layers phase numbers in different time periodsmarginO::Vector{Float64}
: lateral coordinates of continental margin - originmarginE::Vector{Float64}
: lateral coordinates of continental margin - 2nd pointhUp::Float64
: up dip thickness of sediment cover (onshore)hDown::Float64
: down dip thickness of sediment cover (off shore)dTrans::Float64
: half of transition zone
LaMEM.LaMEM_Model.Grid
— TypeStructure that contains the LaMEM grid information
nmark_x::Int64
: number of markers/element in x-directionnmark_y::Int64
: number of markers/element in y-directionnmark_z::Int64
: number of markers/element in x-directionnel_x::Vector{Int64}
: number of elements in x-directionnel_y::Vector{Int64}
: number of elements in y-directionnel_z::Vector{Int64}
: number of elements in z-directioncoord_x::Vector{Float64}
: coordinates in x-directioncoord_y::Vector{Float64}
: coordinates in y-directioncoord_z::Vector{Float64}
: coordinates in z-directionnseg_x::Int64
: number of segments in x-direction (if we employ variable grid spacing in x-direction)nseg_y::Int64
: number of segments in y-direction (if we employ variable grid spacing in y-direction)nseg_z::Int64
: number of segments in z-direction (if we employ variable grid spacing in z-direction)bias_x::Vector{Float64}
: bias in x-direction (if we employ variable grid spacing in x-direction)bias_y::Vector{Float64}
: bias in y-direction (if we employ variable grid spacing in y-direction)bias_z::Vector{Float64}
: bias in z-direction (if we employ variable grid spacing in z-direction)Grid::GeophysicalModelGenerator.LaMEM_grid
: Contains the LaMEM Grid objectPhases::Array{Int32}
: Phases; 3D phase informationTemp::Array{Float64}
: Temp; 3D phase information
Example 1
julia> d=LaMEM.Grid(coord_x=[0.0, 0.7, 0.8, 1.0], bias_x=[0.3,1.0,3.0], nel_x=[10,4,2])
LaMEM grid with 1D refinement:
nel : ([10, 4, 2], [16], [16])
marker/cell : (3, 3, 3)
x ϵ [0.0, 0.7, 0.8, 1.0], bias=[0.3, 1.0, 3.0], nseg=3, Δmin=0.025000000000000022, Δmax=0.1499999999999999
y ϵ [-10.0 : 0.0]
z ϵ [-10.0 : 0.0]
Example 2
julia> d=LaMEM.Grid(nel=(10,20))
LaMEM grid with constant Δ:
nel : ([10], [1], [20])
marker/cell : (3, 3, 3)
x ϵ [-10.0 : 10.0]
y ϵ [-10.0 : 0.0]
z ϵ [-10.0 : 0.0]
LaMEM.LaMEM_Model.Materials
— TypeStructure that contains the material properties in the current simulation
Phases::Vector{Phase}
: Different Materials implementedSofteningLaws::Vector{Softening}
: Softening laws implementedPhaseTransitions::Vector{PhaseTransition}
: Internal Phase Transitions (that change the ID of markers) implementedDikes::Vector{Dike}
: Dikes implemented (mostly for MOR simulations)
LaMEM.LaMEM_Model.Model
— Type; Model
Structure that holds all the information to create a LaMEM input file
Scaling::Scaling
: Scaling parametersGrid::Grid
: LaMEM GridTime::Any
: Time optionsFreeSurface::Any
: Free surface optionsBoundaryConditions::Any
: Boundary conditionsSolutionParams::Any
: Global solution parametersSolver::Any
: Solver options and optional PETSc optionsModelSetup::Any
: Model setupOutput::Any
: Output optionsMaterials::Any
: Material parameters for each of the phases
LaMEM.LaMEM_Model.Model
— MethodModel(args...)
Allow to define a model setup by specifying some of the basic objects
Example
julia> d = Model(Grid(nel=(10,1,20)), Scaling(NO_units()))
LaMEM Model setup
|
|-- Scaling : GeoParams.Units.GeoUnits{GeoParams.Units.NONE}
|-- Grid : nel=(10, 1, 20); xϵ(-10.0, 10.0), yϵ(-10.0, 0.0), zϵ(-10.0, 0.0)
|-- Time : nstep_max=50; nstep_out=1; time_end=1.0; dt=0.05
|-- Boundary conditions : noslip=[0, 0, 0, 0, 0, 0]
|-- Solution parameters :
|-- Solver options : direct solver; superlu_dist; penalty term=10000.0
|-- Model setup options : Type=geom;
|-- Output options : filename=output; pvd=1; avd=0; surf=0
|-- Materials : 1 phases;
LaMEM.LaMEM_Model.Model
— MethodModel(;
Scaling=Scaling(GEO_units()),
Grid=Grid(),
Time=Time(),
FreeSurface=FreeSurface(),
BoundaryConditions=BoundaryConditions(),
SolutionParams=SolutionParams(),
Solver=Solver(),
ModelSetup=ModelSetup(),
Output=Output(),
Materials=Materials()
)
Creates a LaMEM Model setup.
Scaling::Scaling
Grid::Grid
Time::Any
FreeSurface::Any
BoundaryConditions::Any
SolutionParams::Any
Solver::Any
ModelSetup::Any
Output::Any
Materials::Any
LaMEM.LaMEM_Model.ModelSetup
— TypeStructure that contains the LaMEM Model Setup and Advection options
msetup::String
: Setup type - can begeom
(phases are assigned from geometric primitives),files
(from julia input),polygons
(from geomIO input, which requirespoly_file
to be specified)rand_noise::Int64
: add random noise to the particle locationrand_noiseGP::Int64
: random noise flag, subsequently applied to geometric primitivesbg_phase::Int64
: background phase IDsave_mark::Int64
: save marker to disk flagmark_load_file::String
: marker input file (extension is .xxxxxxxx.dat), if usingmsetup
=files
mark_save_file::String
: marker output file (extension is .xxxxxxxx.dat)poly_file::String
: polygon geometry file (redundant), if usingmsetup
=polygons
temp_file::String
: initial temperature file (redundant), if not set on markersadvect::String
: advection scheme; options=none
(no advection);basic
(Euler classical implementation [default]);Euler
(Euler explicit in time);rk2
(Runge-Kutta 2nd order in space)interp::String
: velocity interpolation scheme; options =stag
(trilinear interpolation from FDSTAG points),minmod
( MINMOD interpolation to nodes, trilinear interpolation to markers + correction),stagp
( STAG_P empirical approach by T. Gerya)stagp_a::Float64
: STAG_P velocity interpolation parametermark_ctrl::String
: marker control type; options aresubgrid
(default; marker control enforced over fine scale grid),none
(none),basic
(AVD for cells + corner insertion), andavd
(pure AVD for all control volumes)nmark_lim::Vector{Int64}
: min/max number per cell (marker control)nmark_avd::Vector{Int64}
: x-y-z AVD refinement factors (avd marker control)nmark_sub::Int64
: max number of same phase markers per subcell (subgrid marker control)geom_primitives::Vector
: Different geometric primitives that can be selected if wemsetup
=
geom; see
geom_Sphere`
LaMEM.LaMEM_Model.Output
— TypeStructure that contains the LaMEM output options
out_file_name::Any
: output file nameout_dir::Any
: output directoryparam_file_name::Any
: parameter filenamewrite_VTK_setup::Any
: write VTK initial model setupout_pvd::Any
: activate writing .pvd fileout_phase::Any
: dominant phaseout_density::Any
: densityout_visc_total::Any
: total (viscoelastoplastic) viscosityout_visc_creep::Any
: creep viscosityout_velocity::Any
: velocityout_pressure::Any
: (dynamic) pressureout_tot_press::Any
: total pressureout_eff_press::Any
: effective pressureout_over_press::Any
out_litho_press::Any
out_pore_press::Any
out_temperature::Any
: temperatureout_dev_stress::Any
: deviatoric strain rate tensorout_j2_dev_stress::Any
: second invariant of deviatoric stress tensorout_strain_rate::Any
: deviatoric strain rate tensorout_j2_strain_rate::Any
: second invariant of strain rate tensorout_shmax::Any
out_ehmax::Any
out_yield::Any
out_rel_dif_rate::Any
: relative proportion of diffusion creep strainrateout_rel_dis_rate::Any
: relative proportion of dislocation creep strainrateout_rel_prl_rate::Any
: relative proportion of peierls creep strainrateout_rel_pl_rate::Any
: relative proportion of plastic strainrateout_plast_strain::Any
: accumulated plastic strainout_plast_dissip::Any
: plastic dissipationout_tot_displ::Any
out_moment_res::Any
out_cont_res::Any
out_energ_res::Any
out_melt_fraction::Any
out_fluid_density::Any
out_conductivity::Any
out_vel_gr_tensor::Any
out_surf::Any
: activate surface outputout_surf_pvd::Any
: activate writing .pvd fileout_surf_velocity::Any
: surface velocityout_surf_topography::Any
: surface topographyout_surf_amplitude::Any
: amplitude of topography (=topo-average(topo))out_mark::Any
: activate marker outputout_mark_pvd::Any
: activate writing .pvd fileout_avd::Any
: activate AVD phase outputout_avd_pvd::Any
: activate writing .pvd fileout_avd_ref::Any
: AVD grid refinement factorout_ptr::Any
: activateout_ptr_ID::Any
: ID of the passive tracersout_ptr_phase::Any
: phase of the passive tracersout_ptr_Pressure::Any
: interpolated pressureout_ptr_Temperature::Any
: temperatureout_ptr_MeltFraction::Any
: melt fraction computed using P-T of the markerout_ptr_Active::Any
: option that highlight the marker that are currently activeout_ptr_Grid_Mf::Any
: option that allow to store the melt fraction seen within the cell
LaMEM.LaMEM_Model.Phase
— TypeDefines the material properties for each of the phases
ID::Union{Nothing, Int64}
: Material phase IDName::Union{Nothing, String}
: Description of the phaserho::Union{Nothing, Float64}
: Density [kg/m^3]eta::Union{Nothing, Float64}
: Linear viscosity [Pas]visID::Union{Nothing, Int64}
: material ID for phase visualization (default is ID)diff_prof::Union{Nothing, String}
: DIFFUSION creep profile; example: "DryOlivinediffcreep-HirthKohlstedt_2003"disl_prof::Union{Nothing, String}
: DISLOCATION creep profile; example: "Granite-Tireletal_2008"peir_prof::Union{Nothing, String}
: PEIERLS creep profile; example: "OlivinePeierls-Kameyama1999"rho_n::Union{Nothing, Float64}
: depth-dependent density model parameterrho_c::Union{Nothing, Float64}
: depth-dependent density model parameterbeta::Union{Nothing, Float64}
: pressure-dependent density model parameterG::Union{Nothing, Float64}
: shear modulusKb::Union{Nothing, Float64}
: bulk modulusE::Union{Nothing, Float64}
: Young's modulusnu::Union{Nothing, Float64}
: Poisson's ratioKp::Union{Nothing, Float64}
: pressure dependence parameterBd::Union{Nothing, Float64}
: DIFFUSION creep pre-exponential constantEd::Union{Nothing, Float64}
: activation energyVd::Union{Nothing, Float64}
: activation volumeeta0::Union{Nothing, Float64}
: POWER LAW reference viscositye0::Union{Nothing, Float64}
: reference strain rateBn::Union{Nothing, Float64}
: DISLOCATION creep pre-exponential constantEn::Union{Nothing, Float64}
: activation energyVn::Union{Nothing, Float64}
: activation volumen::Union{Nothing, Float64}
: power law exponentBp::Union{Nothing, Float64}
: PEIERLS creep pre-exponential constantEp::Union{Nothing, Float64}
: activation energyVp::Union{Nothing, Float64}
: activation volumetaup::Union{Nothing, Float64}
: scaling stressgamma::Union{Nothing, Float64}
: approximation parameterq::Union{Nothing, Float64}
: stress-dependence parametereta_fk::Union{Nothing, Float64}
: reference viscosity for Frank-Kamenetzky viscositygamma_fk::Union{Nothing, Float64}
: gamma parameter for Frank-Kamenetzky viscosityTRef_fk::Union{Nothing, Float64}
: reference Temperature for Frank-Kamenetzky viscosity (if not set it is 0°C)ch::Union{Nothing, Float64}
: cohesionfr::Union{Nothing, Float64}
: friction angleeta_st::Union{Nothing, Float64}
: stabilization viscosity (default is eta_min)rp::Union{Nothing, Float64}
: pore-pressure ratiochSoftID::Union{Nothing, Int64}
: friction softening law IDfrSoftID::Union{Nothing, Int64}
: cohesion softening law IDhealID::Union{Nothing, Int64}
: healing ID, points to healTau in Softeningalpha::Union{Nothing, Float64}
: thermal expansivityCp::Union{Nothing, Float64}
: specific heat (capacity), J⋅K−1⋅kg−1k::Union{Nothing, Float64}
: thermal conductivityA::Union{Nothing, Float64}
: radiogenic heat productionT::Union{Nothing, Float64}
: optional temperature to set within the phaseLatent_hx::Union{Nothing, Float64}
: optional, used for dike heating, J/kgT_liq::Union{Nothing, Float64}
: optional, used for dike heating, liquidus temperature of material, celsiusT_sol::Union{Nothing, Float64}
: optional, used for dike heating, solidus temperature of material, celsiusT_Nu::Union{Nothing, Float64}
: default value for thermal conductivity boundarynu_k::Union{Nothing, Float64}
: optional parameter, Nusselt number for use with conductivityrho_ph::Union{Nothing, String}
: name of the phase diagram you want to use (still needs rho to be defined for the initial guess of pressure)rho_ph_dir::Union{Nothing, String}
: in case the phase diagram has a different path provide the path (without the name of the actual PD) heremfc::Union{Nothing, Float64}
: melt fraction viscosity correction factor (positive scalar)
LaMEM.LaMEM_Model.PhaseTransition
— TypeDefines phase transitions on markers (that change the Phase ID of a marker depending on some conditions)
ID::Int64
: Phase_transition law IDType::String
: [Constant, Clapeyron, Box]: Constant - the phase transition occurs only at a fixed value of the parameter; Clapeyron - clapeyron slopeName_Clapeyron::Union{Nothing, Int64}
: Type of predefined Clapeyron slope, such as MantleTransition660kmPTBox_Bounds::Union{Nothing, Vector{Float64}}
: box bound coordinates: [left, right, front, back, bottom, top]BoxVicinity::Union{Nothing, Int64}
: 1: only check particles in the vicinity of the box boundaries (2: in all directions)Parameter_transition::String
: [T = Temperature, P = Pressure, Depth = z-coord, X=x-coord, Y=y-coord, APS = accumulated plastic strain, MeltFraction, t = time] parameter that triggers the phase transitionConstantValue::Union{Nothing, Float64}
: Value of the parameter [unit of T,P,z, APS]number_phases::Union{Nothing, Int64}
: The number of involved phases [default=1]PhaseAbove::Union{Nothing, Vector{Int64}}
: Above the chosen value the phase is 1, below it, the value is PhaseBelowPhaseBelow::Union{Nothing, Vector{Int64}}
PhaseInside::Union{Nothing, Vector{Int64}}
: Phase within the box [use -1 if you don't want to change the phase inside the box]PhaseOutside::Union{Nothing, Vector{Int64}}
: Phase outside the box [use -1 if you don't want to change the phase outside the box. If combined with OutsideToInside, all phases that come in are set to PhaseInside]PhaseDirection::String
: [BothWays=default; BelowToAbove; AboveToBelow] Direction in which transition worksResetParam::String
: [APS] Parameter to reset on particles below PT or within boxPTBox_TempType::String
: # Temperature condition witin the box [none, constant, linear, halfspace]PTBox_topTemp::Union{Nothing, Float64}
: Temp @ top of box [for linear & halfspace]PTBox_botTemp::Union{Nothing, Float64}
: Temp @ bottom of box [for linear & halfspace]PTBox_thermalAge::Union{Nothing, Float64}
: Thermal age, usually in geo-units [Myrs] [only in case of halfspace]PTBox_cstTemp::Union{Nothing, Float64}
: Temp within box [only for constant T]v_box::Union{Nothing, Float64}
: [optional] only for NotInAirBox, velocity with which box moves in cm/yrt0_box::Union{Nothing, Float64}
: [optional] beginning time of movemen in Myrt1_box::Union{Nothing, Float64}
: [optional] end time of movement in Myr
LaMEM.LaMEM_Model.Scaling
— TypeScaling{T} is a structure that contains the scaling info, employed in the current simulation
Scaling::Any
: Scaling object (as in GeoParams), which can beGEO_units()
,NO_units()
, orSI_units()
LaMEM.LaMEM_Model.Softening
— TypeDefines strain softening parameters
ID::Int64
: softening law IDAPS1::Float64
: Begin of softening, in units of accumulated plastic strain (APS)APS2::Float64
: End of softening, in units of accumulated plastic strain (APS)A::Float64
: Reduction ratioLm::Float64
: Material length scale (in selected units, e.g. km in geo)APSheal2::Union{Nothing, Float64}
: APS when healTau2 activateshealTau::Union{Nothing, Float64}
: healing timescale parameter [Myr]healTau2::Union{Nothing, Float64}
: healing timescale parameter [Myr] starting at APS=APSheal2
LaMEM.LaMEM_Model.SolutionParams
— TypeStructure that contains the LaMEM global solution parameters.
gravity::Vector{Float64}
: gravitational acceleration vectorFSSA::Float64
: free surface stabilization parameter [0 - 1]; The value has to be between 0 and 1
shear_heat_eff::Float64
: shear heating efficiency parameter [0 - 1]Adiabatic_Heat::Float64
: Adiabatic Heating activation flag and efficiency. 0.0 - 1.0act_temp_diff::Int64
: temperature diffusion activation flagact_therm_exp::Int64
: thermal expansion activation flagact_steady_temp::Int64
: steady-state temperature initial guess activation flagsteady_temp_t::Float64
: time for (quasi-)steady-state temperature initial guessnstep_steady::Int64
: number of steps for (quasi-)steady-state temperature initial guess (default = 1)act_heat_rech::Int64
: recharge heat in anomalous bodies after (quasi-)steady-state temperature initial guess (=2: recharge after every diffusion step of initial guess)init_lith_pres::Int64
: sets initial pressure to be the lithostatic pressure (stabilizes compressible setups in the first steps)init_guess::Int64
: create an initial guess step (using constant viscosityeta_ref
before starting the simulationp_litho_visc::Int64
: use lithostatic instead of dynamic pressure for creep lawsp_litho_plast::Int64
: use lithostatic pressure for plasticityp_lim_plast::Int64
: limit pressure at first iteration for plasticityp_shift::Int64
: add a constant value [MPa] to the total pressure field, before evaluating plasticity (e.g., when the domain is located @ some depth within the crust)act_p_shift::Int64
: pressure shift activation flag (enforce zero pressure on average in the top cell layer); note: this overwrites p_shift above!eta_min::Float64
: viscosity lower bound [Pas]eta_max::Float64
: viscosity upper limit [Pas]eta_ref::Float64
: Reference viscosity (used for the initial guess) [Pas]T_ref::Float64
: Reference temperature [C]RUGC::Float64
: universal gas constant (you need to change this only for non-dimensional setups)min_cohes::Float64
: cohesion lower bound [Pa]min_fric::Float64
: friction lower bound [degree]tau_ult::Float64
: ultimate yield stress [Pa]rho_fluid::Float64
: fluid density for depth-dependent density modelgw_level_type::String
: ground water level type for pore pressure computation (see below)gw_level::Float64
: ground water level at the free surface (if defined)biot::Float64
: Biot pressure parameterget_permea::Float64
: effective permeability computation activation flagrescal::Float64
: stencil rescaling flag (for internal constraints, for example while computing permeability)mfmax::Float64
: maximum melt fraction affecting viscosity reductionlmaxit::Int64
: maximum number of local rheology iterationslrtol::Float64
: local rheology iterations relative toleranceact_dike::Int64
: dike activation flag (additonal term in divergence)useTk::Int64
: switch to use T-dependent conductivity, 0: not activedikeHeat::Int64
: switch to use Behn & Ito heat source in the dikeCompute_velocity_gradient::Int64
: compute the velocity gradient tensor 1: active, 0: not active. If active, it automatically activates the output in the .pvd file
LaMEM.LaMEM_Model.Solver
— TypeStructure that contains the LaMEM solver options
SolverType::String
: solver employed ["direct"
or"multigrid"
]DirectSolver::String
: mumps/superlu_dist/pastix/umfpack (requires these external PETSc packages to be installed!)DirectPenalty::Float64
: penalty parameter [employed if we use a direct solver]MGLevels::Int64
: number of MG levels [default=3]MGSweeps::Int64
: number of MG smoothening steps per level [default=10]MGSmoother::String
: type of smoothener used [chebyshev or jacobi]MGJacobiDamp::Float64
: Dampening parameter [only employed for Jacobi smoothener; default=0.6]MGCoarseSolver::String
: coarse grid solver if using multigrid ["direct"
/"mumps"
/"superlu_dist"
or"redundant"
- more options specifiable through the command-line options-crs_ksp_type
&-crs_pc_type
]MGRedundantNum::Int64
: How many times do we copy the coarse grid? [only employed for redundant solver; default is 4]MGRedundantSolver::String
: The coarse grid solver for each of the redundant solves [only employed for redundant; options are"mumps"
/"superlu_dist"
with default"superlu_dist"
]PETSc_options::Vector{String}
: List with (optional) PETSc options
LaMEM.LaMEM_Model.Time
— TypeStructure that contains the LaMEM timestepping information. An explanation of the paramneters is given in the struct `Time_info`
time_end::Float64
: simulation end timedt::Float64
: initial time stepdt_min::Float64
: minimum time step (declare divergence if lower value is attempted)dt_max::Float64
: maximum time stepdt_out::Float64
: output step (output at least at fixed time intervals)inc_dt::Float64
: time step increment per time step (fraction of unit)CFL::Float64
: CFL (Courant-Friedrichs-Lewy) criterionCFLMAX::Float64
: CFL criterion for elasticitynstep_max::Int64
: maximum allowed number of steps (lower bound: timeend/dtmax)nstep_out::Int64
: save output every n steps; Set this to -1 to deactivate saving outputnstep_rdb::Int64
: save restart database every n stepsnum_dt_periods::Int64
: number of time stepping periodstime_dt_periods::Vector{Int64}
: timestamps where timestep should be fixed (first entry has to 0)step_dt_periods::Vector{Float64}
: target timesteps ar timestamps abovenstep_ini::Int64
: save output for n initial stepstime_tol::Float64
: relative tolerance for time comparisons
LaMEM.LaMEM_Model.geom_Sphere
— TypeLaMEM geometric primitive `sphere` object
phase::Int64
: phaseradius::Float64
: radius of spherecenter::Vector{Float64}
: center of sphereTemperature::Union{Nothing, String}
: optional: Temperature of the sphere. possibilities: [constant, or nothing]cstTemp::Union{Nothing, Float64}
: required in case of [constant]: temperature value [in Celcius in case of GEO units]
GeophysicalModelGenerator.AddBox!
— MethodAddBox!(model::Model; xlim=Tuple{2}, [ylim=Tuple{2}], zlim=Tuple{2},
Origin=nothing, StrikeAngle=0, DipAngle=0,
phase = ConstantPhase(1),
T=nothing )
Adds a box with phase & temperature structure to a 3D model setup. This simplifies creating model geometries in geodynamic models See the documentation of the GMG routine
GeophysicalModelGenerator.AddCylinder!
— MethodAddCylinder!(model::Model; # required input
base=Tuple{3}, cap=Tuple{3}, radius=Tuple{1}, # center and radius of the sphere
phase = ConstantPhase(1), # Sets the phase number(s) in the sphere
T=nothing ) # Sets the thermal structure (various fucntions are available)
See the documentation of the GMG routine
GeophysicalModelGenerator.AddEllipsoid!
— MethodAddEllipsoid!(model::Model; # required input
cen=Tuple{3}, axes=Tuple{3}, # center and semi-axes of the ellpsoid
Origin=nothing, StrikeAngle=0, DipAngle=0, # origin & dip/strike
phase = ConstantPhase(1), # Sets the phase number(s) in the box
T=nothing )
See the documentation of the GMG routine
GeophysicalModelGenerator.AddSphere!
— MethodAddSphere!(model::Model; cen=Tuple{3}, radius=Tuple{1}, phase = ConstantPhase(1), T=nothing)
See the documentation of the GMG routine
LaMEM.IO_functions.Read_LaMEM_simulation
— MethodTimestep, FileNames, Time = Read_LaMEM_simulation(model::Model; phase=false, surf=false, passive_tracers=false)
Reads a LaMEM simulation as specified in model
and returns the timesteps, times and filenames of that simulation once it is finished.
LaMEM.IO_functions.Read_LaMEM_timestep
— Functiondata, time = Read_LaMEM_timestep(model::Model, TimeStep::Int64=0; fields=nothing, phase=false, surf=false, last=true)
Reads a specific Timestep
from a simulation specified in model
LaMEM.LaMEM_Model.Check_LaMEM_Model
— MethodCheck_LaMEM_Model(m::Model)
Checks the LaMEM Setup Model m
for errors
LaMEM.LaMEM_Model.Create_Grid
— MethodThis creates a LaMEM grid
LaMEM.LaMEM_Model.Write_LaMEM_InputFile
— FunctionWrite_LaMEM_InputFile(d::Model,fname::String; dir=pwd())
Writes a LaMEM input file based on the data stored in Model
LaMEM.LaMEM_Model.Write_LaMEM_InputFile
— MethodWrite_LaMEM_InputFile(io, d::BoundaryConditions)
Writes the boundary conditions related parameters to file
LaMEM.LaMEM_Model.Write_LaMEM_InputFile
— MethodWrite_LaMEM_InputFile(io, d::FreeSurface)
Writes the free surface related parameters to file
LaMEM.LaMEM_Model.Write_LaMEM_InputFile
— MethodWrite_LaMEM_InputFile(io, d::Grid)
This writes grid info to a LaMEM input file
Example
julia> d=LaMEM.Grid(coord_x=[0.0, 0.7, 0.8, 1.0], bias_x=[0.3,1.0,3.0], nel_x=[10,4,2])
julia> io = open("test.dat","w")
julia> LaMEM.Write_LaMEM_InputFile(io, d)
julia> close(io)
LaMEM.LaMEM_Model.Write_LaMEM_InputFile
— MethodWrite_LaMEM_InputFile(io, d::Output)
Writes the free surface related parameters to file
LaMEM.LaMEM_Model.Write_LaMEM_InputFile
— MethodWrite_LaMEM_InputFile(io, d::ModelSetup)
Writes options related to the Model Setup to disk
LaMEM.LaMEM_Model.Write_LaMEM_InputFile
— MethodWrite_LaMEM_InputFile(io, d::Output)
Writes the free surface related parameters to file
LaMEM.LaMEM_Model.Write_LaMEM_InputFile
— MethodWrite_LaMEM_InputFile(io, d::SolutionParams)
Writes the boundary conditions related parameters to file
LaMEM.LaMEM_Model.Write_LaMEM_InputFile
— MethodWrite_LaMEM_InputFile(io, d::Solver)
Writes the free surface related parameters to file
LaMEM.LaMEM_Model.Write_LaMEM_InputFile
— MethodWrites the Time related parameters to file
LaMEM.LaMEM_Model.Write_LaMEM_InputFile
— MethodWrite_LaMEM_InputFile(io, d::geom_Sphere)
LaMEM.LaMEM_Model.Write_LaMEM_InputFile_PETSc
— MethodWrite_LaMEM_InputFile_PETSc(io, d::Solver)
Writes the (optional) PETSc options to file
LaMEM.LaMEM_Model.add_dike!
— Methodadd_dike!(model::Model, dike::Dike)
This adds a phase transition phase_trans
to model
LaMEM.LaMEM_Model.add_geom!
— Methodadd_geom!(model::Model, geom_object)
This adds an internal geometric primitive object geom_object
to the LaMEM Model Setup model
.
Currently available primitive geom objects are:
geom_Sphere
LaMEM.LaMEM_Model.add_petsc!
— Methodadd_petsc!(model::Model, option::String)
Adds one or more PETSc options to the model
Example
julia> d = Model()
julia> add_petsc!(d,"-snes_npicard 3")
LaMEM.LaMEM_Model.add_phase!
— Methodadd_phase!(model::Model, phase::Phase)
This adds a phase
(with material properties) to model
LaMEM.LaMEM_Model.add_phase!
— Methodadd_phase!(model::Model, phases...)
Add several phases @ once.
LaMEM.LaMEM_Model.add_phasetransition!
— Methodadd_phasetransition!(model::Model, phase_trans::PhaseTransition)
This adds a phase transition phase_trans
to model
LaMEM.LaMEM_Model.add_softening!
— Methodadd_softening!(model::Model, soft::Softening)
This adds a plastic softening law soft
to model
LaMEM.LaMEM_Model.create_initialsetup
— Functioncreate_initialsetup(model::Model, cores::Int64=1, args::String="")
Creates the initial model setup of LaMEM from model
, which includes:
- Writing the LaMEM (*.dat) input file
- Write the VTK file (if requested when
model.Output.write_VTK_setup=true
) - Write the marker files to disk (if
model.ModelSetup.msetup="files"
)
LaMEM.LaMEM_Model.get_doc
— Methodhelp_info::String = get_doc(structure, field::Symbol)
This returns a string with the documentation for a parameter field
that is within the structure
. Note that this structure must be a help structure of the current one.
LaMEM.LaMEM_Model.replace_phase!
— Methodreplace_phase!(model::Model, phase_new::Phase; ID::Int64=nothing, Name::String=nothing)
This replaces a phase
within a LaMEM Model Setup model
with phase_new
either based on its Name
or ID
. Note that it is expected that only one such phase is present in the current setup.
LaMEM.LaMEM_Model.rm_last_phase!
— Methodrm_last_phase!(model::Model, phase::Phase)
This removes the last added phase
from model
LaMEM.LaMEM_Model.rm_phase!
— Methodrm_phase!(model::Model, ID::Int64)
This removes a phase with ID
from model
LaMEM.LaMEM_Model.set_geom!
— MethodThis sets the geometry
LaMEM.Run.run_lamem
— Functionrun_lamem(model::Model, cores::Int64=1, args:String=""; wait=true)
Performs a LaMEM run for the parameters that are specified in model