Peano
Loading...
Searching...
No Matches
advection.py
Go to the documentation of this file.
1# This file is part of the ExaHyPE2 project. For conditions of distribution and
2# use, please see the copyright notice at www.peano-framework.org
3import sympy
4
5"""
6We import the required modules for this project.
7We always need the 'peano4' module as this is our core project.
8Since we are creating a SymHyPE/ExaHyPE 2 application, we additionally
9need to import the 'exahype2' and 'symhype' modules.
10"""
11import peano4
12import exahype2
13import symhype
14
15"""
16We specify the space dimensions here.
17We support either 2- or 3-dimensional problems.
18"""
19dimensions = 2
20
21"""
22The number of finite volumes per axis in one patch.
23"""
24patch_size = 16
25
26"""
27The number of levels the mesh is refined.
28"""
29depth = 5
30
31"""
32The simulation end time.
33"""
34end_time = 1.0
35
36"""
37Define the advection speed.
38"""
39advection_speed = 1.0
40
41"""
42Choose domain size and offset.
43"""
44size = [1.0, 1.0, 1.0]
45offset = [0.0, 0.0, 0.0]
46
47"""
48Choose how often a snapshot is written.
49"""
50time_in_between_two_snapshots = end_time / 10
51# time_in_between_two_snapshots = 0 # Set to 0 to disable I/O
52
53"""
54Switch between 'Release', 'Debug', 'Asserts', 'Trace', 'Stats'.
55"""
56compile_mode = "Release"
57
58"""
59We first create a new ExaHyPE 2 project.
60For this, we specify the (nested) namespaces, the name of our main file and our executable name.
61"""
62my_project = exahype2.Project(
63 namespace=["tutorials", "symhype", "advection"],
64 project_name="Advection",
65 directory=".",
66 executable="Advection",
67)
68
69"""
70Add the Finite Volumes solver using named arguments.
71This is the way you can add further PDE terms.
72This requires the 'BlockStructured' toolbox and 'ExaHyPE' to be built.
73"""
74my_solver = exahype2.solvers.fv.godunov.GlobalAdaptiveTimeStep(
75 name="FVSolver",
76 patch_size=patch_size,
77 unknowns=dimensions,
78 auxiliary_variables=0,
79 min_volume_h=(1.1 * min(size[0:dimensions]) / (3.0**depth)),
80 max_volume_h=(1.1 * min(size[0:dimensions]) / (3.0**depth)),
81 time_step_relaxation=0.5,
82)
83
84"""
85We want to define our PDE symbolically.
86"""
87my_pde = symhype.FirstOrderConservativePDEFormulation(
88 unknowns=dimensions, auxiliary_variables=0, dimensions=dimensions
89)
90v = my_pde.name_Q_entries(0, dimensions, "v")
91c = sympy.symbols("c")
92
93"""
94Define the equation system
95"""
96# Flux [unknowns, dimensions]
97if dimensions == 2:
98 my_pde.F[0, 0] = v[0]
99 my_pde.F[1, 0] = 0
100
101 my_pde.F[0, 1] = 0
102 my_pde.F[1, 1] = v[1]
103else:
104 my_pde.F[0, 0] = v[0]
105 my_pde.F[1, 0] = 0
106 my_pde.F[2, 0] = 0
107
108 my_pde.F[0, 1] = 0
109 my_pde.F[1, 1] = v[1]
110 my_pde.F[2, 1] = 0
111
112 my_pde.F[0, 2] = 0
113 my_pde.F[1, 2] = 0
114 my_pde.F[2, 2] = v[2]
115
116# Eigenvalues [unknowns, dimensions]
117if dimensions == 2:
118 my_pde.eigenvalues[0] = c
119 my_pde.eigenvalues[1] = c
120else:
121 my_pde.eigenvalues[0] = c
122 my_pde.eigenvalues[1] = c
123 my_pde.eigenvalues[2] = c
124
125my_pde.substitute_expression(c, advection_speed)
126
127"""
128Since 'my_pde' only holds the PDE without initial- or boundary conditions,
129we still need to properly define initial- and boundary conditions.
130This gives us then a complete description of a 'scenario'.
131"""
132
133# Initial conditions
134my_pde.initial_values[0] = my_pde.x[0]
135my_pde.initial_values[1] = my_pde.x[1]
136if dimensions == 3:
137 my_pde.initial_values[2] = my_pde.x[2]
138
139"""
140Specify which implementation our solvers uses.
141Here we want to set the implementation we get from our symbolically defined PDE,
142i.e., we get the C++ implementation which is generated by SymHyPE.
143"""
144my_solver.set_implementation(
145 initial_conditions=my_pde.implementation_of_initial_conditions(),
146 boundary_conditions=my_pde.implementation_of_homogeneous_Neumann_BC(),
147 flux=my_pde.implementation_of_flux(),
148 max_eigenvalue=my_pde.implementation_of_max_eigenvalue(),
149)
150
151"""
152To see which variables (unknowns + auxiliary variables) we can visualise,
153let's add a plot description for the variables to our solver.
154"""
155my_solver.plot_description = my_pde.unknown_identifier_for_plotter()
156
157"""
158Add the solver to our project
159"""
160my_project.add_solver(my_solver)
161
162"""
163Configure some global parameters
164"""
165my_project.set_global_simulation_parameters(
166 dimensions=dimensions,
167 size=size[0:dimensions],
168 offset=offset[0:dimensions],
169 min_end_time=end_time,
170 max_end_time=end_time,
171 first_plot_time_stamp=0.0,
172 time_in_between_plots=time_in_between_two_snapshots,
173 periodic_BC=[False, False, False],
174)
175
176"""
177This defines where the output files should go.
178If you omit this, output files are automatically put into the application's folder.
179"""
180my_project.set_output_path("solutions")
181
182"""
183Configure load balancer for parallel execution.
184"""
185my_project.set_load_balancer("new ::exahype2::LoadBalancingConfiguration")
186
187"""
188We need to set the location of our core libraries ('Peano4').
189This helps us to resolve any dependencies.
190Additionally, we specify the build mode which you can also change to a different mode.
191"""
192my_project.set_Peano4_installation(
193 "../../../", mode=peano4.output.string_to_mode(compile_mode)
194)
195
196"""
197We generate and grab the underlying core project of 'Peano4'.
198This gives us access to some functions we want to use to finalise and build this project.
199"""
200my_project = my_project.generate_Peano4_project(verbose=False)
201
202"""
203Finally, we want to build our project.
204First, all of the necessary glue code is generated in the application folder,
205then 'make' is invoked automatically which compiles the generated code and links against our core libraries
206and toolboxes which have been built before.
207You can also always invoke 'make' yourself to compile, or cleanup with 'make clean'.
208"""
209my_project.build(make=True, make_clean_first=True, throw_away_data_after_build=True)
210
211print("\nPDE is: ")
212print(my_pde.__str__())
213print(my_solver)
214print(my_project)