diff --git a/src/LaMEM_ModelGeneration/ModelSetup.jl b/src/LaMEM_ModelGeneration/ModelSetup.jl index 5ba538a9..da631b87 100644 --- a/src/LaMEM_ModelGeneration/ModelSetup.jl +++ b/src/LaMEM_ModelGeneration/ModelSetup.jl @@ -1,6 +1,8 @@ # Model Setup -export ModelSetup, Write_LaMEM_InputFile, geom_Sphere, set_geom! +export ModelSetup, Write_LaMEM_InputFile, + geom_Sphere, geom_Ellipsoid, geom_Box, geom_RidgeSeg, geom_Hex, geom_Layer, geom_Cylinder, + set_geom! """ Structure that contains the LaMEM Model Setup and Advection options @@ -98,6 +100,8 @@ end # Geometric primitives====== + +# ------- """ LaMEM geometric primitive `sphere` object @@ -144,6 +148,328 @@ function Write_LaMEM_InputFile(io, d::geom_Sphere) println(io, " ") return nothing end +# ------- + +# ------- +""" + LaMEM geometric primitive `Ellipsoid` object + + $(TYPEDFIELDS) +""" +Base.@kwdef struct geom_Ellipsoid + "phase" + phase::Int64 = 1 + + "semi-axes of ellipsoid in `x`, `y` and `z` " + axes::Vector{Float64} = [2.0, 1.5, 1.0] + + "center of sphere" + center::Vector{Float64} = [1.0, 2.0, 3.0] + + "optional: Temperature of the sphere. possibilities: [constant, or nothing]" + Temperature::Union{String,Nothing} = nothing + + "required in case of [constant]: temperature value [in Celcius in case of GEO units]" + cstTemp::Union{Float64,Nothing} = nothing +end + + +function show(io::IO, d::geom_Ellipsoid) + println(io, "Ellipsoid(ph=$(d.phase), axes=$(d.axes), cen=$(d.center), T=$(d.Temperature)=$(d.cstTemp))") + return nothing +end + +""" + Write_LaMEM_InputFile(io, d::geom_Ellipsoid) +""" +function Write_LaMEM_InputFile(io, d::geom_Ellipsoid) + fields = fieldnames(typeof(d)) + println(io, " ") + for f in fields + if !isnothing(getfield(d,f)) + name = rpad(String(f),15) + comment = get_doc(geom_Ellipsoid, f) + data = getfield(d,f) + println(io," $name = $(write_vec(data)) # $(comment)") + end + end + println(io, " ") + return nothing +end +# ------- + +# ------- +""" + LaMEM geometric primitive `Box` object + + $(TYPEDFIELDS) +""" +Base.@kwdef struct geom_Box + "phase" + phase::Int64 = 1 + + "box bound coordinates: `left`, `right`, `front`, `back`, `bottom`, `top` " + bounds::Vector{Float64} = [1.0 2.0 1.0 2.0 1.0 2.0] + + "optional: Temperature structure. possibilities: [constant, linear, halfspace]" + Temperature::Union{String,Nothing} = nothing + + "required in case of [`constant`]: temperature value [in Celcius in case of GEO units]" + cstTemp::Union{Float64,Nothing} = nothing + + "required in case of [`linear,halfspace`]: temperature @ top [in Celcius in case of GEO units]" + topTemp::Union{Float64,Nothing} = nothing + + "required in case of [`linear,halfspace`]: temperature @ top [in Celcius in case of GEO units]" + botTemp::Union{Float64,Nothing} = nothing + + "required in case of [`halfspace`]: thermal age of lithosphere [in Myrs if GEO units are used]" + thermalAge::Union{Float64,Nothing} = nothing +end + +function show(io::IO, d::geom_Box) + println(io, "Box(ph=$(d.phase), bounds=$(d.bounds), T=$(d.Temperature)=$(d.cstTemp), [top,box]=[$(d.topTemp), $(d.botTemp)], thermalAge=$(d.thermalAge))") + return nothing +end + +""" + Write_LaMEM_InputFile(io, d::geom_Box) +""" +function Write_LaMEM_InputFile(io, d::geom_Box) + fields = fieldnames(typeof(d)) + println(io, " ") + for f in fields + if !isnothing(getfield(d,f)) + name = rpad(String(f),15) + comment = get_doc(geom_Box, f) + data = getfield(d,f) + println(io," $name = $(write_vec(data)) # $(comment)") + end + end + println(io, " ") + return nothing +end +# ------- + +# ------- +""" + LaMEM geometric primitive `RidgeSeg` object + + $(TYPEDFIELDS) +""" +Base.@kwdef struct geom_RidgeSeg + "phase" + phase::Int64 = 1 + + "box bound coordinates: `left`, `right`, `front`, `back`, `bottom`, `top` " + bounds::Vector{Float64} = [1.0, 2.0, 1.0, 2.0, 1.0, 2.0] + + "coordinate order: left, right [can be different for oblique ridge]" + ridgeseg_x::Vector{Float64} = [1.5, 1.5] + + "coordinate order: front, back [can be different for oblique ridge]" + ridgeseg_y::Vector{Float64} = [1.0, 2.0] + + "initial temperature structure [ridge must be set to `halfspace_age` --> setTemp=4]" + Temperature::String = "halfspace_age" + + "required in case of [`linear,halfspace`]: temperature @ top [in Celcius in case of GEO units]" + topTemp::Float64 = 0.0 + + "required in case of [`linear,halfspace`]: temperature @ top [in Celcius in case of GEO units]" + botTemp::Float64 = 1300.0 + + "minimum age of seafloor at ridge [in `Myr` in case of GEO units]" + age0::Float64 = 0.001 + + "[optional] parameter that indicates the maximum thermal age of a plate " + maxAge::Union{Float64,Nothing} = nothing + + "[optional] parameter that indicates the spreading velocity of the plate; if not defined it uses bvel_velin specified elsewhere" + v_spread::Union{Float64,Nothing} = nothing +end + +function show(io::IO, d::geom_RidgeSeg) + println(io, "RidgeSeg(ph=$(d.phase), bounds=$(d.bounds), ridgeseg_x=$(d.ridgeseg_x),ridgeseg_y=$(d.ridgeseg_y), T=$(d.Temperature) = [top,box]=[$(d.topTemp), $(d.botTemp)], age0=$(d.age0)), maxAge=$(d.maxAge), v_spread=$(d.v_spread))") + return nothing +end + +""" + Write_LaMEM_InputFile(io, d::geom_RidgeSeg) +""" +function Write_LaMEM_InputFile(io, d::geom_RidgeSeg) + fields = fieldnames(typeof(d)) + println(io, " ") + for f in fields + if !isnothing(getfield(d,f)) + name = rpad(String(f),15) + comment = get_doc(geom_RidgeSeg, f) + data = getfield(d,f) + println(io," $name = $(write_vec(data)) # $(comment)") + end + end + println(io, " ") + return nothing +end +# ------- + +# ------- +""" + LaMEM geometric primitive `Hex` object to define hexahedral elements + + $(TYPEDFIELDS) +""" +Base.@kwdef struct geom_Hex + "phase" + phase::Int64 = 1 + + """ + `x`-`y`-`z` coordinates for each of 8 nodes (24 parameters) + (counter)-clockwise for an arbitrary face, followed by the opposite face + """ + coord::Vector{Float64} = [0.25, 0.25, 0.25, 0.5, 0.2, 0.2, 0.6, 0.7, 0.25, 0.3, 0.5, 0.3, 0.2, 0.3, 0.75, 0.6, 0.15, 0.75, 0.5, 0.6, 0.80, 0.2, 0.4, 0.75] +end + +function show(io::IO, d::geom_Hex) + println(io, "Hex(ph=$(d.phase), coord=$(d.coord))") + return nothing +end + +""" + Write_LaMEM_InputFile(io, d::geom_Hex) +""" +function Write_LaMEM_InputFile(io, d::geom_Hex) + fields = fieldnames(typeof(d)) + println(io, " ") + for f in fields + if !isnothing(getfield(d,f)) + name = rpad(String(f),15) + comment = get_doc(geom_Hex, f) + data = getfield(d,f) + println(io," $name = $(write_vec(data)) # $(comment)") + end + end + println(io, " ") + return nothing +end +# ------- + +# ------- +""" + LaMEM geometric primitive `Layer` object + + $(TYPEDFIELDS) +""" +Base.@kwdef struct geom_Layer + "phase" + phase::Int64 = 1 + + "top of layer" + top::Float64 = 5.0 + + "bottom of layer" + bottom::Float64 = 3.0 + + "optional: add a cosine perturbation on top of the interface (if 1)" + cosine::Union{Int64,Nothing} = nothing + + "required if cosine: wavelength in x-direction" + wavelength::Union{Float64,Nothing} = nothing + + "required if cosine: amplitude of perturbation" + amplitude::Union{Float64,Nothing} = nothing + + "optional: Temperature structure. possibilities: [constant, linear, halfspace]" + Temperature::Union{String,Nothing} = nothing + + "required in case of [`constant`]: temperature value [in Celcius in case of GEO units]" + cstTemp::Union{Float64,Nothing} = nothing + + "required in case of [`linear,halfspace`]: temperature @ top [in Celcius in case of GEO units]" + topTemp::Union{Float64,Nothing} = nothing + + "required in case of [`linear,halfspace`]: temperature @ top [in Celcius in case of GEO units]" + botTemp::Union{Float64,Nothing} = nothing + + "required in case of [`halfspace`]: thermal age of lithosphere [in Myrs if GEO units are used]" + thermalAge::Union{Float64,Nothing} = nothing +end + +function show(io::IO, d::geom_Layer) + println(io, "Layer(ph=$(d.phase), bot/top=[$(d.bot),$(d.top)], T=$(d.Temperature)=$(d.cstTemp), [Ttop,Tbot]=[$(d.topTemp), $(d.botTemp)], thermalAge=$(d.thermalAge), cosine perturbation=$(d.cosine), wavelength=$(d.wavelength), amplitude=$(d.amplitude))") + return nothing +end + +""" + Write_LaMEM_InputFile(io, d::geom_Layer) +""" +function Write_LaMEM_InputFile(io, d::geom_Layer) + fields = fieldnames(typeof(d)) + println(io, " ") + for f in fields + if !isnothing(getfield(d,f)) + name = rpad(String(f),15) + comment = get_doc(geom_Layer, f) + data = getfield(d,f) + println(io," $name = $(write_vec(data)) # $(comment)") + end + end + println(io, " ") + return nothing +end +# ------- + +# ------- +""" + LaMEM geometric primitive `Cylinder` object + + $(TYPEDFIELDS) +""" +Base.@kwdef struct geom_Cylinder + "phase" + phase::Int64 = 1 + + "radius of cylinder " + radius::Float64 = 1.5 + + "center of base of cylinder" + base::Vector{Float64} = [1.0, 2.0, 3.0] + + "center of cap of cylinder" + cap::Vector{Float64} = [3.0, 5.0, 7.0] + + "optional: Temperature structure. possibilities: [constant]" + Temperature::Union{String,Nothing} = nothing + + "required in case of [`constant`]: temperature value [in Celcius in case of GEO units]" + cstTemp::Union{Float64,Nothing} = nothing + +end + +function show(io::IO, d::geom_Cylinder) + println(io, "Cylinder(ph=$(d.phase), radius=$(d.radius) base=$(d.base), cap=$(d.cap), T=$(d.Temperature)=$(d.cstTemp) )") + return nothing +end + +""" + Write_LaMEM_InputFile(io, d::geom_Cylinder) +""" +function Write_LaMEM_InputFile(io, d::geom_Cylinder) + fields = fieldnames(typeof(d)) + println(io, " ") + for f in fields + if !isnothing(getfield(d,f)) + name = rpad(String(f),15) + comment = get_doc(geom_Cylinder, f) + data = getfield(d,f) + println(io," $name = $(write_vec(data)) # $(comment)") + end + end + println(io, " ") + return nothing +end +# ------- + # ========================== @@ -193,7 +519,6 @@ function Write_LaMEM_InputFile(io, d::ModelSetup) Write_LaMEM_InputFile(io, object) end - println("WARNING! Still need to implement geometric primitives!!") end