Peano 4
Loading...
Searching...
No Matches
noh.py
Go to the documentation of this file.
1import os
2import argparse
3import numpy as np
4
5import peano4
6import swift2
7import swift2.sphtools as sphtools
8import dastgen2
9
10import shutil
11
12
13
16
17parser = argparse.ArgumentParser(description="SPH test script")
18parser.add_argument(
19 "-j", "--parallel-builds", dest="j", type=int, default=-1, help="Parallel builds"
20)
21parser.add_argument(
22 "-v",
23 "--verbose",
24 dest="verbose",
25 action="store_true",
26 default=False,
27 help="Verbose",
28)
29parser.add_argument(
30 "-c",
31 "--clean",
32 dest="clean_temporary_directories",
33 action="store_true",
34 help="Remove temporary directories before script starts",
35)
36parser.add_argument(
37 "-np",
38 "--particle-number",
39 dest="particle_number",
40 required=True,
41 type=int,
42 help="Particle number (1D)",
43)
44parser.add_argument(
45 "-et",
46 "--end-time",
47 dest="end_time",
48 type=float,
49 default=0.2,
50 help="End of simulation",
51)
52parser.add_argument(
53 "-plot",
54 "--plot-delta",
55 dest="plot_delta",
56 type=float,
57 default=0.01,
58 help="Plot delta (0 switches off)",
59)
60parser.add_argument(
61 "-dt",
62 "--timestep-size",
63 dest="timestep_size",
64 type=float,
65 default=1e-4,
66 help="Global timestep size (initial value)",
67)
68parser.add_argument(
69 "-cfl",
70 "--cfl-factor",
71 dest="cfl_factor",
72 type=float,
73 default=0.1,
74 help="Value of the CFL constant",
75)
76parser.add_argument(
77 "-cs",
78 "--cell-size",
79 dest="cell_size",
80 type=float,
81 default=0.3,
82 help="Cell size (default 0.3)",
83)
84parser.add_argument(
85 "-pd",
86 "--peano-dir",
87 dest="peanodir",
88 default="../../../..",
89 help="Peano4 directory",
90)
91parser.add_argument(
92 "-asserts",
93 "--asserts",
94 dest="asserts",
95 action="store_true",
96 default=False,
97 help="Switch on assertions",
98)
99parser.add_argument(
100 "-k",
101 "--kernel",
102 dest="kernel",
103 type=str,
104 default="quartic_spline",
105 choices=sphtools.sph_kernel_list,
106 help="SPH kernel function to use.",
107)
108args = parser.parse_args()
109
110
111
114
115project_namespace = ["tests", "swift2", "noh2D"]
116
117project = swift2.Project(
118 namespace=project_namespace, project_name="Noh2D", executable="noh2D"
119)
120
121# Initialise the SPH particle --------------------------------------------------
122# If multiple particle species are present, we can specify the parameters
123# independently for each case.
124
125name = "HydroPart"
126
127dimensions = 2
128hydro_dimensions = 2 # dimensions for the hydro
129
131 name,
132 dimensions_hydro=hydro_dimensions,
133 cfl_factor=args.cfl_factor,
134 initial_time_step_size=args.timestep_size,
135 constant_time_step_size=True,
136 swift_project_namespace=project_namespace,
137 particles_per_cell=20,
138 min_h=args.cell_size,
139 max_h=args.cell_size,
140)
141
142
143
149
150offset = [0, 0, 0]
151domain_size = [2, 2, 2]
152periodic_boundary_conditions = [False, False, False]
153
154# HYDRO_DIMENSION can be smaller than the Peano4 dimensions
155HYDRO_DIMENSIONS = hydro_dimensions
156
157# Particle number along 1 dimension
158HYDRO_PART_NUMBER = args.particle_number
159
160# Time integration parameters
161GLOBAL_TIME_STEP_SIZE = args.timestep_size
162CFL_CONSTANT = args.cfl_factor
163
164# EoS parameter ---------------------------------------------------------------
165gamma_hydro_list = [str(5 / 3), str(7 / 5), str(4 / 3), str(2 / 1)]
166gamma_hydro_list_symbols = [
167 "HYDRO_GAMMA_5_3",
168 "HYDRO_GAMMA_7_5",
169 "HYDRO_GAMMA_4_3",
170 "HYDRO_GAMMA_2_1",
171]
172GAMMA = gamma_hydro_list[0]
173GAMMA_HYDRO_SYMBOL = gamma_hydro_list_symbols[0]
174
175# SPH kernel -------------------------------------------------------
176
177kernel = args.kernel
178KERNEL_SYMBOL = sphtools.sph_kernel_macro_name[kernel]
179KERNEL_SUPPORT_RAD = sphtools.sph_kernel_H_over_h[HYDRO_DIMENSIONS - 1][kernel]
180
181
182# Additional Particle Parameters -----------------------------------------------
183dx = domain_size[0] / args.particle_number
184hydro_H_max = 10.0 * dx
185hydro_h_max = swift2.sphtools.sph_kernel_get_smoothing_length_from_search_radius(
186 hydro_H_max, kernel, hydro_dimensions
187)
188
189particle.h_hydro_max = hydro_h_max
190
191# particle.h_hydro_min = 1e-6
192# particle.h_max_iterations = 10
193# particle.h_tolerance = 1e-4
194# particle.eta_factor = 1.2348
195# particle.beta_av = 3.0
196# particle.alpha_av = 0.8
197
198# Optional: Reduce the mantissa size for particle attributes -------------------
199# Native size is 53 (double precision)
200# particle.mantissa_size = 23
201
202
203# deploy particle properties
204particle.set_parameters()
205
206
207# Add particle to the SWIFT project --------------------------------------------
208# Notice that you any change to the particle will not take effect after adding
209# the particle_set to the project
210particle_set = project.add_particle_species(particle)
211
212if args.asserts:
213 build_mode = peano4.output.CompileMode.Asserts
214 # build_mode = peano4.output.CompileMode.Debug
215else:
216 build_mode = peano4.output.CompileMode.Release
217
218# Set Peano4 simulation parameters ---------------------------------------------
219project.set_global_simulation_parameters(
220 dimensions=dimensions,
221 offset=offset,
222 domain_size=domain_size,
223 min_end_time=args.end_time,
224 max_end_time=0.0,
225 first_plot_time_stamp=0.0,
226 time_in_between_plots=args.plot_delta,
227 periodic_BC=periodic_boundary_conditions,
228 plotter_precision=8,
229)
230
231project.set_Peano4_installation(args.peanodir, build_mode)
232
233project.set_load_balancing(
234 # "toolbox::loadbalancing::strategies::RecursiveBipartition",
235 "toolbox::loadbalancing::strategies::SpreadOutOnceGridStagnates",
236 "new toolbox::loadbalancing::DefaultConfiguration()",
237)
238
239
240
243
244# Noh problem (2D) snippet --------------------------------------------------
245# The particle attributes are set via an initalisation snippet that sets
246# the initial conditions specific for this problem, i.e. the velocity field
247# that converges towards the center of the domain and that generates the
248# collapse.
249
250ICs = swift2.scenario.NohProblemIC(dimensions=HYDRO_DIMENSIONS)
251initialisation_snippet = ICs.initialisation_snippet
252
253project.algorithm_step_initial_conditions.add_action_set(
255 particle_set=particle_set,
256 initialisation_call=initialisation_snippet,
257 distance_between_particles_along_one_axis=1.0 / args.particle_number,
258 additional_includes="""
259#include "swift2/kernels/legacy/equation_of_state.h"
260""",
261 )
262)
263
264
288 "particles",
289 particle_set,
290 "globaldata::{}::getSpecies().getMinTimeStamp()".format(particle.name),
291)
292particle_plotter.add_attribute_to_plot(particle.part_id, 1)
293particle_plotter.add_attribute_to_plot(particle.velocity, dimensions)
294particle_plotter.add_attribute_to_plot(particle.acceleration, dimensions)
295particle_plotter.add_attribute_to_plot(particle.mass, 1)
296particle_plotter.add_attribute_to_plot(particle.density, 1)
297particle_plotter.add_attribute_to_plot(particle.pressure, 1)
298particle_plotter.add_attribute_to_plot(particle.u, 1)
299particle_plotter.add_attribute_to_plot(particle.smoothL, 1)
300particle_plotter.add_attribute_to_plot(particle.is_boundary_part, 1)
301if args.asserts:
302 particle_plotter.add_attribute_to_plot(particle.density_neighbourcount, 1)
303 particle_plotter.add_attribute_to_plot(particle.force_neighbourcount, 1)
304project.algorithm_step_plot.add_action_set(particle_plotter)
305
306
307
318
319peano4_project = project.generate_Peano4_project(verbose=args.verbose)
320
321# Some sanity checks ----------------------------------------------------------
322if GAMMA not in gamma_hydro_list:
323 print(f"Please check the value of GAMMA. You have chosen: {GAMMA}")
324
325print("-----------------------------------------------------------------------")
326print("SPH parameters summary:")
327print(f"Effective dimensions of the Simulation: {HYDRO_DIMENSIONS}")
328print(f"Global time-step size: {GLOBAL_TIME_STEP_SIZE}")
329print(f"Number of particles (1D): {HYDRO_PART_NUMBER}")
330print(f"Adiabatic index: {GAMMA}")
331print(f"Smoothing length free factor: {particle.eta_factor}")
332print("-----------------------------------------------------------------------")
333
334
335# Export Constants.h file ------------------------------------------------------
336peano4_project.constants.export_const_with_type(
337 "HYDRO_DIMENSIONS", str(HYDRO_DIMENSIONS), "double"
338)
339peano4_project.constants.export_const_with_type(
340 "HYDRO_PART_NUMBER", str(HYDRO_PART_NUMBER), "double"
341)
342peano4_project.constants.export_const_with_type(
343 "GLOBAL_TIME_STEP_SIZE", str(GLOBAL_TIME_STEP_SIZE), "double"
344)
345peano4_project.constants.export_const_with_type("GAMMA", GAMMA, "double")
346
347# Export symbols for compilation (e.g. for SPH kernel) -----------------------
348peano4_project.output.makefile.add_CXX_flag(
349 "-D" + "HYDRO_DIMENSION=" + str(HYDRO_DIMENSIONS)
350)
351peano4_project.output.makefile.add_CXX_flag("-D" + GAMMA_HYDRO_SYMBOL)
352peano4_project.output.makefile.add_CXX_flag("-D" + KERNEL_SYMBOL)
353
354
355
358
359# If we want to start a project from scratch, we may want to remove some of the generated code
360if args.clean_temporary_directories:
361 dirs = ["vertexdata/", "repositories/", "globaldata/", "observers/"]
362 for dir in dirs:
363 shutil.rmtree(dir)
364
365# It is safe (and recommended) to overwrite any existing main file, since this
366# is generated automatically by the Swift2 graph compiler.
367peano4_project.generate(
368 overwrite=peano4.output.Overwrite.Always, throw_away_data_after_generation=True
369)
370peano4_project.build(make_clean_first=True, number_of_parallel_builds=args.j)
371
372# remove output files from previous runs
373output_files = [
374 f
375 for f in os.listdir(".")
376 if f.endswith(".peano-patch-file") or f.endswith(".vtu") or f.endswith(".pvd")
377]
378for f in output_files:
379 os.remove(f)
Swift2 project.
Definition Project.py:23
Initialisation snippet for the Noh problem in different spatial dimensions.
Definition NohProblem.py:5
This file is part of the SWIFT 2 project.
Definition __init__.py:1