From bf89d251be0abf59b6ebc254b209f67c311b1871 Mon Sep 17 00:00:00 2001 From: Boris Kaus Date: Thu, 8 Feb 2024 15:24:02 +0100 Subject: [PATCH 01/10] add AddLayer and first set of default parameters --- src/LaMEM.jl | 2 +- src/LaMEM_ModelGeneration/DefaultParams.jl | 23 ++++++++++++++++++++++ src/LaMEM_ModelGeneration/GMG_interface.jl | 15 ++++++++++++-- src/LaMEM_ModelGeneration/Utils.jl | 22 ++++++++++++++++++++- 4 files changed, 58 insertions(+), 4 deletions(-) diff --git a/src/LaMEM.jl b/src/LaMEM.jl index b9fc50b..fa48f5c 100644 --- a/src/LaMEM.jl +++ b/src/LaMEM.jl @@ -36,7 +36,7 @@ export LaMEM_Model, Model, Write_LaMEM_InputFile, create_initialsetup, add_phase!, rm_phase!, rm_last_phase!, replace_phase!, add_petsc!, add_softening!, add_phaseaggregate!, add_phasetransition!, add_dike!, add_geom!, set_air, copy_phase, add_topography!, AboveSurface!, BelowSurface!, - prepare_lamem + prepare_lamem, isdefault using .Run.LaMEM_jll diff --git a/src/LaMEM_ModelGeneration/DefaultParams.jl b/src/LaMEM_ModelGeneration/DefaultParams.jl index 3fdc213..ba3051d 100644 --- a/src/LaMEM_ModelGeneration/DefaultParams.jl +++ b/src/LaMEM_ModelGeneration/DefaultParams.jl @@ -35,6 +35,29 @@ function UpdateDefaultParameters(model::Model) model.SolutionParams.act_p_shift = 1 end + # output additional fields at all times; stress, strainrate, density, pressure, velocity, temperature + + # is using MG and 2D , set da_refine_y to 1 + + # if we have a free surface, you'll generally want output + if model.FreeSurface.surf_use==1 + model.Output.out_surf=1 + model.Output.out_surf_pvd = 1 + model.Output.out_surf_velocity = 1 + model.Output.out_surf_topography = 1 + model.Output.out_surf_amplitude = 1 + end + + # Scaling: if we use default values, employ smarter default values based on the model setup + if isdefault(model.Scaling,Scaling()) + le = abs(diff(model.Grid.coord_z)[1])*km # length + η = model.SolutionParams.eta_ref*Pas # viscosity + model.Scaling = Scaling(GEO_units(length=le, viscosity = η)) + end + + # exx_strain_rates: no need to specify exx_num_periods + + # if surf_use=1 and surf_level==nothing, set it to zero return model end diff --git a/src/LaMEM_ModelGeneration/GMG_interface.jl b/src/LaMEM_ModelGeneration/GMG_interface.jl index 215a9bb..df2213e 100644 --- a/src/LaMEM_ModelGeneration/GMG_interface.jl +++ b/src/LaMEM_ModelGeneration/GMG_interface.jl @@ -2,7 +2,7 @@ # # Some wrappers around GMG routines -import GeophysicalModelGenerator: AddBox!, AddSphere!, AddEllipsoid!, AddCylinder!, AboveSurface, BelowSurface +import GeophysicalModelGenerator: AddBox!, AddLayer!, AddSphere!, AddEllipsoid!, AddCylinder!, AboveSurface, BelowSurface export AboveSurface!, BelowSurface! """ @@ -12,11 +12,22 @@ export AboveSurface!, BelowSurface! 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 +See the documentation of the GMG routine for the full options. """ AddBox!(model::Model; kwargs...) = AddBox!(model.Grid.Phases, model.Grid.Temp, model.Grid.Grid; kwargs...) +""" + AddLayer!(model::Model; xlim, ylim, zlim=Tuple{2}, + phase = ConstantPhase(1), + T=nothing ) + +Adds a layer with phase & temperature structure to a 3D model setup. This simplifies creating model geometries in geodynamic models +See the documentation of the GMG routine for the full options. + +""" +AddLayer!(model::Model; kwargs...) = AddLayer!(model.Grid.Phases, model.Grid.Temp, model.Grid.Grid; kwargs...) + """ AddSphere!(model::Model; cen=Tuple{3}, radius=Tuple{1}, phase = ConstantPhase(1), T=nothing) diff --git a/src/LaMEM_ModelGeneration/Utils.jl b/src/LaMEM_ModelGeneration/Utils.jl index ef68cd6..a5fd693 100644 --- a/src/LaMEM_ModelGeneration/Utils.jl +++ b/src/LaMEM_ModelGeneration/Utils.jl @@ -5,7 +5,9 @@ export add_phase!, rm_phase!, rm_last_phase!, replace_phase!, add_vbox!, rm_last_vbox!, rm_vbox!, add_softening!, add_phasetransition!, add_phaseaggregate!, add_dike!, add_geom! , cross_section, - set_air, copy_phase + set_air, copy_phase, + isdefault + """ @@ -370,4 +372,22 @@ function add_topography!(model::Model, topography::CartData; surf_air_phase=0, s model.BoundaryConditions.open_top_bound=open_top_bound return nothing +end + + +""" + isdefault(s1::S, s_default::S) + +Checks whether a struct `s1` has default parameters `s_default` +""" +function isdefault(s1, s_default) + + default = true + for f in fieldnames(typeof(s1)) + if getfield(s1,f) != getfield(s_default,f) + default = false + end + end + + return default end \ No newline at end of file From b943736a22541878122fd1cf20168e40cbefc791 Mon Sep 17 00:00:00 2001 From: Boris Kaus Date: Thu, 8 Feb 2024 16:04:54 +0100 Subject: [PATCH 02/10] more defaults & dir issues --- src/LaMEM_ModelGeneration/DefaultParams.jl | 29 +++++++++++++++++++--- src/LaMEM_ModelGeneration/Model.jl | 10 ++++---- src/LaMEM_ModelGeneration/ModelSetup.jl | 5 ++-- 3 files changed, 33 insertions(+), 11 deletions(-) diff --git a/src/LaMEM_ModelGeneration/DefaultParams.jl b/src/LaMEM_ModelGeneration/DefaultParams.jl index ba3051d..8e6636d 100644 --- a/src/LaMEM_ModelGeneration/DefaultParams.jl +++ b/src/LaMEM_ModelGeneration/DefaultParams.jl @@ -35,9 +35,17 @@ function UpdateDefaultParameters(model::Model) model.SolutionParams.act_p_shift = 1 end - # output additional fields at all times; stress, strainrate, density, pressure, velocity, temperature + # Output some fields @ all times + model.Output.out_density = 1 + model.Output.out_j2_dev_stress = 1 + model.Output.out_j2_strain_rate = 1 + model.Output.out_pressure = 1 + model.Output.out_temperature = 1 - # is using MG and 2D , set da_refine_y to 1 + # if using multigrid in a 2D setup + if model.Solver.SolverType=="multigrid" && model.Grid.nel_y[1]==1 + push!(model.Solver.PETSc_options,"-da_refine_y 1") + end # if we have a free surface, you'll generally want output if model.FreeSurface.surf_use==1 @@ -56,8 +64,23 @@ function UpdateDefaultParameters(model::Model) end # exx_strain_rates: no need to specify exx_num_periods + if length(model.BoundaryConditions.exx_strain_rates)==1 & model.BoundaryConditions.exx_num_periods>1 + # we only specified a single strainrate so it'll be constant with time + model.BoundaryConditions.exx_num_periods = 1 + model.BoundaryConditions.exx_time_delims = 1e20; + end + + if model.FreeSurface.surf_use==1 + # if we use a free surface & surf_level==nothing, set it to zero + if isnothing(model.FreeSurface.surf_level) + model.FreeSurface.surf_level = 0.0; + end + if isnothing(model.FreeSurface.surf_air_phase) + model.FreeSurface.surf_air_phase = 0 + end - # if surf_use=1 and surf_level==nothing, set it to zero + + end return model end diff --git a/src/LaMEM_ModelGeneration/Model.jl b/src/LaMEM_ModelGeneration/Model.jl index 6f2824f..1f1077b 100644 --- a/src/LaMEM_ModelGeneration/Model.jl +++ b/src/LaMEM_ModelGeneration/Model.jl @@ -188,13 +188,15 @@ Performs a LaMEM run for the parameters that are specified in `model` """ function run_lamem(model::Model, cores::Int64=1, args::String=""; wait=true) + cur_dir = pwd(); + + #if !isdir(model.Output.out_dir); mkdir(model.Output.out_dir); end # create directory if needed create_initialsetup(model, cores, args); - cur_dir = pwd(); if !isempty(model.Output.out_dir) cd(model.Output.out_dir) end - + run_lamem(model.Output.param_file_name, cores, args; wait=wait) cd(cur_dir) @@ -219,9 +221,7 @@ function prepare_lamem(model::Model, cores::Int64=1, args::String=""; verbose=fa println("Creating LaMEM input files in the directory: $(model.Output.out_dir)") cur_dir = pwd(); - if !isempty(model.Output.out_dir) - cd(model.Output.out_dir) - end + create_initialsetup(model, cores, args, verbose=verbose); cd(cur_dir) diff --git a/src/LaMEM_ModelGeneration/ModelSetup.jl b/src/LaMEM_ModelGeneration/ModelSetup.jl index 902ec5d..5ba538a 100644 --- a/src/LaMEM_ModelGeneration/ModelSetup.jl +++ b/src/LaMEM_ModelGeneration/ModelSetup.jl @@ -54,7 +54,7 @@ Base.@kwdef mutable struct ModelSetup nmark_avd::Vector{Int64} = [3, 3, 3] "max number of same phase markers per subcell (subgrid marker control)" - nmark_sub::Int64 = 1 + nmark_sub::Int64 = 3 "Different geometric primitives that can be selected if we `msetup``=`geom`; see `geom_Sphere`" geom_primitives::Vector = [] @@ -169,12 +169,11 @@ function Write_LaMEM_InputFile(io, d::ModelSetup) (f == :msetup) || (f == :rand_noise) || (f == :nmark_lim) || + (f == :nmark_sub) || (f == :advect) || (f == :interp) || (f == :mark_ctrl) - - if (f != :geom_primitives) # only print if value differs from reference value From c7711c5cfaaaca1edbd4456608d8d2ff89b12120 Mon Sep 17 00:00:00 2001 From: Boris Kaus Date: Thu, 8 Feb 2024 16:15:30 +0100 Subject: [PATCH 03/10] update tests --- test/mesh_refinement_test.jl | 2 +- test/test_julia_setups.jl | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/test/mesh_refinement_test.jl b/test/mesh_refinement_test.jl index c890c6d..aac9512 100644 --- a/test/mesh_refinement_test.jl +++ b/test/mesh_refinement_test.jl @@ -31,7 +31,7 @@ using GeophysicalModelGenerator # read last timestep data,time = Read_LaMEM_timestep(model,last=true); - @test sum(data.fields.velocity[3][:,:,:]) ≈ 0.36892682f0 # check Vz + @test sum(data.fields.velocity[3][:,:,:]) ≈ 0.3680135f0 # check Vz # cleanup the directory rm(model.Output.out_dir, force=true, recursive=true) diff --git a/test/test_julia_setups.jl b/test/test_julia_setups.jl index d0164d9..375c075 100644 --- a/test/test_julia_setups.jl +++ b/test/test_julia_setups.jl @@ -151,7 +151,7 @@ using GeophysicalModelGenerator # read last timestep data,time = Read_LaMEM_timestep(model,last=true); - @test sum(data.fields.phase) ≈ 29047.736f0 + @test sum(data.fields.phase) ≈ 29060.664f0 # cleanup the directory rm(model.Output.out_dir, force=true, recursive=true) From d4d12a14f4364e1fe1a5adaab6bbf245d5f8a673 Mon Sep 17 00:00:00 2001 From: Boris Kaus Date: Thu, 8 Feb 2024 16:15:44 +0100 Subject: [PATCH 04/10] fix strainrate default --- src/LaMEM_ModelGeneration/DefaultParams.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/LaMEM_ModelGeneration/DefaultParams.jl b/src/LaMEM_ModelGeneration/DefaultParams.jl index 8e6636d..1bc4160 100644 --- a/src/LaMEM_ModelGeneration/DefaultParams.jl +++ b/src/LaMEM_ModelGeneration/DefaultParams.jl @@ -64,10 +64,10 @@ function UpdateDefaultParameters(model::Model) end # exx_strain_rates: no need to specify exx_num_periods - if length(model.BoundaryConditions.exx_strain_rates)==1 & model.BoundaryConditions.exx_num_periods>1 + if length(model.BoundaryConditions.exx_strain_rates)==1 && model.BoundaryConditions.exx_num_periods>1 # we only specified a single strainrate so it'll be constant with time model.BoundaryConditions.exx_num_periods = 1 - model.BoundaryConditions.exx_time_delims = 1e20; + model.BoundaryConditions.exx_time_delims = [1e20]; end if model.FreeSurface.surf_use==1 From 24379d01a96ef7f695b41c243928c62e49b56764 Mon Sep 17 00:00:00 2001 From: Boris Kaus Date: Mon, 12 Feb 2024 12:46:30 +0100 Subject: [PATCH 05/10] check for plasticity --- src/LaMEM.jl | 2 +- src/LaMEM_ModelGeneration/FreeSurface.jl | 3 ++- src/LaMEM_ModelGeneration/Materials.jl | 2 +- src/LaMEM_ModelGeneration/Model.jl | 8 +++++++- src/LaMEM_ModelGeneration/SolutionParams.jl | 4 +++- src/LaMEM_ModelGeneration/Utils.jl | 15 ++++++++++++++- 6 files changed, 28 insertions(+), 6 deletions(-) diff --git a/src/LaMEM.jl b/src/LaMEM.jl index df97647..b088c4b 100644 --- a/src/LaMEM.jl +++ b/src/LaMEM.jl @@ -36,7 +36,7 @@ export LaMEM_Model, Model, Write_LaMEM_InputFile, create_initialsetup, add_phase!, rm_phase!, rm_last_phase!, replace_phase!, add_petsc!, add_softening!, add_phaseaggregate!, add_phasetransition!, add_dike!, add_geom!, set_air, copy_phase, add_topography!, AboveSurface!, BelowSurface!, - prepare_lamem, isdefault + prepare_lamem, isdefault, hasplasticity using .Run.LaMEM_jll diff --git a/src/LaMEM_ModelGeneration/FreeSurface.jl b/src/LaMEM_ModelGeneration/FreeSurface.jl index fb2ebbf..4427f4f 100644 --- a/src/LaMEM_ModelGeneration/FreeSurface.jl +++ b/src/LaMEM_ModelGeneration/FreeSurface.jl @@ -129,7 +129,8 @@ function Write_LaMEM_InputFile(io, d::FreeSurface) if surf_use==1 for f in fields - if getfield(d,f) != getfield(Reference,f) && (f != :Topography) + if (getfield(d,f) != getfield(Reference,f) && (f != :Topography)) + # only print if value differs from reference value name = rpad(String(f),15) comment = get_doc(FreeSurface, f) diff --git a/src/LaMEM_ModelGeneration/Materials.jl b/src/LaMEM_ModelGeneration/Materials.jl index 08e75df..85c0d76 100644 --- a/src/LaMEM_ModelGeneration/Materials.jl +++ b/src/LaMEM_ModelGeneration/Materials.jl @@ -221,7 +221,7 @@ Base.@kwdef mutable struct Softening A::Float64 = 0.7 "Material length scale (in selected units, e.g. km in geo)" - Lm::Float64 = 0.2 + Lm::Union{Nothing,Float64} = nothing # healing parameters "APS when healTau2 activates" diff --git a/src/LaMEM_ModelGeneration/Model.jl b/src/LaMEM_ModelGeneration/Model.jl index 1f1077b..f906765 100644 --- a/src/LaMEM_ModelGeneration/Model.jl +++ b/src/LaMEM_ModelGeneration/Model.jl @@ -160,7 +160,13 @@ function Write_LaMEM_InputFile(d::Model, fname::String="input.dat"; dir=pwd()) # If we want to write an input file Write_Paraview(CartData(d.Grid.Grid, (Phases=d.Grid.Phases,Temp=d.Grid.Temp)),"Model3D") end - + + if any(hasplasticity.(d.Materials.Phases)) + # We have plasticity, so we likely want to see that + d.Output.out_plast_strain = 1 # accumulated plastic strain + d.Output.out_plast_dissip = 1 # plastic dissipation + end + io = open(fname,"w") Write_LaMEM_InputFile(io, d.Scaling) diff --git a/src/LaMEM_ModelGeneration/SolutionParams.jl b/src/LaMEM_ModelGeneration/SolutionParams.jl index 0256e2e..9abce65 100644 --- a/src/LaMEM_ModelGeneration/SolutionParams.jl +++ b/src/LaMEM_ModelGeneration/SolutionParams.jl @@ -183,7 +183,9 @@ function Write_LaMEM_InputFile(io, d::SolutionParams) if getfield(d,f) != getfield(Reference,f) || (f == :eta_ref) || (f == :gravity) || - (f == :act_p_shift) + (f == :act_p_shift) || + (f == :init_guess) || + (f == :p_lim_plast) # only print if value differs from reference value name = rpad(String(f),15) diff --git a/src/LaMEM_ModelGeneration/Utils.jl b/src/LaMEM_ModelGeneration/Utils.jl index a5fd693..46a3532 100644 --- a/src/LaMEM_ModelGeneration/Utils.jl +++ b/src/LaMEM_ModelGeneration/Utils.jl @@ -6,7 +6,7 @@ export add_phase!, rm_phase!, rm_last_phase!, replace_phase!, add_softening!, add_phasetransition!, add_phaseaggregate!, add_dike!, add_geom! , cross_section, set_air, copy_phase, - isdefault + isdefault, hasplasticity @@ -390,4 +390,17 @@ function isdefault(s1, s_default) end return default +end + +""" + hasplasticity(p::Phase) + +`true` if `p` contains plastic parameters (cohesion or friction angle) +""" +function hasplasticity(p::Phase) + plastic = false + if !isnothing(p.ch) || !isnothing(p.fr) + plastic = true + end + return plastic end \ No newline at end of file From a5b3ac9e56cf0960ce61170ccc47e1cf76d443a0 Mon Sep 17 00:00:00 2001 From: Boris Kaus Date: Tue, 13 Feb 2024 10:04:11 +0100 Subject: [PATCH 06/10] upcoming addition to LaMEM --- src/LaMEM_ModelGeneration/Materials.jl | 3 +++ 1 file changed, 3 insertions(+) diff --git a/src/LaMEM_ModelGeneration/Materials.jl b/src/LaMEM_ModelGeneration/Materials.jl index 85c0d76..0f87946 100644 --- a/src/LaMEM_ModelGeneration/Materials.jl +++ b/src/LaMEM_ModelGeneration/Materials.jl @@ -120,6 +120,9 @@ Base.@kwdef mutable struct Phase "stabilization viscosity (default is eta_min)" eta_st::Union{Nothing,Float64} = nothing + "viscoplastic plasticity regularisation viscosity" + eta_vp::Union{Nothing,Float64} = nothing + "pore-pressure ratio" rp::Union{Nothing,Float64} = nothing From 33166ae01ea47e2e3a04b8e805f98121498647e8 Mon Sep 17 00:00:00 2001 From: Boris Kaus Date: Wed, 14 Feb 2024 09:19:54 +0100 Subject: [PATCH 07/10] set mumps as default; add :all to Plots --- src/LaMEM_ModelGeneration/DefaultParams.jl | 5 ++++ src/PlotsExt.jl | 35 ++++++++++++++++++++-- 2 files changed, 37 insertions(+), 3 deletions(-) diff --git a/src/LaMEM_ModelGeneration/DefaultParams.jl b/src/LaMEM_ModelGeneration/DefaultParams.jl index 1bc4160..9fc4d11 100644 --- a/src/LaMEM_ModelGeneration/DefaultParams.jl +++ b/src/LaMEM_ModelGeneration/DefaultParams.jl @@ -42,6 +42,11 @@ function UpdateDefaultParameters(model::Model) model.Output.out_pressure = 1 model.Output.out_temperature = 1 + if isdefault(model.Solver,Solver()) + # the LaMEM default (superlu_dist) currently doesn't work with PETSc_jll + model.Solver.DirectSolver = "mumps"; + end + # if using multigrid in a 2D setup if model.Solver.SolverType=="multigrid" && model.Grid.nel_y[1]==1 push!(model.Solver.PETSc_options,"-da_refine_y 1") diff --git a/src/PlotsExt.jl b/src/PlotsExt.jl index ce39d08..a12db4c 100644 --- a/src/PlotsExt.jl +++ b/src/PlotsExt.jl @@ -7,13 +7,18 @@ import .Plots: heatmap export plot_topo, plot_cross_section, plot_phasediagram """ - plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol=:phase, dim=1, x=nothing, y=nothing, z=nothing, aspect_ratio::Union{Real, Symbol}=:equal) + plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol=:phase, + timestep::Union{Nothing, Int64}=nothing, + dim=1, x=nothing, y=nothing, z=nothing, aspect_ratio::Union{Real, Symbol}=:equal) This plots a cross-section through the LaMEM `model`, of the field `field`. If that is a vector or tensor field, specify `dim` to indicate which of the fields you want to see. If the keyword `timestep` is specified, it will load a timestep """ -function plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol=:phase, timestep=nothing, dim=1, x=nothing, y=nothing, z=nothing, aspect_ratio::Union{Real, Symbol}=:equal) - +function plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol=:phase, + timestep::Union{Nothing, Int64}=nothing, dim=1, + x=nothing, y=nothing, z=nothing, + aspect_ratio::Union{Real, Symbol}=:equal) + if !isnothing(timestep) # load a particular timestep data_cart, time = Read_LaMEM_timestep(model,timestep) @@ -50,6 +55,30 @@ function plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol return hm end +""" + plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol=:phase, + timestep::Symbol=:all, + dim=1, x=nothing, y=nothing, z=nothing, aspect_ratio::Union{Real, Symbol}=:equal) + +This allows plotting all timesteps if you specify `timestep=:all` +""" +function plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol=:phase, + timestep::Symbol, dim=1, x=nothing, y=nothing, z=nothing, aspect_ratio::Union{Real, Symbol}=:equal) + + if timestep==:all + Timesteps,_,_ = Read_LaMEM_simulation(model); + for timestep_val in Timesteps + plot_cross_section(model, args, field=field, timestep=timestep_val, dim=dim, x=x, y=y, z=z, aspect_ratio=aspect_ratio) + end + + else + error("unknown timestep: $timestep") + end + + + return nothing +end + """ plot_topo(topo::CartData, args...) From 90d98f15acd7f7c112665eaaa5ef9f032093f20a Mon Sep 17 00:00:00 2001 From: Boris Kaus Date: Wed, 14 Feb 2024 09:23:48 +0100 Subject: [PATCH 08/10] addition to plot --- src/PlotsExt.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/PlotsExt.jl b/src/PlotsExt.jl index a12db4c..5fe5362 100644 --- a/src/PlotsExt.jl +++ b/src/PlotsExt.jl @@ -63,7 +63,7 @@ end This allows plotting all timesteps if you specify `timestep=:all` """ function plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol=:phase, - timestep::Symbol, dim=1, x=nothing, y=nothing, z=nothing, aspect_ratio::Union{Real, Symbol}=:equal) + timestep::Symbol=:all, dim=1, x=nothing, y=nothing, z=nothing, aspect_ratio::Union{Real, Symbol}=:equal) if timestep==:all Timesteps,_,_ = Read_LaMEM_simulation(model); From 1b04a8e22e78f7ecba2df9e6c5509aaf23cc218d Mon Sep 17 00:00:00 2001 From: Boris Kaus Date: Wed, 14 Feb 2024 09:27:16 +0100 Subject: [PATCH 09/10] deactivatinga s not working --- src/PlotsExt.jl | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/PlotsExt.jl b/src/PlotsExt.jl index 5fe5362..0b30db0 100644 --- a/src/PlotsExt.jl +++ b/src/PlotsExt.jl @@ -55,6 +55,7 @@ function plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol return hm end +#= """ plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol=:phase, timestep::Symbol=:all, @@ -78,7 +79,7 @@ function plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol return nothing end - +=# """ plot_topo(topo::CartData, args...) From ad0392f5e915fa204bb047992e400eb2e99d65fc Mon Sep 17 00:00:00 2001 From: Boris Kaus Date: Fri, 16 Feb 2024 15:23:12 +0100 Subject: [PATCH 10/10] add plot_cross_section_simulation to see the full simulations --- src/PlotsExt.jl | 31 ++++++++++++------------------- 1 file changed, 12 insertions(+), 19 deletions(-) diff --git a/src/PlotsExt.jl b/src/PlotsExt.jl index 0b30db0..bfeaf52 100644 --- a/src/PlotsExt.jl +++ b/src/PlotsExt.jl @@ -4,7 +4,7 @@ println("Adding Plots.jl plotting extensions for LaMEM") using LaMEM, GeophysicalModelGenerator, Statistics, DelimitedFiles using .Plots import .Plots: heatmap -export plot_topo, plot_cross_section, plot_phasediagram +export plot_topo, plot_cross_section, plot_phasediagram, plot_cross_section_simulation """ plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol=:phase, @@ -55,31 +55,24 @@ function plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol return hm end -#= + """ - plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol=:phase, - timestep::Symbol=:all, + plot_cross_section_simulation(model::Union{Model,CartData}, args...; field::Symbol=:phase, dim=1, x=nothing, y=nothing, z=nothing, aspect_ratio::Union{Real, Symbol}=:equal) -This allows plotting all timesteps if you specify `timestep=:all` +As `plot_cross_section`, but for the entire simulation instead of a single timestep. """ -function plot_cross_section(model::Union{Model,CartData}, args...; field::Symbol=:phase, - timestep::Symbol=:all, dim=1, x=nothing, y=nothing, z=nothing, aspect_ratio::Union{Real, Symbol}=:equal) - - if timestep==:all - Timesteps,_,_ = Read_LaMEM_simulation(model); - for timestep_val in Timesteps - plot_cross_section(model, args, field=field, timestep=timestep_val, dim=dim, x=x, y=y, z=z, aspect_ratio=aspect_ratio) - end - - else - error("unknown timestep: $timestep") - end +function plot_cross_section_simulation(model::Union{Model,CartData}, args...; field::Symbol=:phase, + dim=1, x=nothing, y=nothing, z=nothing, aspect_ratio::Union{Real, Symbol}=:equal) + Timesteps,_,_ = Read_LaMEM_simulation(model); + for timestep_val in Timesteps + plot_cross_section(model, args, field=field, timestep=timestep_val, dim=dim, x=x, y=y, z=z, aspect_ratio=aspect_ratio) + end - return nothing + return nothing end -=# + """ plot_topo(topo::CartData, args...)