JuliaGeodynamics / LaMEM.jl

Julia interface to LaMEM (Lithosphere and Mantle Evolution Model)
GNU General Public License v3.0
24 stars 12 forks source link

Passive Tracers support #31

Closed boriskaus closed 5 months ago

boriskaus commented 5 months ago

On popular request, this adds support for PassiveTracers:

julia> model  = Model(Grid(nel=(16,16,16), x=[-1,1], y=[-1,1], z=[-1,1]), PassiveTracers(Passive_Tracer=1, PassiveTracer_Box=[-1,1,-1,1,-1,1]) )
LaMEM Model setup
|
|-- Scaling             :  GeoParams.Units.GeoUnits{GEO}
|-- Grid                :  nel=(16, 16, 16); xϵ(-1.0, 1.0), yϵ(-1.0, 1.0), zϵ(-1.0, 1.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 :  eta_min=1.0e18; eta_max=1.0e25; eta_ref=1.0e20; act_temp_diff=0
|-- Solver options      :  direct solver; superlu_dist; penalty term=10000.0
|-- Model setup options :  Type=files; 
|-- Output options      :  filename=output; pvd=1; avd=0; surf=0
|-- Passive Tracers     :  Type=Always; Res=[100,1,100], Box=[-1.0:1.0,-1.0:1.0,-1.0:1.0]
|-- Materials           :  0 phases; 
julia> matrix = Phase(ID=0,Name="matrix",eta=1e20,rho=3000);
julia> sphere = Phase(ID=1,Name="sphere",eta=1e23,rho=3200)
Phase 1 (sphere): 
  rho    = 3200.0 
  eta    = 1.0e23 
julia> add_phase!(model, sphere, matrix)
julia> AddSphere!(model,cen=(0.0,0.0,0.0), radius=(0.5, ))
julia> run_lamem(model,1)

This will activate paraview output as well if passive tracers are activated (view them in Paraview).

Addresses issue #27

cc @NicolasRiel, @wenrongcao