Note
Go to the end to download the full example code.
Testing Orszag-Tang vortex (3D thin box) with SPMHD#
CI test for the 3D (thin box) Orszag-Tang vortex, following the setup described in the Phantom code paper (Price et al. 2018, Section 5.7.3): a uniform density, periodic box with a sinusoidal velocity and magnetic field perturbation that develops into a system of interacting MHD shocks.
12 import os
13
14 import matplotlib.pyplot as plt
15 import numpy as np
16
17 import shamrock
18
19 # If we use the shamrock executable to run this script instead of the python interpreter,
20 # we should not initialize the system as the shamrock executable needs to handle specific MPI logic
21 if not shamrock.sys.is_initialized():
22 shamrock.change_loglevel(1)
23 shamrock.sys.init("0:0")
Parameters (Price et al. 2018, Sec 5.7.3)
28 kernel = "M4"
29
30 C_cour = 0.3
31 C_force = 0.25
32
33 nx = 32
34 ny = 32
35 nz = 3
36
37 xymin = -0.5
38 xmin = xymin
39 ymin = xymin
40
41 gamma = 5.0 / 3.0
42 betazero = 10.0 / 3.0
43 machzero = 1.0
44 vzero = 1.0
45 bzero = 1.0 / np.sqrt(4.0 * np.pi)
46
47 przero = 0.5 * bzero**2 * betazero
48 rhozero = gamma * przero * machzero
49 gam1 = gamma - 1.0
50 uuzero = przero / (gam1 * rhozero)
51
52 t_target = 0.5
53
54
55 def vel_func(r):
56 x, y, z = r
57 vx = -vzero * np.sin(2.0 * np.pi * (y - ymin))
58 vy = vzero * np.sin(2.0 * np.pi * (x - xmin))
59 return (vx, vy, 0.0)
60
61
62 def mag_func(r):
63 x, y, z = r
64 Bx = -bzero * np.sin(2.0 * np.pi * (y - ymin)) / rhozero
65 By = bzero * np.sin(4.0 * np.pi * (x - xmin)) / rhozero
66 return (Bx, By, 0.0)
67
68
69 def run_orszag_tang(nx, ny, nz=3):
70 ctx = shamrock.Context()
71 ctx.pdata_layout_new()
72
73 codeu = shamrock.UnitSystem(unit_time=1.0, unit_length=1.0, unit_mass=1.2566370621219e-06)
74 ucte = shamrock.Constants(codeu)
75 model = shamrock.get_Model_SPH(context=ctx, vector_type="f64_3", sph_kernel=kernel)
76
77 cfg = model.gen_default_config()
78 cfg.set_units(codeu)
79
80 cfg.set_artif_viscosity_None()
81 cfg.set_IdealMHD(sigma_mhd=1, sigma_u=1)
82 cfg.set_boundary_periodic()
83 cfg.set_eos_adiabatic(gamma)
84 model.set_solver_config(cfg)
85
86 crit_split = int(1e7)
87 crit_merge = 1
88 model.init_scheduler(crit_split, crit_merge)
89
90 (xs, ys, zs) = model.get_box_dim_fcc_3d(1, nx, ny, nz)
91 dr = 1 / xs
92 (xs, ys, zs) = model.get_box_dim_fcc_3d(dr, nx, ny, nz)
93
94 model.resize_simulation_box((-xs / 2, -ys / 2, -zs / 2), (xs / 2, ys / 2, zs / 2))
95 model.add_cube_fcc_3d(dr, (-xs / 2, -ys / 2, -zs / 2), (xs / 2, ys / 2, zs / 2))
96 model.set_value_in_a_box(
97 "uint", "f64", uuzero, (-xs / 2, -ys / 2, -zs / 2), (xs / 2, ys / 2, zs / 2)
98 )
99
100 model.set_field_value_lambda_f64_3("vxyz", vel_func)
101 model.set_field_value_lambda_f64_3("B/rho", mag_func)
102
103 vol_b = xs * ys * zs
104 totmass = rhozero * vol_b
105 pmass = model.total_mass_to_part_mass(totmass)
106 model.set_particle_mass(pmass)
107
108 print(f"[nx={nx}] Total mass :", totmass)
109 print(f"[nx={nx}] Current part mass :", pmass)
110
111 model.set_cfl_cour(C_cour)
112 model.set_cfl_force(C_force)
113
114 model.timestep()
115
116 return ctx, model, xs, ys, zs, pmass, ucte.mu_0()
117
118
119 ctx, model, xs, ys, zs, pmass, mu_0 = run_orszag_tang(nx, ny, nz)
120
121 model.evolve_until(t_target)
[nx=32] Total mass : 0.01435415559948049
[nx=32] Current part mass : 4.505384682825013e-06
---------------- t = 0, dt = 0 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 0.016923076923076926 unconverged cnt = 3186
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 0.01861538461538462 unconverged cnt = 3186
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 0.020476923076923084 unconverged cnt = 3186
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 0.022524615384615393 unconverged cnt = 3186
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 0.024777076923076933 unconverged cnt = 3186
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 0.02725478461538463 unconverged cnt = 3186
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 0.029980263076923093 unconverged cnt = 3178
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
Warning: smoothing length is not converged, rerunning the iterator ... [Smoothinglength][rank=0]
largest h = 0.0329782893846154 unconverged cnt = 58
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0, dt = 3.033731433136635e-05 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 3.033731433136635e-05, dt = 0.001031437962581657 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.0010617752769130234, dt = 0.0016976264484685334 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.0027594017253815568, dt = 0.0021356366249373863 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.004895038350318943, dt = 0.002418806853735241 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.007313845204054184, dt = 0.0026700013750796284 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.009983846579133813, dt = 0.002773642037475183 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.012757488616608996, dt = 0.0028010142809859677 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.015558502897594964, dt = 0.0028577746372138112 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.018416277534808775, dt = 0.002941129869256174 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.02135740740406495, dt = 0.0028564063306287346 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.024213813734693685, dt = 0.0028792313183021627 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.027093045052995847, dt = 0.002902213759073423 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.02999525881206927, dt = 0.0029291361613329285 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.0329243949734022, dt = 0.0029625739452987837 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.03588696891870099, dt = 0.0030030588818953344 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.038890027800596326, dt = 0.003052236188930685 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.04194226398952701, dt = 0.002742368941219658 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.04468463293074667, dt = 0.0027901580401330723 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.04747479097087974, dt = 0.0028326095093300853 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.050307400480209825, dt = 0.0028554783210417113 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.053162878801251534, dt = 0.002835728348659805 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.05599860714991134, dt = 0.0027927986915133934 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.05879140584142473, dt = 0.0028486935768232153 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.06164009941824795, dt = 0.0025706011675115104 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.06421070058575945, dt = 0.0026281393955539935 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.06683883998131344, dt = 0.002686366459414145 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.06952520644072759, dt = 0.002750488301862423 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.07227569474259, dt = 0.00282599512954057 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.07510168987213058, dt = 0.002896294785169768 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.07799798465730035, dt = 0.0028991343909787227 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.08089711904827906, dt = 0.002938456287501252 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.08383557533578032, dt = 0.0029913168773705767 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.08682689221315089, dt = 0.0029571832280422625 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.08978407544119316, dt = 0.0030218771869806042 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.09280595262817376, dt = 0.003085624431122113 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.09589157705929588, dt = 0.0031441906204539853 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.09903576767974986, dt = 0.0032009987795236245 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.10223676645927349, dt = 0.003282075102662769 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.10551884156193625, dt = 0.0033784125545342812 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.10889725411647053, dt = 0.003323999560157249 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.11222125367662777, dt = 0.003296141542512871 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.11551739521914065, dt = 0.0033018687533237726 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.11881926397246442, dt = 0.00328991782363642 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.12210918179610085, dt = 0.003242292423378229 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.12535147421947906, dt = 0.0032121565900727457 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.1285636308095518, dt = 0.0032553347553047 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.1318189655648565, dt = 0.0031677212387990193 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.13498668680365553, dt = 0.003153877934640941 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.13814056473829647, dt = 0.0030214671694061934 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.14116203190770266, dt = 0.002971659631602173 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.14413369153930483, dt = 0.0029489830329003383 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.14708267457220517, dt = 0.0029606050896834606 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.15004327966188863, dt = 0.0028421227584198338 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.15288540242030846, dt = 0.00284432153159224 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.1557297239519007, dt = 0.002831902509442114 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.15856162646134284, dt = 0.002832357112130397 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.16139398357347323, dt = 0.0027830177962468224 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.16417700136972005, dt = 0.0027843453347914057 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.16696134670451146, dt = 0.0026955123916991095 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.16965685909621056, dt = 0.0026406642458044814 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.17229752334201504, dt = 0.002612749455775484 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.1749102727977905, dt = 0.00264450642498602 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.17755477922277652, dt = 0.0026593957242491114 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.18021417494702563, dt = 0.0026873882192648734 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.1829015631662905, dt = 0.0026082572489186514 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.18550982041520914, dt = 0.0025773743193899133 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.18808719473459906, dt = 0.0026127883656065394 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.1906999831002056, dt = 0.002658191343635494 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.1933581744438411, dt = 0.002671156310461559 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.19602933075430265, dt = 0.0025820322522970124 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.19861136300659965, dt = 0.0026349938126811467 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.2012463568192808, dt = 0.002682325791323772 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.20392868261060457, dt = 0.0025133197113263142 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.20644200232193088, dt = 0.002541578607633636 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.2089835809295645, dt = 0.0025869081186845654 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.21157048904824907, dt = 0.0026188534590828854 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.21418934250733196, dt = 0.002660852595712909 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.21685019510304487, dt = 0.0026038109250228997 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.21945400602806778, dt = 0.0026432774343134913 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.22209728346238128, dt = 0.002510128157296585 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.22460741161967787, dt = 0.0023911730958748826 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.22699858471555276, dt = 0.0024156216438691786 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.22941420635942195, dt = 0.0023927813338326703 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.2318069876932546, dt = 0.0023120878804907758 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.2341190755737454, dt = 0.0023313415047011158 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.23645041707844652, dt = 0.002365539822237645 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.23881595690068416, dt = 0.002336595718729727 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.2411525526194139, dt = 0.002416285436429048 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.24356883805584295, dt = 0.002358883032914871 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.24592772108875782, dt = 0.0023733814629670357 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.24830110255172486, dt = 0.002404665516944625 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.25070576806866945, dt = 0.0024555495802389194 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.25316131764890837, dt = 0.0024377006691237947 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.25559901831803217, dt = 0.002538332802226522 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.2581373511202587, dt = 0.0025776357362716222 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.26071498685653033, dt = 0.0025931183039620102 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.26330810516049236, dt = 0.0025805021488309455 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.2658886073093233, dt = 0.0024896839895845145 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.2683782912989078, dt = 0.0025233634383153404 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.2709016547372231, dt = 0.0025797502368433253 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.2734814049740664, dt = 0.002623379466377052 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.27610478444044345, dt = 0.002688311088951301 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.27879309552939474, dt = 0.002704022753635474 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.2814971182830302, dt = 0.00271036982186821 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.28420748810489843, dt = 0.002660715756450602 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.28686820386134904, dt = 0.002707707923686889 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.28957591178503594, dt = 0.0026928636680525375 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.2922687754530885, dt = 0.0027384484119079735 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.29500722386499645, dt = 0.0027880510531527337 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.29779527491814917, dt = 0.002844778442487287 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.30064005336063643, dt = 0.0027946993649719605 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3034347527256084, dt = 0.0028852857509845974 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.306320038476593, dt = 0.002934193447490817 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3092542319240838, dt = 0.002805127096990373 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3120593590210742, dt = 0.0028885768692379532 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3149479358903121, dt = 0.0028438220935711493 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.31779175798388326, dt = 0.0027826996416755864 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.32057445762555886, dt = 0.002862843765714108 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.32343730139127297, dt = 0.0027971870523728825 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.32623448844364583, dt = 0.0028555860961777635 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3290900745398236, dt = 0.0027891373559386892 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.33187921189576225, dt = 0.0028374530308694157 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.33471666492663166, dt = 0.002910277362997476 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.33762694228962914, dt = 0.0029080522420572735 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3405349945316864, dt = 0.0026503662647422285 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.34318536079642864, dt = 0.0026251733035954406 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3458105341000241, dt = 0.002667145819351389 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3484776799193755, dt = 0.0027274570658869226 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3512051369852624, dt = 0.0028098955514387004 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3540150325367011, dt = 0.0029189913858782114 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3569340239225793, dt = 0.002887284302166794 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.35982130822474606, dt = 0.002925846623724978 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.36274715484847103, dt = 0.002924587631044734 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3656717424795158, dt = 0.0029985772979750592 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3686703197774908, dt = 0.0030065698663788816 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3716768896438697, dt = 0.002889750802606165 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.37456664044647586, dt = 0.0029599366644587383 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3775265771109346, dt = 0.002985873189398533 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.38051245030033315, dt = 0.0029973351969724282 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3835097854973056, dt = 0.002960214870350956 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.38647000036765655, dt = 0.0029819795510055264 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3894519799186621, dt = 0.003037154163210768 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.39248913408187286, dt = 0.003101989741769651 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.39559112382364253, dt = 0.0030902577660956963 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.3986813815897382, dt = 0.0029923244121638003 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.401673706001902, dt = 0.0030563478905641104 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.40473005389246614, dt = 0.0031428891716875646 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.4078729430641537, dt = 0.0031927419052565718 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.4110656849694103, dt = 0.0030678299848718475 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.4141335149542822, dt = 0.003098914585759383 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.4172324295400416, dt = 0.0031525424609631417 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.4203849720010047, dt = 0.003210725851575855 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.42359569785258055, dt = 0.0032892657971469 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.42688496364972744, dt = 0.003357149761513026 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.43024211341124047, dt = 0.003378844428219875 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.43362095783946036, dt = 0.0032369895947596773 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.43685794743422, dt = 0.003348427072921212 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.44020637450714123, dt = 0.003369462124969044 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.44357583663211025, dt = 0.003477608970171239 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.4470534456022815, dt = 0.003552916803054394 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.4506063624053359, dt = 0.0035442874703774836 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.4541506498757134, dt = 0.003544946330021918 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.4576955962057353, dt = 0.0035783505304329295 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.4612739467361682, dt = 0.003617092430427224 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.46489103916659547, dt = 0.00370003851308489 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.46859107767968033, dt = 0.003774781852861722 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.47236585953254207, dt = 0.0037186488750704004 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.47608450840761246, dt = 0.003616173015480347 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.4797006814230928, dt = 0.003684027487475568 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.4833847089105684, dt = 0.0037591441321218663 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.48714385304269026, dt = 0.0036814092107613334 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.4908252622534516, dt = 0.0037493908910763364 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.4945746531445279, dt = 0.0038492532129886367 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
---------------- t = 0.49842390635751654, dt = 0.001576093642483456 ----------------
Warning: High interface/patch volume ratio. [InterfaceGen][rank=0]
This can lead to high mpi overhead, try to increase the patch split crit
EvolveUntilResults(reach_target_time=true, reach_niter_max=false, reach_max_walltime=false, iter_count=175)
Diagnostics
126 kinetic_energy = shamrock.model_sph.analysisEnergyKinetic(model=model).get_kinetic_energy()
127 total_momentum = shamrock.model_sph.analysisTotalMomentum(model=model).get_total_momentum()
128
129 data = ctx.collect_data()
130 hpart = data["hpart"]
131 hfact = model.get_hfact()
132 rho = pmass * (hfact / hpart) ** 3
133 B_on_rho = data["B/rho"]
134
135 magnetic_energy = 0.5 * pmass * np.sum(rho * np.sum(B_on_rho**2, axis=1)) / mu_0
136
137 rho_max = np.max(rho)
138 rho_min = np.min(rho)
139
140 momentum_norm = np.sqrt(total_momentum[0] ** 2 + total_momentum[1] ** 2 + total_momentum[2] ** 2)
141
142 print(f"kinetic_energy = {kinetic_energy}")
143 print(f"magnetic_energy = {magnetic_energy}")
144 print(f"momentum_norm = {momentum_norm}")
145 print(f"rho_max = {rho_max}")
146 print(f"rho_min = {rho_min}")
147
148
149 plot_extra_resols = [] # 128, 148
150
151 do_plot = False
152 dump_folder = "_to_trash"
153 if do_plot and shamrock.sys.world_rank() == 0:
154 os.makedirs(dump_folder, exist_ok=True)
155
156
157 def render_density_z0(model, xs, ys, render_res=300, min_normalization=1e-9):
158 kwargs = dict(
159 center=(0.0, 0.0, 0.0),
160 delta_x=(xs, 0.0, 0.0),
161 delta_y=(0.0, ys, 0.0),
162 nx=render_res,
163 ny=render_res,
164 )
165 raw = np.asarray(model.render_cartesian_slice("rho", "f64", **kwargs))
166 unity = np.asarray(model.render_cartesian_slice("unity", "f64", **kwargs))
167 return np.where(unity < min_normalization, np.nan, raw / unity)
168
169
170 def pressure_cut(model, y0, xs, n_pts=400, min_normalization=1e-9):
171 x_arr = np.linspace(-xs / 2, xs / 2, n_pts)
172 positions = [(x, y0, 0.0) for x in x_arr]
173
174 rho_field = model.compute_field("rho", "f64")
175 uint_field = model.compute_field("uint", "f64")
176
177 def compute_pressure(size, rho, uint):
178 return (gamma - 1.0) * rho * uint
179
180 P_field = shamrock.map_fields_f64(compute_pressure, rho=rho_field, uint=uint_field)
181
182 P_raw = np.asarray(model.render_slice(P_field, positions))
183 unity = np.asarray(model.render_slice("unity", "f64", positions))
184 P_cut = np.where(unity < min_normalization, np.nan, P_raw / unity)
185 return x_arr, P_cut
186
187
188 if do_plot:
189 y_cuts = [0.3125, 0.4277]
190
191 # snapshot at t=0.5 (already reached above) for the primary resolution
192 resols = [(nx, ny)]
193 density_t05 = [render_density_z0(model, xs, ys)]
194 cuts_t05 = [[pressure_cut(model, y0, xs) for y0 in y_cuts]]
195
196 # continue to t=1 for the Fig. 32 bottom row
197 model.evolve_until(1.0)
198 density_t10 = [render_density_z0(model, xs, ys)]
199
200 for nx_e, ny_e in plot_extra_resols:
201 ctx_e, model_e, xs_e, ys_e, zs_e, pmass_e, mu_0_e = run_orszag_tang(nx_e, ny_e, nz)
202
203 model_e.evolve_until(0.5)
204 resols.append((nx_e, ny_e))
205 density_t05.append(render_density_z0(model_e, xs_e, ys_e))
206 cuts_t05.append([pressure_cut(model_e, y0, xs_e) for y0 in y_cuts])
207
208 model_e.evolve_until(1.0)
209 density_t10.append(render_density_z0(model_e, xs_e, ys_e))
210
211 if shamrock.sys.world_rank() == 0:
212 # Figure 32: density z=0, rows = t=0.5 / t=1, columns = resolutions
213 n_res = len(resols)
214 fig32, axs32 = plt.subplots(2, n_res, figsize=(max(4 * n_res, 6), 8), squeeze=False)
215 for icol, (nxr, nyr) in enumerate(resols):
216 extent = [-xs / 2, xs / 2, -ys / 2, ys / 2]
217 axs32[0, icol].imshow(
218 density_t05[icol], origin="lower", extent=extent, cmap="gist_heat"
219 )
220 axs32[0, icol].set_title(f"nx={nxr}, t=0.5")
221 axs32[1, icol].imshow(
222 density_t10[icol], origin="lower", extent=extent, cmap="gist_heat"
223 )
224 axs32[1, icol].set_title(f"nx={nxr}, t=1")
225 for ax in (axs32[0, icol], axs32[1, icol]):
226 ax.set_xlabel("x")
227 ax.set_ylabel("y")
228 fig32.suptitle("Orszag-Tang vortex: density in a z=0 cross section")
229 fig32.tight_layout()
230 fig32.savefig(os.path.join(dump_folder, "mhd_orszag_tang_density_z0.png"), dpi=150)
231 plt.close(fig32)
232
233 # Figure 33: horizontal pressure cuts at t=0.5, one panel per y0,
234 # one line per resolution
235 fig33, axs33 = plt.subplots(len(y_cuts), 1, figsize=(7, 8), squeeze=False)
236 axs33 = axs33[:, 0]
237 for irow, y0 in enumerate(y_cuts):
238 for (nxr, nyr), cuts in zip(resols, cuts_t05):
239 x_arr, P_cut = cuts[irow]
240 axs33[irow].plot(x_arr, P_cut, label=f"nx={nxr}")
241 axs33[irow].set_title(f"y={y0}")
242 axs33[irow].set_xlabel("x")
243 axs33[irow].set_ylabel("P")
244 axs33[irow].legend()
245 fig33.suptitle("Orszag-Tang vortex: horizontal pressure cuts (z=0, t=0.5)")
246 fig33.tight_layout()
247 fig33.savefig(os.path.join(dump_folder, "mhd_orszag_tang_pressure_cuts.png"), dpi=150)
248 plt.close(fig33)
249
250 test_pass = True
251 err_log = ""
252
253
254 expect_kinetic_energy = 0.0021986783020805675
255 expect_magnetic_energy = 0.0029601504573239087
256 expect_rho_max = 0.3400457561374047
257 expect_rho_min = 0.14169354827263864
258
259 tol = 1e-4 # relative tolerance
260 max_momentum_norm = 1e-3
261
262
263 def float_equal(val1, val2, prec):
264 return abs(val1 - val2) < prec
265
266
267 if not float_equal(kinetic_energy, expect_kinetic_energy, tol * expect_kinetic_energy):
268 err_log += "error on the kinetic energy is outside of tolerances:\n"
269 err_log += f" expected error = {expect_kinetic_energy} +- {tol * expect_kinetic_energy}\n"
270 err_log += f" obtained error = {kinetic_energy}\n"
271 test_pass = False
272
273 if not float_equal(magnetic_energy, expect_magnetic_energy, tol * expect_magnetic_energy):
274 err_log += "error on the magnetic energy is outside of tolerances:\n"
275 err_log += f" expected error = {expect_magnetic_energy} +- {tol * expect_magnetic_energy}\n"
276 err_log += f" obtained error = {magnetic_energy}\n"
277 test_pass = False
278
279 if momentum_norm > max_momentum_norm:
280 err_log += "total momentum norm is outside of tolerances:\n"
281 err_log += f" expected |p| < {max_momentum_norm}\n"
282 err_log += f" obtained |p| = {momentum_norm}\n"
283 test_pass = False
284
285 if not float_equal(rho_max, expect_rho_max, tol * expect_rho_max):
286 err_log += "error on rho_max is outside of tolerances:\n"
287 err_log += f" expected error = {expect_rho_max} +- {tol * expect_rho_max}\n"
288 err_log += f" obtained error = {rho_max}\n"
289 test_pass = False
290
291 if not float_equal(rho_min, expect_rho_min, tol * expect_rho_min):
292 err_log += "error on rho_min is outside of tolerances:\n"
293 err_log += f" expected error = {expect_rho_min} +- {tol * expect_rho_min}\n"
294 err_log += f" obtained error = {rho_min}\n"
295 test_pass = False
296
297 if test_pass == False:
298 exit("Test did not pass L2 margins : \n" + err_log)
kinetic_energy = 0.0021986783020805675
magnetic_energy = 0.0029601504573239087
momentum_norm = 7.683768519158053e-05
rho_max = 0.3400457561374047
rho_min = 0.14169354827263864
Total running time of the script: (0 minutes 10.636 seconds)
Estimated memory usage: 162 MB