Skip to content

Commit 306c097

Browse files
buchanankerswellbuchanankerswell
committed
Add coordinate system to particle/property/function
1 parent a69b574 commit 306c097

File tree

12 files changed

+438
-7
lines changed

12 files changed

+438
-7
lines changed

‎include/aspect/particle/property/function.h‎

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -88,6 +88,12 @@ namespace aspect
8888
*/
8989
std::unique_ptr<Functions::ParsedFunction<dim>> function;
9090

91+
/**
92+
* The coordinate representation to evaluate the function. Possible
93+
* choices are depth, cartesian and spherical.
94+
*/
95+
Utilities::Coordinates::CoordinateSystem coordinate_system;
96+
9197
/**
9298
* A private variable that stores the number of particle property
9399
* function components.

‎source/particle/property/function.cc‎

Lines changed: 25 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -19,6 +19,7 @@
1919
*/
2020

2121
#include <aspect/particle/property/function.h>
22+
#include <aspect/geometry_model/interface.h>
2223

2324
namespace aspect
2425
{
@@ -37,8 +38,11 @@ namespace aspect
3738
Function<dim>::initialize_one_particle_property(const Point<dim> &position,
3839
std::vector<double> &data) const
3940
{
41+
// convert the position into the selected coordinate system
42+
const Utilities::NaturalCoordinate<dim> point = this->get_geometry_model().cartesian_to_other_coordinates(position, coordinate_system);
43+
4044
for (unsigned int i = 0; i < n_components; ++i)
41-
data.push_back(function->value(position, i));
45+
data.push_back(function->value(Utilities::convert_array_to_point<dim>(point.get_coordinates()), i));
4246
}
4347

4448
template <int dim>
@@ -56,6 +60,16 @@ namespace aspect
5660
{
5761
prm.enter_subsection("Function");
5862
{
63+
prm.declare_entry ("Coordinate system", "cartesian",
64+
Patterns::Selection ("cartesian|spherical|depth"),
65+
"A selection that determines the assumed coordinate "
66+
"system for the function variables. Allowed values "
67+
"are `cartesian', `spherical', and `depth'. `spherical' coordinates "
68+
"are interpreted as r,phi or r,phi,theta in 2d/3d "
69+
"respectively with theta being the polar angle. `depth' "
70+
"will create a function, in which only the first "
71+
"parameter is non-zero, which is interpreted to "
72+
"be the depth of the point.");
5973
prm.declare_entry ("Number of components", "1",
6074
Patterns::Integer (0),
6175
"The number of function components where each component is described "
@@ -71,20 +85,24 @@ namespace aspect
7185
Function<dim>::parse_parameters (ParameterHandler &prm)
7286
{
7387
prm.enter_subsection("Function");
74-
n_components = prm.get_integer ("Number of components");
75-
try
88+
{
89+
n_components = prm.get_integer ("Number of components");
90+
try
7691
{
7792
function = std::make_unique<Functions::ParsedFunction<dim>>(n_components);
7893
function->parse_parameters (prm);
7994
}
80-
catch (...)
95+
catch (...)
8196
{
8297
std::cerr << "ERROR: FunctionParser failed to parse\n"
83-
<< "\t'Particles.Function'\n"
84-
<< "with expression\n"
85-
<< "\t'" << prm.get("Function expression") << "'";
98+
<< "\t'Particles.Function'\n"
99+
<< "with expression\n"
100+
<< "\t'" << prm.get("Function expression") << "'";
86101
throw;
87102
}
103+
104+
coordinate_system = Utilities::Coordinates::string_to_coordinate_system(prm.get("Coordinate system"));
105+
}
88106
prm.leave_subsection();
89107
}
90108
}

‎tests/particle_property_multiple_functions_spherical_coordinates.prm‎

Lines changed: 118 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,118 @@
1+
# This test is identical to particle_property_multiple_functions.prm
2+
# but with spherical coordinates.
3+
#
4+
# A description of the van Keken et al. benchmark using smooth initial
5+
# conditions instead of the discontinuous ones in van-keken-discontinuous.prm.
6+
# See the manual for more information. This test adds particles to the cookbook
7+
# parameter file with two particle properties defined by two different functions.
8+
# One is a discontinuous function and the other is a smooth function.
9+
10+
# MPI: 2
11+
12+
set Dimension = 2
13+
set End time = 70
14+
set Use years in output instead of seconds = false
15+
16+
subsection Geometry model
17+
set Model name = box
18+
19+
subsection Box
20+
set X extent = 0.9142
21+
set Y extent = 1.0000
22+
end
23+
end
24+
25+
subsection Boundary velocity model
26+
set Tangential velocity boundary indicators = left, right
27+
set Zero velocity boundary indicators = bottom, top
28+
end
29+
30+
subsection Material model
31+
set Model name = simple
32+
33+
subsection Simple model
34+
set Reference density = 1010
35+
set Viscosity = 1e2
36+
set Thermal expansion coefficient = 0
37+
set Density differential for compositional field 1 = -10
38+
end
39+
end
40+
41+
subsection Gravity model
42+
set Model name = vertical
43+
44+
subsection Vertical
45+
set Magnitude = 10
46+
end
47+
end
48+
49+
############### Parameters describing the temperature field
50+
# Note: The temperature plays no role in this model
51+
52+
subsection Initial temperature model
53+
set Model name = function
54+
55+
subsection Function
56+
set Function expression = 0
57+
end
58+
end
59+
60+
############### Parameters describing the compositional field
61+
# Note: The compositional field is what drives the flow
62+
# in this example
63+
64+
subsection Compositional fields
65+
set Number of fields = 1
66+
end
67+
68+
subsection Initial composition model
69+
set Model name = function
70+
71+
subsection Function
72+
set Variable names = x,z
73+
set Function constants = pi=3.1415926
74+
set Function expression = 0.5*(1+tanh((0.2+0.02*cos(pi*x/0.9142)-z)/0.02))
75+
end
76+
end
77+
78+
############### Parameters describing the discretization
79+
80+
subsection Mesh refinement
81+
set Initial adaptive refinement = 0
82+
set Strategy = composition
83+
set Initial global refinement = 4
84+
set Time steps between mesh refinement = 0
85+
set Coarsening fraction = 0.05
86+
set Refinement fraction = 0.3
87+
end
88+
89+
############### Parameters describing what to do with the solution
90+
91+
subsection Postprocess
92+
set List of postprocessors = velocity statistics, composition statistics, particles
93+
94+
subsection Particles
95+
set Time between data output = 70
96+
set Data output format = vtu
97+
end
98+
end
99+
100+
subsection Particles
101+
set List of particle properties = velocity, function, initial composition, initial position
102+
set Particle generator name = probability density function
103+
104+
subsection Function
105+
set Coordinate system = spherical
106+
set Number of components = 2
107+
set Variable names = x,z
108+
set Function expression = if( (z>0.2+0.02*cos(pi*x/0.9142)) , 0 , 1 ); 0.5*(1+tanh((0.2+0.02*cos(pi*x/0.9142)-z)/0.02))
109+
end
110+
111+
subsection Generator
112+
subsection Probability density function
113+
set Variable names = x,z
114+
set Function expression = x*x*z
115+
set Number of particles = 10
116+
end
117+
end
118+
end

‎tests/particle_property_multiple_functions_spherical_coordinates/particles.pvd‎

Lines changed: 10 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,10 @@
1+
<?xml version="1.0"?>
2+
<!--
3+
#This file was generated by the deal.II library on 2016/9/20 at 14:39:51
4+
-->
5+
<VTKFile type="Collection" version="0.1" ByteOrder="LittleEndian">
6+
<Collection>
7+
<DataSet timestep="0" group="" part="0" file="particles/particles-00000.pvtu"/>
8+
<DataSet timestep="70" group="" part="0" file="particles/particles-00001.pvtu"/>
9+
</Collection>
10+
</VTKFile>

‎tests/particle_property_multiple_functions_spherical_coordinates/particles/particles-00000.0000.vtu‎

Lines changed: 47 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,47 @@
1+
<?xml version="1.0" ?>
2+
<!--
3+
# vtk DataFile Version 3.0
4+
#This file was generated
5+
<VTKFile type="UnstructuredGrid" version="0.1" compressor="vtkZLibDataCompressor" byte_order="LittleEndian">
6+
<UnstructuredGrid>
7+
<FieldData>
8+
<DataArray type="Float32" Name="TIME" NumberOfTuples="1" format="ascii">0</DataArray>
9+
</FieldData>
10+
<Piece NumberOfPoints="3" NumberOfCells="3" >
11+
<Points>
12+
<DataArray type="Float32" NumberOfComponents="3" format="binary">
13+
AQAAACQAAAAkAAAAJQAAAA==eAH7d0PWXrX8gh0DENzxSref7LUUzJ630cB+SfwzMBsA6zYLTg==
14+
</DataArray>
15+
</Points>
16+
17+
<Cells>
18+
<DataArray type="Int32" Name="connectivity" format="binary">
19+
AQAAAAwAAAAMAAAAEQAAAA==eAFjYGBgYARiJiAGAAAcAAQ=
20+
</DataArray>
21+
<DataArray type="Int32" Name="offsets" format="binary">
22+
AQAAAAwAAAAMAAAAEQAAAA==eAFjZGBgYAJiZiAGAAA0AAc=
23+
</DataArray>
24+
<DataArray type="UInt8" Name="types" format="binary">
25+
AQAAAAMAAAADAAAACwAAAA==eAFjZGQEAAAJAAQ=
26+
</DataArray>
27+
</Cells>
28+
<PointData Scalars="scalars">
29+
<DataArray type="Float32" Name="velocity" NumberOfComponents="3" format="binary">
30+
AQAAACQAAAAkAAAACwAAAA==eAFjYCAMAAAkAAE=
31+
</DataArray>
32+
<DataArray type="Float32" Name="function" NumberOfComponents="3" format="binary">
33+
AQAAACQAAAAkAAAAHwAAAA==eAFjYGBgmBN2WB9IgUFZhIYpjM2wMFMHxAYAWAkERg==
34+
</DataArray>
35+
<DataArray type="Float32" Name="initial position" NumberOfComponents="3" format="binary">
36+
AQAAACQAAAAkAAAAJQAAAA==eAH7d0PWXrX8gh0DENzxSref7LUUzJ630cB+SfwzMBsA6zYLTg==
37+
</DataArray>
38+
<DataArray type="Float32" Name="id" format="binary">
39+
AQAAAAwAAAAMAAAAEAAAAA==eAFjYACBBnsg4QAABIcBAA==
40+
</DataArray>
41+
<DataArray type="Float32" Name="initial C_1" format="binary">
42+
AQAAAAwAAAAMAAAAFAAAAA==eAGbE3ZYvyxCw5RhYaYOAB9FBEY=
43+
</DataArray>
44+
</PointData>
45+
</Piece>
46+
</UnstructuredGrid>
47+
</VTKFile>

‎tests/particle_property_multiple_functions_spherical_coordinates/particles/particles-00000.0001.vtu‎

Lines changed: 47 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,47 @@
1+
<?xml version="1.0" ?>
2+
<!--
3+
# vtk DataFile Version 3.0
4+
#This file was generated
5+
<VTKFile type="UnstructuredGrid" version="0.1" compressor="vtkZLibDataCompressor" byte_order="LittleEndian">
6+
<UnstructuredGrid>
7+
<FieldData>
8+
<DataArray type="Float32" Name="TIME" NumberOfTuples="1" format="ascii">0</DataArray>
9+
</FieldData>
10+
<Piece NumberOfPoints="7" NumberOfCells="7" >
11+
<Points>
12+
<DataArray type="Float32" NumberOfComponents="3" format="binary">
13+
AQAAAFQAAABUAAAASAAAAA==eAFjn6Rm77LF2p4BCFiffbNbrl0EZnflcNjbs5WA2cXxCfYr2z3B7MioSHvRc6Fg9rKgAHvuUylg9qK+KPt3E+rAbACDXRU2
14+
</DataArray>
15+
</Points>
16+
17+
<Cells>
18+
<DataArray type="Int32" Name="connectivity" format="binary">
19+
AQAAABwAAAAcAAAAGwAAAA==eAFjYGBgYARiJiBmBmIWIGYFYjYgBgAA/AAW
20+
</DataArray>
21+
<DataArray type="Int32" Name="offsets" format="binary">
22+
AQAAABwAAAAcAAAAGwAAAA==eAFjZGBgYAJiZiBmAWJWIGYDYnYgBgABbAAd
23+
</DataArray>
24+
<DataArray type="UInt8" Name="types" format="binary">
25+
AQAAAAcAAAAHAAAACwAAAA==eAFjZAQDAAAjAAg=
26+
</DataArray>
27+
</Cells>
28+
<PointData Scalars="scalars">
29+
<DataArray type="Float32" Name="velocity" NumberOfComponents="3" format="binary">
30+
AQAAAFQAAABUAAAADAAAAA==eAFjYKA+AAAAVAAB
31+
</DataArray>
32+
<DataArray type="Float32" Name="function" NumberOfComponents="3" format="binary">
33+
AQAAAFQAAABUAAAADAAAAA==eAFjYKA+AAAAVAAB
34+
</DataArray>
35+
<DataArray type="Float32" Name="initial position" NumberOfComponents="3" format="binary">
36+
AQAAAFQAAABUAAAASAAAAA==eAFjn6Rm77LF2p4BCFiffbNbrl0EZnflcNjbs5WA2cXxCfYr2z3B7MioSHvRc6Fg9rKgAHvuUylg9qK+KPt3E+rAbACDXRU2
37+
</DataArray>
38+
<DataArray type="Float32" Name="id" format="binary">
39+
AQAAABwAAAAcAAAAHgAAAA==eAFjYHBwYGBoAOIFQHwAiB8AMYMjA4OAIwBHAgTT
40+
</DataArray>
41+
<DataArray type="Float32" Name="initial C_1" format="binary">
42+
AQAAABwAAAAcAAAACwAAAA==eAFjYMANAAAcAAE=
43+
</DataArray>
44+
</PointData>
45+
</Piece>
46+
</UnstructuredGrid>
47+
</VTKFile>

‎tests/particle_property_multiple_functions_spherical_coordinates/particles/particles-00000.pvtu‎

Lines changed: 19 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,19 @@
1+
<?xml version="1.0"?>
2+
<!--
3+
#This file was generated
4+
<VTKFile type="PUnstructuredGrid" version="0.1" byte_order="LittleEndian">
5+
<PUnstructuredGrid GhostLevel="0">
6+
<PPointData Scalars="scalars">
7+
<PDataArray type="Float32" Name="velocity" NumberOfComponents="3" format="ascii"/>
8+
<PDataArray type="Float32" Name="function" NumberOfComponents="3" format="ascii"/>
9+
<PDataArray type="Float32" Name="initial position" NumberOfComponents="3" format="ascii"/>
10+
<PDataArray type="Float32" Name="id" format="ascii"/>
11+
<PDataArray type="Float32" Name="initial C_1" format="ascii"/>
12+
</PPointData>
13+
<PPoints>
14+
<PDataArray type="Float32" NumberOfComponents="3"/>
15+
</PPoints>
16+
<Piece Source="particles-00000.0000.vtu"/>
17+
<Piece Source="particles-00000.0001.vtu"/>
18+
</PUnstructuredGrid>
19+
</VTKFile>

‎tests/particle_property_multiple_functions_spherical_coordinates/particles/particles-00001.0000.vtu‎

Lines changed: 48 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,48 @@
1+
<?xml version="1.0" ?>
2+
<!--
3+
# vtk DataFile Version 3.0
4+
#This file was generated
5+
<VTKFile type="UnstructuredGrid" version="0.1" compressor="vtkZLibDataCompressor" byte_order="LittleEndian">
6+
<UnstructuredGrid>
7+
<FieldData>
8+
<DataArray type="Float32" Name="CYCLE" NumberOfTuples="1" format="ascii">1</DataArray>
9+
<DataArray type="Float32" Name="TIME" NumberOfTuples="1" format="ascii">70</DataArray>
10+
</FieldData>
11+
<Piece NumberOfPoints="3" NumberOfCells="3" >
12+
<Points>
13+
<DataArray type="Float32" NumberOfComponents="3" format="binary">
14+
AQAAACQAAAAkAAAAJQAAAA==eAHTuS5nP0v1tB0DEDxyTrcva5wDZk9/bGh/5Ot9MBsA2bMLhQ==
15+
</DataArray>
16+
</Points>
17+
18+
<Cells>
19+
<DataArray type="Int32" Name="connectivity" format="binary">
20+
AQAAAAwAAAAMAAAAEQAAAA==eAFjYGBgYARiJiAGAAAcAAQ=
21+
</DataArray>
22+
<DataArray type="Int32" Name="offsets" format="binary">
23+
AQAAAAwAAAAMAAAAEQAAAA==eAFjZGBgYAJiZiAGAAA0AAc=
24+
</DataArray>
25+
<DataArray type="UInt8" Name="types" format="binary">
26+
AQAAAAMAAAADAAAACwAAAA==eAFjZGQEAAAJAAQ=
27+
</DataArray>
28+
</Cells>
29+
<PointData Scalars="scalars">
30+
<DataArray type="Float32" Name="velocity" NumberOfComponents="3" format="binary">
31+
AQAAACQAAAAkAAAAJQAAAA==eAGbbB1qMWWB7E4GIPDKZ91WMrsezF4m1GKxKN0OzAYA1EsKlQ==
32+
</DataArray>
33+
<DataArray type="Float32" Name="function" NumberOfComponents="3" format="binary">
34+
AQAAACQAAAAkAAAAHwAAAA==eAFjYGBgmBN2WB9IgUFZhIYpjM2wMFMHxAYAWAkERg==
35+
</DataArray>
36+
<DataArray type="Float32" Name="initial position" NumberOfComponents="3" format="binary">
37+
AQAAACQAAAAkAAAAJQAAAA==eAH7d0PWXrX8gh0DENzxSref7LUUzJ630cB+SfwzMBsA6zYLTg==
38+
</DataArray>
39+
<DataArray type="Float32" Name="id" format="binary">
40+
AQAAAAwAAAAMAAAAEAAAAA==eAFjYACBBnsg4QAABIcBAA==
41+
</DataArray>
42+
<DataArray type="Float32" Name="initial C_1" format="binary">
43+
AQAAAAwAAAAMAAAAFAAAAA==eAGbE3ZYvyxCw5RhYaYOAB9FBEY=
44+
</DataArray>
45+
</PointData>
46+
</Piece>
47+
</UnstructuredGrid>
48+
</VTKFile>

0 commit comments

Comments
 (0)