Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Changed input style for add_polygon! #119

Merged
merged 8 commits into from
Oct 9, 2024
Merged

Changed input style for add_polygon! #119

Show file tree
Hide file tree
Changes from 2 commits
Commits
Show all changes
8 commits
Select commit Hold shift + click to select a range
f933def
Easier handling
TatjanaWeiler May 28, 2024
41b06b6
Update Setup_geometry.jl
TatjanaWeiler May 28, 2024
fd39d65
Making it time stable
TatjanaWeiler Jun 3, 2024
3c5027d
Update Setup_geometry.jl
TatjanaWeiler Jun 18, 2024
e0a2be1
change nicer handling
TatjanaWeiler Jun 19, 2024
db96858
Update Tutorial_polygon_geometry.jl
TatjanaWeiler Oct 9, 2024
c82028d
Update .typos.toml
TatjanaWeiler Oct 9, 2024
6c597cf
Merge branch 'main' into main
boriskaus Oct 9, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 3 additions & 3 deletions docs/src/man/tutorial_Polygon_structures.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -47,13 +47,13 @@ To include the geological structures of a passive margin into the model, we use
# unlimited number of points possible to create the polygon

# add sediment basin
add_polygon!(Phase, Temp, Cart; xlim=[0.0,0.0, 160.0, 200.0],ylim=[100.0,300.0], zlim=[0.0,-10.0,-20.0,0.0], phase = ConstantPhase(8), T=LinearTemp(Ttop=20, Tbot=30));
add_polygon!(Phase, Temp, Cart; xlim=(0.0,0.0, 160.0, 200.0),ylim=(100.0,300.0), zlim=(0.0,-10.0,-20.0,0.0), phase = ConstantPhase(8), T=LinearTemp(Ttop=20, Tbot=30));

# add thinning of the continental crust attached to the slab and its thickness
add_polygon!(Phase, Temp, Cart; xlim=[0.0, 200.0, 0.0],ylim=[500.0,800.0], zlim=[-100.0,-150.0,-150.0], phase = ConstantPhase(5), T=LinearTemp(Ttop=1000, Tbot=1100));
add_polygon!(Phase, Temp, Cart; xlim=(0.0, 200.0, 0.0),ylim=(500.0,800.0), zlim=(-100.0,-150.0,-150.0), phase = ConstantPhase(5), T=LinearTemp(Ttop=1000, Tbot=1100));

# add accretionary prism
add_polygon!(Phase, Temp, Cart; xlim=[800.0, 600.0, 800.0],ylim=[100.0,800.0], zlim=[0.0,0.0,-60.0], phase = ConstantPhase(8), T=LinearTemp(Ttop=20, Tbot=30));
add_polygon!(Phase, Temp, Cart; xlim=(800.0, 600.0, 800.0),ylim=(100.0,800.0), zlim=(0.0,0.0,-60.0), phase = ConstantPhase(8), T=LinearTemp(Ttop=20, Tbot=30));
```

#### 3. Export final model setup to a paraview file
Expand Down
16 changes: 13 additions & 3 deletions src/Setup_geometry.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -593,7 +593,7 @@ end


"""
add_polygon!(Phase, Temp, Grid::AbstractGeneralGrid; xlim::Vector(), ylim=Vector(2), zlim=Vector(), phase = ConstantPhase(1), T=nothing, cell=false )
add_polygon!(Phase, Temp, Grid::AbstractGeneralGrid; xlim=Tuple{}, ylim=Tuple{2}, zlim=Tuple{}, phase = ConstantPhase(1), T=nothing, cell=false )
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What is the point of Tuple{} and Tuple{2}? If you want to make empty tuples you can to just xlim=()

julia> Tuple{} == ()
false


Adds a polygon with phase & temperature structure to a 3D model setup. This simplifies creating model geometries in geodynamic models

Expand Down Expand Up @@ -625,7 +625,7 @@ LaMEM Grid:
z ϵ [-2.0 : 0.0]
julia> Phases = zeros(Int32, size(Grid.X));
julia> Temp = zeros(Float64, size(Grid.X));
julia> add_polygon!(Phase, Temp, Cart; xlim=[0.0,0.0, 1.6, 2.0],ylim=[0.0,0.8], zlim=[0.0,-1.0,-2.0,0.0], phase = ConstantPhase(8), T=ConstantTemp(30))
julia> add_polygon!(Phase, Temp, Cart; xlim=(0,0, 1.6, 2.0),ylim=(0,0.8), zlim=(0,-1,-2,0), phase = ConstantPhase(8), T=ConstantTemp(30))
julia> Model3D = ParaviewData(Grid, (Phases=Phases,Temp=Temp)); # Create Cartesian model
julia> write_paraview(Model3D,"LaMEM_ModelSetup") # Save model to paraview
1-element Vector{String}:
Expand All @@ -634,10 +634,20 @@ julia> write_paraview(Model3D,"LaMEM_ModelSetup") # Save model to para

"""
function add_polygon!(Phase, Temp, Grid::AbstractGeneralGrid; # required input
xlim::Vector=[], ylim::Vector=[], zlim::Vector=[], # limits of the box
xlim=Tuple{}, ylim=Tuple{2}, zlim=Tuple{}, # limits of the box
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Same as above, these are not empty tuples.

phase = ConstantPhase(1), # Sets the phase number(s) in the box
T=nothing, cell=false ) # Sets the thermal structure (various functions are available)


xlim = collect(xlim)
Copy link
Member

@albert-de-montserrat albert-de-montserrat Jun 1, 2024 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Here you are changing the type of xlim from whatever it is originally to a Array. This yields a type instability, you could fix it just by changing the name of the xlim either in the left hand side or the right hind side of the expression. Same for ylim and zlim

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for your comments. I changed the names of xlim, ylim and zlim.
I have not changed the input style because I used the same as for add_box and all the other add function. Should I change that everywhere? And should I also change that in the help functions?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I have not changed the input style because I used the same as for add_box and all the other add function. Should I change that everywhere? And should I also change that in the help functions?

It would be a good idea to change everything that looks like that indeed. I have the suspicion that all/most those functions will crash if the keyword arguments are not provided

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Many functions need a Tuple with a specific length and not an empty one. Would e. g. Tuple{Float64, Float64} be better even though at the moment Int64 and Float64 are possible or is it better to place everywhere empty tuples and define it only in the help function what kind of input is needed?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The problem is the following. If you do, for example,

x=Tuple{2}

you can see that x is not a actual usable tuple. x is actually a type if you define it that way.

julia> typeof(x)
DataType

This means that if the user calls any of these functions without explicitly passing x into the function, the function will crash because x is not an actual tuple. I guess this is what you want to do:

julia> function foo(; x::NTuple{2, T} = (1,1)) where T
           return x
       end
foo (generic function with 1 method)

Now you can see how this works depending on what I pass as x

julia> foo(; x=(1,1))
(1, 1)

julia> foo()
(1, 1)

julia> foo(; x=(1.0, 1.0))
(1.0, 1.0)

julia> foo(; x=(1.0, 1.0, 1.0))
ERROR: MethodError: no method matching var"#foo#5"(::Tuple{Float64, Float64, Float64}, ::typeof(foo))

Closest candidates are:
  var"#foo#5"(::Tuple{T, T}, ::typeof(foo)) where T
   @ Main ~/Desktop/DevPkgs/GeophysicalModelGenerator.jl/volcano.jl:70

Stacktrace:
 [1] top-level scope
   @ ~/Desktop/DevPkgs/GeophysicalModelGenerator.jl/volcano.jl:77

if you define x like that in the function call, this means that it can be a Tuple{T, T} where T Any, but in the right hand side you need to provide a default value for x; in this case x=(1,1), but could be whatever other tuple.

Hope this makes any sense!

ylim = collect(ylim)
zlim = collect(zlim)

xlim = Float64.(xlim)
ylim = Float64.(ylim)
zlim = Float64.(zlim)


# Retrieve 3D data arrays for the grid
X,Y,Z = coordinate_grids(Grid, cell=cell)

Expand Down
2 changes: 1 addition & 1 deletion test/test_setup_geometry.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -230,7 +230,7 @@ Temp = ones(Float64,(length(x),length(y),length(z)))*1350;
add_box!(Phase, Temp, Cart; xlim=(0.0,600.0),ylim=(0.0,600.0), zlim=(-80.0, 0.0), phase = ConstantPhase(5), T=T=ConstantTemp(120.0));

# add accretionary prism
add_polygon!(Phase, Temp, Cart; xlim=[500.0, 200.0, 500.0],ylim=[100.0,400.0], zlim=[0.0,0.0,-60.0], phase = ConstantPhase(8), T=LinearTemp(Ttop=20, Tbot=30))
add_polygon!(Phase, Temp, Cart; xlim=(500.0, 200.0, 500.0),ylim=(100.0,400.0), zlim=(0.0,0.0,-60.0), phase = ConstantPhase(8), T=LinearTemp(Ttop=20, Tbot=30))

@test maximum(Phase) == 8
@test minimum(Temp) == 21.40845070422536
Expand Down
Loading