Note
Go to the end to download the full example code.
Custom warp disc simulation#
This simple example shows how to run a simulation of a disc with a custom initial warp
9 import numpy as np
10
11 import shamrock
12
13 # If we use the shamrock executable to run this script instead of the python interpreter,
14 # we should not initialize the system as the shamrock executable needs to handle specific MPI logic
15 if not shamrock.sys.is_initialized():
16 shamrock.change_loglevel(1)
17 shamrock.sys.init("0:0")
Setup units
23 si = shamrock.UnitSystem()
24 sicte = shamrock.Constants(si)
25 codeu = shamrock.UnitSystem(
26 unit_time=3600 * 24 * 365,
27 unit_length=sicte.au(),
28 unit_mass=sicte.sol_mass(),
29 )
30 ucte = shamrock.Constants(codeu)
31 G = ucte.G()
List parameters
37 # Resolution
38 Npart = 1000000
39
40 # Sink parameters
41 center_mass = 1.0
42 center_racc = 0.1
43
44 # Disc parameter
45 disc_mass = 0.01 # sol mass
46 rout = 60.0 # au
47 rin = 1.0 # au
48 H_r_0 = 0.05
49 q = 0.5
50 p = 3.0 / 2.0
51 r0 = 1.0
52
53 # Viscosity parameter
54 alpha_AV = 1.0e-3 / 0.08
55 alpha_u = 1.0
56 beta_AV = 2.0
57
58 # Integrator parameters
59 C_cour = 0.3
60 C_force = 0.25
61
62 # Warp parameters
63 inclination = 30.0
64
65
66 # Disc profiles
67 def sigma_profile(r):
68 sigma_0 = 1.0 # We do not care as it will be renormalized
69 return sigma_0 * (r / r0) ** (-p)
70
71
72 def kep_profile(r):
73 return (G * center_mass / r) ** 0.5
74
75
76 def omega_k(r):
77 return kep_profile(r) / r
78
79
80 def cs_profile(r):
81 cs_in = (H_r_0 * r0) * omega_k(r0)
82 return ((r / r0) ** (-q)) * cs_in
83
84
85 def inc_profile(r): # profile of the inclination angle
86 if r < 10.0:
87 effective_inc = 0.0
88 elif r < 30.0 and r > 10:
89 lx = 0.1 * (1.0 + np.sin((r - 20) * np.pi / 20.0))
90 lz = np.sqrt(1 - lx * lx)
91 effective_inc = np.arccos(lz)
92
93 else:
94 lx = 0.2
95 lz = np.sqrt(1 - lx * lx)
96 effective_inc = np.arccos(lz)
97
98 return effective_inc
99
100
101 def psi_profile(r):
102 return 0.0
103
104
105 def k_profile(r): # profile of the warp direction
106 incl_rad = inc_profile(r) * np.pi / 180.0
107 return (0.0, np.cos(incl_rad), np.sin(incl_rad))
Utility functions and quantities deduced from the base one
113 # Deduced quantities
114 pmass = disc_mass / Npart
115 bmin = (-rout * 2, -rout * 2, -rout * 2)
116 bmax = (rout * 2, rout * 2, rout * 2)
117
118 cs0 = cs_profile(r0)
119
120
121 def rot_profile(r):
122 return ((kep_profile(r) ** 2) - (2 * p + q) * cs_profile(r) ** 2) ** 0.5
123
124
125 def H_profile(r):
126 H = cs_profile(r) / omega_k(r)
127 # fact = (2.**0.5) * 3. # factor taken from phantom, to fasten thermalizing
128 fact = 1.0
129 return fact * H
Start the context The context holds the data of the code We then init the layout of the field (e.g. the list of fields used by the solver)
137 ctx = shamrock.Context()
138 ctx.pdata_layout_new()
Attach a SPH model to the data and configure it
143 model = shamrock.get_Model_SPH(context=ctx, vector_type="f64_3", sph_kernel="M4")
144
145 # Generate the default config
146 cfg = model.gen_default_config()
147 # Use disc alpha model viscosity
148 cfg.set_artif_viscosity_ConstantDisc(alpha_u=alpha_u, alpha_AV=alpha_AV, beta_AV=beta_AV)
149 # use the Lodato Price 2007 equation of state
150 cfg.set_eos_locally_isothermalLP07(cs0=cs0, q=q, r0=r0)
151 # Use the given code units
152 cfg.set_units(codeu)
153 # Change particle mass
154 cfg.set_particle_mass(pmass)
155 # Set the CFL
156 cfg.set_cfl_cour(C_cour)
157 cfg.set_cfl_force(C_force)
158
159 # Set the solver config to be the one stored in cfg
160 model.set_solver_config(cfg)
161
162 # Print the solver config
163 model.get_current_config().print_status()
164
165 # We want the patches to split above 10^8 part and merge if smaller than 1 part (e.g. disable patch)
166 model.init_scheduler(int(1e8), 1)
167
168 # Set the simulation box size
169 model.resize_simulation_box(bmin, bmax)
----- SPH Solver configuration -----
units :
unit_length : 149597870700
unit_mass : 1.98847e+30
unit_current : 1
unit_temperature : 1
unit_qte : 1
unit_lumint : 1
part mass 1e-08 ( can be changed using .set_part_mass() )
cfl force 0.25
cfl courant 0.3
--- artificial viscosity config
Config Type : constant disc
alpha_AV = 0.0125
alpha_u = 1
beta_AV = 2
--- artificial viscosity config (deduced)
-------------
EOS config f64_3 :
locally isothermal (Lodato Price 2007) :
--- Bondaries config
Config Type : Free boundaries
--- Bondaries config config (deduced)
-------------
------------------------------------
Add the sink particle
174 # null position and velocity
175 model.add_sink(center_mass, (0, 0, 0), (0, 0, 0), center_racc)
Create the setup
180 setup = model.get_setup()
181 gen_disc = setup.make_generator_disc_mc(
182 part_mass=pmass,
183 disc_mass=disc_mass,
184 r_in=rin,
185 r_out=rout,
186 sigma_profile=sigma_profile,
187 H_profile=H_profile,
188 rot_profile=rot_profile,
189 cs_profile=cs_profile,
190 random_seed=666,
191 )
192
193 # apply the custom warp
194 warp = setup.make_modifier_custom_warp(
195 parent=gen_disc,
196 inc_profile=inc_profile,
197 psi_profile=psi_profile,
198 k_profile=k_profile,
199 )
200 # Print the dot graph of the setup
201 print(warp.get_dot())
digraph G {
rankdir=LR;
node_0 [label="GeneratorMCDisc"];
node_1 [label="ApplyCustomWarp"];
node_0 -> node_1;
node_3 [label="Simulation"];
node_1 -> node_3;
}
Apply the setup
205 setup.apply_setup(warp)
Info: pushing data in scheduler, N = 1000000 [DataInserterUtility][rank=0]
Info: reattributing data ... [DataInserterUtility][rank=0]
Info: reattributing data done in 53.53 ms [DataInserterUtility][rank=0]
Info: Compute load ... [DataInserterUtility][rank=0]
Info: run scheduler step ... [DataInserterUtility][rank=0]
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 17.04 us (36.8%)
Info: summary : [LoadBalance][rank=0]
Info: - strategy "psweep" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: - strategy "round robin" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: Loadbalance stats : [LoadBalance][rank=0]
npatch = 1
min = 1000000
max = 1000000
avg = 1000000
efficiency = 100.00%
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 1914.00 ns (0.3%)
patch tree reduce : 1523.00 ns (0.3%)
gen split merge : 771.00 ns (0.1%)
split / merge op : 0/0
apply split merge : 982.00 ns (0.2%)
LB compute : 581.56 us (97.9%)
LB move op cnt : 0
LB apply : 3.64 us (0.6%)
Info: the setup took : 8.115601336000001 s [SPH setup][rank=0]
Run a single step to init the integrator and smoothing length of the particles
Here the htolerance is the maximum factor of evolution of the smoothing length in each
Smoothing length iterations, increasing it affect the performance negatively but increse the
convergence rate of the smoothing length
this is why we increase it temporely to 1.3 before lowering it back to 1.1 (default value)
Note that both change_htolerance can be removed and it will work the same but would converge
more slowly at the first timestep
216 model.change_htolerance(1.3)
217 model.timestep()
218 model.change_htolerance(1.1)
Warning: .change_htolerance(val) is deprecated, [SPH][rank=0]
-> calling this is replaced internally by .change_htolerances(coarse=val, fine=min(val, 1.1))
see: https://shamrock-code.github.io/Shamrock/mkdocs/models/sph/smoothing_length_tolerance
---------------- t = 0, dt = 0 ----------------
Info: summary : [LoadBalance][rank=0]
Info: - strategy "psweep" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: - strategy "round robin" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: Loadbalance stats : [LoadBalance][rank=0]
npatch = 1
min = 1000000
max = 1000000
avg = 1000000
efficiency = 100.00%
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 31.47 us (3.2%)
patch tree reduce : 8.52 us (0.9%)
gen split merge : 741.00 ns (0.1%)
split / merge op : 0/0
apply split merge : 1142.00 ns (0.1%)
LB compute : 871.53 us (90.0%)
LB move op cnt : 0
LB apply : 3.46 us (0.4%)
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 2.84 us (68.2%)
Info: free boundaries skipping geometry update [PositionUpdated][rank=0]
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 0.6759353366995453 unconverged cnt = 1000000
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 0.8787159377094089 unconverged cnt = 1000000
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 1.1423307190222316 unconverged cnt = 1000000
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 1.4721686506842913 unconverged cnt = 999999
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 1.9043892021115991 unconverged cnt = 999999
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 2.475705962745079 unconverged cnt = 999998
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 3.2184177515686025 unconverged cnt = 999991
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 4.183943077039183 unconverged cnt = 999968
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 4.328497738553627 unconverged cnt = 999908
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 4.534227032083806 unconverged cnt = 999610
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 4.534227032084451 unconverged cnt = 996367
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 4.534227032084451 unconverged cnt = 952305
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 4.534227032084451 unconverged cnt = 723608
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 4.534227032084451 unconverged cnt = 135319
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 4.534227032084451 unconverged cnt = 679
Info: conservation infos : [sph::Model][rank=0]
sum v = (1.0280937103436009e-05,-1.2210939977754751e-05,1.5490724643677914e-06)
sum a = (-2.236166980751353e-19,1.8295911660692887e-19,-2.0265925007549178e-21)
sum e = 0.013574918344084199
sum de = 1.0594529285579178e-20
Info: cfl dt = 6.670106715351201e-05 cfl multiplier : 0.01 [sph::Model][rank=0]
Info: processing rate infos : [sph::Model][rank=0]
---------------------------------------------------------------------------------------
| rank | rate (N.s^-1) | Nobj | t compute (s) | MPI | alloc | mem (max) |
---------------------------------------------------------------------------------------
| 0 | 3.3281e+04 | 1000000 | 3.005e+01 | 0 % | 0 % | 1.00 GB |
---------------------------------------------------------------------------------------
Info: estimated rate : 0 (tsim/hr) [sph::Model][rank=0]
Warning: .change_htolerance(val) is deprecated, [SPH][rank=0]
-> calling this is replaced internally by .change_htolerances(coarse=val, fine=min(val, 1.1))
see: https://shamrock-code.github.io/Shamrock/mkdocs/models/sph/smoothing_length_tolerance
Manipulating the simulation#
Dump files (path relative to where you have started shamrock)
225 dump_folder = "_to_trash"
226 import os
227
228 os.system("mkdir -p " + dump_folder)
229
230 # VTK dump
231 model.do_vtk_dump(dump_folder + "/init_disc.vtk", True)
232
233 # Shamrock restart dump files
234 model.dump(dump_folder + "/init_disc.sham")
235
236 # Phantom dump
237 dump = model.make_phantom_dump()
238 dump.save_dump(dump_folder + "/init_disc.phdump")
Info: dump to _to_trash/init_disc.vtk [VTK Dump][rank=0]
- took 143.84 ms, bandwidth = 371.28 MB/s
Info: Dumping state to _to_trash/init_disc.sham [SPH][rank=0]
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 14.92 us (31.4%)
Info: dump to _to_trash/init_disc.sham [Shamrock Dump][rank=0]
- took 129.41 ms, bandwidth = 884.40 MB/s
Single timestep
242 model.evolve_once()
---------------- t = 0, dt = 6.670106715351201e-05 ----------------
Info: summary : [LoadBalance][rank=0]
Info: - strategy "psweep" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: - strategy "round robin" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: Loadbalance stats : [LoadBalance][rank=0]
npatch = 1
min = 1000000
max = 1000000
avg = 1000000
efficiency = 100.00%
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 22.10 us (2.3%)
patch tree reduce : 3.71 us (0.4%)
gen split merge : 551.00 ns (0.1%)
split / merge op : 0/0
apply split merge : 1122.00 ns (0.1%)
LB compute : 880.77 us (93.0%)
LB move op cnt : 0
LB apply : 3.43 us (0.4%)
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 3.06 us (71.9%)
Info: free boundaries skipping geometry update [PositionUpdated][rank=0]
Info: conservation infos : [sph::Model][rank=0]
sum v = (1.0279926159979139e-05,-1.2210250385901987e-05,1.5490729563133356e-06)
sum a = (-1.0299920638612292e-18,-6.606856988583543e-19,3.871536138993191e-20)
sum e = 0.013574918539838787
sum de = 4.368034739352948e-08
Info: cfl dt = 0.002267972213422564 cfl multiplier : 0.34 [sph::Model][rank=0]
Info: processing rate infos : [sph::Model][rank=0]
---------------------------------------------------------------------------------------
| rank | rate (N.s^-1) | Nobj | t compute (s) | MPI | alloc | mem (max) |
---------------------------------------------------------------------------------------
| 0 | 1.8587e+05 | 1000000 | 5.380e+00 | 0 % | 0 % | 1.00 GB |
---------------------------------------------------------------------------------------
Info: estimated rate : 0.04463134510163855 (tsim/hr) [sph::Model][rank=0]
Info: iteration since start : 2 [SPH][rank=0]
Info: time since start : 166.02977888400002 (s) [SPH][rank=0]
Evolve until a given time (code units)
246 model.evolve_until(0.001)
---------------- t = 6.670106715351201e-05, dt = 0.000933298932846488 ----------------
Info: summary : [LoadBalance][rank=0]
Info: - strategy "psweep" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: - strategy "round robin" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: Loadbalance stats : [LoadBalance][rank=0]
npatch = 1
min = 1000000
max = 1000000
avg = 1000000
efficiency = 100.00%
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 26.17 us (2.6%)
patch tree reduce : 4.97 us (0.5%)
gen split merge : 521.00 ns (0.1%)
split / merge op : 0/0
apply split merge : 1323.00 ns (0.1%)
LB compute : 950.49 us (93.9%)
LB move op cnt : 0
LB apply : 3.41 us (0.3%)
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 3.59 us (75.4%)
Info: free boundaries skipping geometry update [PositionUpdated][rank=0]
Info: conservation infos : [sph::Model][rank=0]
sum v = (1.0266790796449059e-05,-1.2201295869683673e-05,1.5490793478006986e-06)
sum a = (-1.5924219408380846e-19,1.5009423825346202e-18,-8.43939291763352e-20)
sum e = 0.0135749567116035
sum de = 6.545898343482271e-07
Info: cfl dt = 0.003738871556162255 cfl multiplier : 0.56 [sph::Model][rank=0]
Info: processing rate infos : [sph::Model][rank=0]
---------------------------------------------------------------------------------------
| rank | rate (N.s^-1) | Nobj | t compute (s) | MPI | alloc | mem (max) |
---------------------------------------------------------------------------------------
| 0 | 1.7557e+05 | 1000000 | 5.696e+00 | 0 % | 0 % | 1.00 GB |
---------------------------------------------------------------------------------------
Info: estimated rate : 0.5899025131217444 (tsim/hr) [sph::Model][rank=0]
Info: iteration since start : 3 [SPH][rank=0]
Info: time since start : 171.983172661 (s) [SPH][rank=0]
True
Get the sinks positions
250 print(model.get_sinks())
[{'pos': (1.515725044321675e-11, -1.0333699921380228e-11, -7.375378607462444e-15), 'velocity': (3.0327272966099705e-08, -2.0599291006627177e-08, -1.4750954840409288e-11), 'sph_acceleration': (3.034186124576687e-05, -2.0521496431259187e-05, -1.4751180698426988e-08), 'ext_acceleration': (0.0, 0.0, 0.0), 'mass': 1.0, 'angular_momentum': (0.0, 0.0, 0.0), 'accretion_radius': 0.1}]
Get the fields as python dictionary of numpy arrays
Warning
Do not do this on a large distributed simulation as this gather all data on MPI rank 0 and will use a lot of memory (and crash if the simulation is too large)
260 print(ctx.collect_data())
Info: collected : 1 patches [PatchScheduler][rank=0]
adding -> xyz
adding -> vxyz
adding -> axyz
adding -> axyz_ext
adding -> hpart
adding -> uint
adding -> duint
{'xyz': array([[-47.63546795, -36.78576494, -0.28681528],
[-34.77056074, -49.22238153, -1.92203334],
[-37.11732682, -46.91115851, -0.87803673],
...,
[ 34.60829161, 46.90640238, 2.01249056],
[ 36.49163377, 47.83425436, 0.38794924],
[ 36.11418584, 46.5831198 , 4.33455869]], shape=(1000000, 3)), 'vxyz': array([[ 0.49239165, -0.63683688, -0.10046196],
[ 0.65451934, -0.45713655, -0.13356885],
[ 0.63003086, -0.49609147, -0.12856823],
...,
[-0.65645734, 0.47859752, 0.13396347],
[-0.63558971, 0.48382671, 0.12970318],
[-0.64631718, 0.48879407, 0.13189137]], shape=(1000000, 3)), 'axyz': array([[ 8.18978867e-03, 6.22404543e-03, -1.40943368e-03],
[ 5.48414381e-03, 8.06454778e-03, -7.60167061e-04],
[ 6.60823869e-03, 7.75994499e-03, -1.84601310e-03],
...,
[-6.31798362e-03, -9.40736223e-03, 1.51289461e-03],
[-5.41499475e-03, -7.45433611e-03, 9.81209327e-04],
[-6.55516559e-03, -8.86301318e-03, 3.56114421e-05]],
shape=(1000000, 3)), 'axyz_ext': array([[ 8.61381241e-03, 6.65188550e-03, 5.18641492e-05],
[ 6.25348519e-03, 8.85264510e-03, 3.45677687e-04],
[ 6.83397636e-03, 8.63719927e-03, 1.61662565e-04],
...,
[-6.87609049e-03, -9.31952004e-03, -3.99848319e-04],
[-6.60555418e-03, -8.65874520e-03, -7.02248566e-05],
[-6.89667804e-03, -8.89591643e-03, -8.27764908e-04]],
shape=(1000000, 3)), 'hpart': array([2.8826018 , 2.76533743, 2.20585668, ..., 2.29494582, 2.18852028,
3.9705829 ], shape=(1000000,)), 'uint': array([-2.93846067e-09, 1.70135720e-10, 1.52133429e-09, ...,
2.79567483e-09, -1.70947890e-09, 1.90116364e-09],
shape=(1000000,)), 'duint': array([-2.93865947e-06, 1.69970420e-07, 1.52031272e-06, ...,
2.79501149e-06, -1.71007070e-06, 1.90162506e-06],
shape=(1000000,))}
Performing a timestep loop
264 dt_stop = 0.001
265 for i in range(10):
266 t_target = i * dt_stop
267 # skip if the model is already past the target
268 if model.get_time() > t_target:
269 continue
270
271 model.evolve_until(i * dt_stop)
272
273 # Dump name is "dump_xxxx.sham" where xxxx is the timestep
274 model.dump(dump_folder + f"/dump_{i:04}.sham")
Info: iteration since start : 3 [SPH][rank=0]
Info: time since start : 172.92844971800002 (s) [SPH][rank=0]
Info: Dumping state to _to_trash/dump_0001.sham [SPH][rank=0]
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 46.88 us (62.6%)
Info: dump to _to_trash/dump_0001.sham [Shamrock Dump][rank=0]
- took 129.55 ms, bandwidth = 883.42 MB/s
---------------- t = 0.001, dt = 0.001 ----------------
Info: summary : [LoadBalance][rank=0]
Info: - strategy "psweep" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: - strategy "round robin" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: Loadbalance stats : [LoadBalance][rank=0]
npatch = 1
min = 1000000
max = 1000000
avg = 1000000
efficiency = 100.00%
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 7.80 us (1.2%)
patch tree reduce : 3.43 us (0.5%)
gen split merge : 712.00 ns (0.1%)
split / merge op : 0/0
apply split merge : 1122.00 ns (0.2%)
LB compute : 603.62 us (96.3%)
LB move op cnt : 0
LB apply : 3.73 us (0.6%)
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 3.03 us (72.1%)
Info: free boundaries skipping geometry update [PositionUpdated][rank=0]
Info: conservation infos : [sph::Model][rank=0]
sum v = (1.0265766172814262e-05,-1.2200679229537624e-05,1.5490798399582005e-06)
sum a = (-1.2705494208814505e-18,-8.199278929421627e-19,-3.970301504132022e-20)
sum e = 0.013574963320910426
sum de = 1.3096156907043732e-06
Info: cfl dt = 0.004723266764344826 cfl multiplier : 0.7066666666666667 [sph::Model][rank=0]
Info: processing rate infos : [sph::Model][rank=0]
---------------------------------------------------------------------------------------
| rank | rate (N.s^-1) | Nobj | t compute (s) | MPI | alloc | mem (max) |
---------------------------------------------------------------------------------------
| 0 | 1.7030e+05 | 1000000 | 5.872e+00 | 0 % | 0 % | 1.00 GB |
---------------------------------------------------------------------------------------
Info: estimated rate : 0.6130946054795301 (tsim/hr) [sph::Model][rank=0]
Info: iteration since start : 4 [SPH][rank=0]
Info: time since start : 178.94952252200002 (s) [SPH][rank=0]
Info: Dumping state to _to_trash/dump_0002.sham [SPH][rank=0]
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 7.11 us (57.1%)
Info: dump to _to_trash/dump_0002.sham [Shamrock Dump][rank=0]
- took 126.73 ms, bandwidth = 903.09 MB/s
---------------- t = 0.002, dt = 0.001 ----------------
Info: summary : [LoadBalance][rank=0]
Info: - strategy "psweep" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: - strategy "round robin" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: Loadbalance stats : [LoadBalance][rank=0]
npatch = 1
min = 1000000
max = 1000000
avg = 1000000
efficiency = 100.00%
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 7.78 us (1.4%)
patch tree reduce : 1202.00 ns (0.2%)
gen split merge : 711.00 ns (0.1%)
split / merge op : 0/0
apply split merge : 1162.00 ns (0.2%)
LB compute : 515.06 us (96.0%)
LB move op cnt : 0
LB apply : 3.37 us (0.6%)
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 2.81 us (72.2%)
Info: free boundaries skipping geometry update [PositionUpdated][rank=0]
Info: conservation infos : [sph::Model][rank=0]
sum v = (1.0265751934751674e-05,-1.2200757015636397e-05,1.549079840565865e-06)
sum a = (-1.4568966692773966e-19,1.8295911660692887e-19,-3.336681154916403e-20)
sum e = 0.013574964954389044
sum de = 1.9651311760317477e-06
Info: cfl dt = 0.005383384395416449 cfl multiplier : 0.8044444444444444 [sph::Model][rank=0]
Info: processing rate infos : [sph::Model][rank=0]
---------------------------------------------------------------------------------------
| rank | rate (N.s^-1) | Nobj | t compute (s) | MPI | alloc | mem (max) |
---------------------------------------------------------------------------------------
| 0 | 1.7855e+05 | 1000000 | 5.601e+00 | 0 % | 0 % | 1.00 GB |
---------------------------------------------------------------------------------------
Info: estimated rate : 0.6427793121536948 (tsim/hr) [sph::Model][rank=0]
Info: iteration since start : 5 [SPH][rank=0]
Info: time since start : 184.70098212300002 (s) [SPH][rank=0]
Info: Dumping state to _to_trash/dump_0003.sham [SPH][rank=0]
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 8.58 us (59.8%)
Info: dump to _to_trash/dump_0003.sham [Shamrock Dump][rank=0]
- took 126.48 ms, bandwidth = 904.84 MB/s
---------------- t = 0.003, dt = 0.001 ----------------
Info: summary : [LoadBalance][rank=0]
Info: - strategy "psweep" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: - strategy "round robin" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: Loadbalance stats : [LoadBalance][rank=0]
npatch = 1
min = 1000000
max = 1000000
avg = 1000000
efficiency = 100.00%
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 7.32 us (0.6%)
patch tree reduce : 1273.00 ns (0.1%)
gen split merge : 691.00 ns (0.1%)
split / merge op : 0/0
apply split merge : 1091.00 ns (0.1%)
LB compute : 1192.24 us (98.2%)
LB move op cnt : 0
LB apply : 4.20 us (0.3%)
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 3.12 us (71.3%)
Info: free boundaries skipping geometry update [PositionUpdated][rank=0]
Info: conservation infos : [sph::Model][rank=0]
sum v = (1.0265738057867366e-05,-1.2200834791363192e-05,1.5490798415554962e-06)
sum a = (1.2095630486791409e-18,-1.7618285302889447e-19,-1.6095280359056812e-20)
sum e = 0.013574967243892004
sum de = 2.6210752343887347e-06
Info: cfl dt = 0.005827444918441602 cfl multiplier : 0.8696296296296296 [sph::Model][rank=0]
Info: processing rate infos : [sph::Model][rank=0]
---------------------------------------------------------------------------------------
| rank | rate (N.s^-1) | Nobj | t compute (s) | MPI | alloc | mem (max) |
---------------------------------------------------------------------------------------
| 0 | 1.8038e+05 | 1000000 | 5.544e+00 | 0 % | 0 % | 1.00 GB |
---------------------------------------------------------------------------------------
Info: estimated rate : 0.6493765323466556 (tsim/hr) [sph::Model][rank=0]
Info: iteration since start : 6 [SPH][rank=0]
Info: time since start : 190.390227366 (s) [SPH][rank=0]
Info: Dumping state to _to_trash/dump_0004.sham [SPH][rank=0]
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 9.07 us (61.3%)
Info: dump to _to_trash/dump_0004.sham [Shamrock Dump][rank=0]
- took 126.50 ms, bandwidth = 904.74 MB/s
---------------- t = 0.004, dt = 0.001 ----------------
Info: summary : [LoadBalance][rank=0]
Info: - strategy "psweep" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: - strategy "round robin" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: Loadbalance stats : [LoadBalance][rank=0]
npatch = 1
min = 1000000
max = 1000000
avg = 1000000
efficiency = 100.00%
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 8.05 us (1.7%)
patch tree reduce : 1202.00 ns (0.2%)
gen split merge : 671.00 ns (0.1%)
split / merge op : 0/0
apply split merge : 1353.00 ns (0.3%)
LB compute : 443.82 us (91.9%)
LB move op cnt : 0
LB apply : 21.66 us (4.5%)
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 2.62 us (67.8%)
Info: free boundaries skipping geometry update [PositionUpdated][rank=0]
Info: conservation infos : [sph::Model][rank=0]
sum v = (1.0265724541865484e-05,-1.2200912555107853e-05,1.5490798429141267e-06)
sum a = (-5.319366908757006e-19,-2.0328790734103208e-19,-1.9513190650125923e-20)
sum e = 0.013574970189721664
sum de = 3.2774348542771947e-06
Info: cfl dt = 0.006072868038638684 cfl multiplier : 0.9130864197530864 [sph::Model][rank=0]
Info: processing rate infos : [sph::Model][rank=0]
---------------------------------------------------------------------------------------
| rank | rate (N.s^-1) | Nobj | t compute (s) | MPI | alloc | mem (max) |
---------------------------------------------------------------------------------------
| 0 | 1.8069e+05 | 1000000 | 5.534e+00 | 0 % | 0 % | 1.00 GB |
---------------------------------------------------------------------------------------
Info: estimated rate : 0.6504769394943632 (tsim/hr) [sph::Model][rank=0]
Info: iteration since start : 7 [SPH][rank=0]
Info: time since start : 196.07406788300003 (s) [SPH][rank=0]
Info: Dumping state to _to_trash/dump_0005.sham [SPH][rank=0]
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 9.04 us (61.8%)
Info: dump to _to_trash/dump_0005.sham [Shamrock Dump][rank=0]
- took 127.26 ms, bandwidth = 899.34 MB/s
---------------- t = 0.005, dt = 0.001 ----------------
Info: summary : [LoadBalance][rank=0]
Info: - strategy "psweep" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: - strategy "round robin" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: Loadbalance stats : [LoadBalance][rank=0]
npatch = 1
min = 1000000
max = 1000000
avg = 1000000
efficiency = 100.00%
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 8.18 us (0.7%)
patch tree reduce : 1423.00 ns (0.1%)
gen split merge : 781.00 ns (0.1%)
split / merge op : 0/0
apply split merge : 1142.00 ns (0.1%)
LB compute : 1181.30 us (98.1%)
LB move op cnt : 0
LB apply : 3.51 us (0.3%)
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 2.94 us (74.7%)
Info: free boundaries skipping geometry update [PositionUpdated][rank=0]
Info: conservation infos : [sph::Model][rank=0]
sum v = (1.0265711386447873e-05,-1.2200990305259818e-05,1.5490798446295448e-06)
sum a = (-2.236166980751353e-18,1.6296913905172739e-18,-5.221164026434711e-21)
sum e = 0.01357497379226242
sum de = 3.934196921285916e-06
Info: cfl dt = 0.006202050139896026 cfl multiplier : 0.9420576131687243 [sph::Model][rank=0]
Info: processing rate infos : [sph::Model][rank=0]
---------------------------------------------------------------------------------------
| rank | rate (N.s^-1) | Nobj | t compute (s) | MPI | alloc | mem (max) |
---------------------------------------------------------------------------------------
| 0 | 1.8182e+05 | 1000000 | 5.500e+00 | 0 % | 0 % | 1.00 GB |
---------------------------------------------------------------------------------------
Info: estimated rate : 0.654563108231495 (tsim/hr) [sph::Model][rank=0]
Info: iteration since start : 8 [SPH][rank=0]
Info: time since start : 201.71865622 (s) [SPH][rank=0]
Info: Dumping state to _to_trash/dump_0006.sham [SPH][rank=0]
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 6.60 us (55.1%)
Info: dump to _to_trash/dump_0006.sham [Shamrock Dump][rank=0]
- took 129.02 ms, bandwidth = 887.03 MB/s
---------------- t = 0.006, dt = 0.001 ----------------
Info: summary : [LoadBalance][rank=0]
Info: - strategy "psweep" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: - strategy "round robin" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: Loadbalance stats : [LoadBalance][rank=0]
npatch = 1
min = 1000000
max = 1000000
avg = 1000000
efficiency = 100.00%
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 8.03 us (1.9%)
patch tree reduce : 1132.00 ns (0.3%)
gen split merge : 671.00 ns (0.2%)
split / merge op : 0/0
apply split merge : 1082.00 ns (0.3%)
LB compute : 411.08 us (95.1%)
LB move op cnt : 0
LB apply : 3.27 us (0.8%)
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 2.73 us (66.3%)
Info: free boundaries skipping geometry update [PositionUpdated][rank=0]
Info: conservation infos : [sph::Model][rank=0]
sum v = (1.0265698591307466e-05,-1.2201068040200577e-05,1.5490798466892382e-06)
sum a = (-5.55653613398821e-19,-7.453889935837843e-20,8.763151409386775e-21)
sum e = 0.013574978051887832
sum de = 4.591348607139268e-06
Info: cfl dt = 0.006268212111452408 cfl multiplier : 0.9613717421124829 [sph::Model][rank=0]
Info: processing rate infos : [sph::Model][rank=0]
---------------------------------------------------------------------------------------
| rank | rate (N.s^-1) | Nobj | t compute (s) | MPI | alloc | mem (max) |
---------------------------------------------------------------------------------------
| 0 | 1.8040e+05 | 1000000 | 5.543e+00 | 0 % | 0 % | 1.00 GB |
---------------------------------------------------------------------------------------
Info: estimated rate : 0.6494499119827504 (tsim/hr) [sph::Model][rank=0]
Info: iteration since start : 9 [SPH][rank=0]
Info: time since start : 207.409355766 (s) [SPH][rank=0]
Info: Dumping state to _to_trash/dump_0007.sham [SPH][rank=0]
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 8.45 us (61.0%)
Info: dump to _to_trash/dump_0007.sham [Shamrock Dump][rank=0]
- took 128.83 ms, bandwidth = 888.36 MB/s
---------------- t = 0.007, dt = 0.001 ----------------
Info: summary : [LoadBalance][rank=0]
Info: - strategy "psweep" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: - strategy "round robin" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: Loadbalance stats : [LoadBalance][rank=0]
npatch = 1
min = 1000000
max = 1000000
avg = 1000000
efficiency = 100.00%
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 7.94 us (0.6%)
patch tree reduce : 1282.00 ns (0.1%)
gen split merge : 701.00 ns (0.1%)
split / merge op : 0/0
apply split merge : 1092.00 ns (0.1%)
LB compute : 1255.88 us (98.2%)
LB move op cnt : 0
LB apply : 3.75 us (0.3%)
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 3.49 us (78.5%)
Info: free boundaries skipping geometry update [PositionUpdated][rank=0]
Info: conservation infos : [sph::Model][rank=0]
sum v = (1.0265686156132855e-05,-1.2201145758317953e-05,1.5490798490814083e-06)
sum a = (-8.131516293641283e-20,8.131516293641283e-20,3.181005763633923e-20)
sum e = 0.013574982968963246
sum de = 5.248876911121978e-06
Info: cfl dt = 0.0062940770393870886 cfl multiplier : 0.9742478280749886 [sph::Model][rank=0]
Info: processing rate infos : [sph::Model][rank=0]
---------------------------------------------------------------------------------------
| rank | rate (N.s^-1) | Nobj | t compute (s) | MPI | alloc | mem (max) |
---------------------------------------------------------------------------------------
| 0 | 1.8049e+05 | 1000000 | 5.540e+00 | 0 % | 0 % | 1.00 GB |
---------------------------------------------------------------------------------------
Info: estimated rate : 0.6497754237426541 (tsim/hr) [sph::Model][rank=0]
Info: iteration since start : 10 [SPH][rank=0]
Info: time since start : 213.096873046 (s) [SPH][rank=0]
Info: Dumping state to _to_trash/dump_0008.sham [SPH][rank=0]
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 7.22 us (57.2%)
Info: dump to _to_trash/dump_0008.sham [Shamrock Dump][rank=0]
- took 128.73 ms, bandwidth = 889.01 MB/s
---------------- t = 0.008, dt = 0.0010000000000000009 ----------------
Info: summary : [LoadBalance][rank=0]
Info: - strategy "psweep" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: - strategy "round robin" : max = 1000000 min = 1000000 [LoadBalance][rank=0]
Info: Loadbalance stats : [LoadBalance][rank=0]
npatch = 1
min = 1000000
max = 1000000
avg = 1000000
efficiency = 100.00%
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 7.41 us (0.6%)
patch tree reduce : 1653.00 ns (0.1%)
gen split merge : 671.00 ns (0.1%)
split / merge op : 0/0
apply split merge : 1042.00 ns (0.1%)
LB compute : 1160.57 us (98.1%)
LB move op cnt : 0
LB apply : 3.77 us (0.3%)
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 2.88 us (74.3%)
Info: free boundaries skipping geometry update [PositionUpdated][rank=0]
Info: conservation infos : [sph::Model][rank=0]
sum v = (1.0265674080609692e-05,-1.2201223457996384e-05,1.5490798517944498e-06)
sum a = (-4.743384504624082e-19,-5.692061405548898e-19,2.321730543313838e-20)
sum e = 0.01357498854384286
sum de = 5.906768970121784e-06
Info: cfl dt = 0.006294531947055687 cfl multiplier : 0.9828318853833258 [sph::Model][rank=0]
Info: processing rate infos : [sph::Model][rank=0]
---------------------------------------------------------------------------------------
| rank | rate (N.s^-1) | Nobj | t compute (s) | MPI | alloc | mem (max) |
---------------------------------------------------------------------------------------
| 0 | 1.7913e+05 | 1000000 | 5.583e+00 | 0 % | 0 % | 1.00 GB |
---------------------------------------------------------------------------------------
Info: estimated rate : 0.6448713447247272 (tsim/hr) [sph::Model][rank=0]
Info: iteration since start : 11 [SPH][rank=0]
Info: time since start : 218.82656540500003 (s) [SPH][rank=0]
Info: Dumping state to _to_trash/dump_0009.sham [SPH][rank=0]
Info: Scheduler step timings : [Scheduler][rank=0]
metadata sync : 7.86 us (58.5%)
Info: dump to _to_trash/dump_0009.sham [Shamrock Dump][rank=0]
- took 128.47 ms, bandwidth = 890.84 MB/s
Plot column integrated density
278 import matplotlib.pyplot as plt
279
280 pixel_x = 1200
281 pixel_y = 1080
282 radius = 30
283 center = (0.0, 0.0, 0.0)
284
285 aspect = pixel_x / pixel_y
286 pic_range = [-radius * aspect, radius * aspect, -radius, radius]
287 delta_x = (radius * 2 * aspect, 0.0, 0.0)
288 delta_y = (0.0, 0.0, radius * 2)
289
290 arr_rho = model.render_cartesian_column_integ(
291 "rho", "f64", center=(0.0, 0.0, 0.0), delta_x=delta_x, delta_y=delta_y, nx=pixel_x, ny=pixel_y
292 )
293
294 import copy
295
296 import matplotlib
297
298 my_cmap = copy.copy(matplotlib.colormaps.get_cmap("gist_heat")) # copy the default cmap
299 my_cmap.set_bad(color="black")
300
301 fig_width = 6
302 fig_height = fig_width / aspect
303 plt.figure(figsize=(fig_width, fig_height))
304 res = plt.imshow(arr_rho, cmap=my_cmap, origin="lower", extent=pic_range, norm="log", vmin=1e-9)
305
306 cbar = plt.colorbar(res, extend="both")
307 cbar.set_label(r"$\int \rho \, \mathrm{d} z$ [code unit]")
308 # or r"$\rho$ [code unit]" for slices
309
310 plt.title("t = {:0.3f} [code unit]".format(model.get_time()))
311 plt.xlabel("x")
312 plt.ylabel("z")
313 plt.show()
![t = 0.009 [code unit]](../../_images/sphx_glr_run_sph_custom_warp_profile_001.png)
Info: compute_column_integ field_name: rho, center: (0,0,0), delta_x: (66.66666666666667,0,0), delta_y: (0,0,60), nx: 1200, ny: 1080 [sph::CartesianRender][rank=0]
Info: compute_column_integ took 28.18 s [sph::CartesianRender][rank=0]
Total running time of the script: (2 minutes 8.333 seconds)
Estimated memory usage: 1091 MB