.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "_as_gen/tests_ci/run_mhd_rotor.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr__as_gen_tests_ci_run_mhd_rotor.py: Testing the MHD rotor problem (3D) with SPMHD =============================================== CI test for the 3D MHD rotor problem (Balsara & Spicer 1999; Toth 2000's 'first rotor problem'), following the setup described in the Phantom code paper (Price et al. 2018, Section 5.6.4): a rapidly rotating dense disc of material immersed in a uniform ambient medium threaded by a uniform magnetic field, testing the propagation of rotational discontinuities. .. GENERATED FROM PYTHON SOURCE LINES 12-26 .. code-block:: Python :lineno-start: 13 import os import matplotlib.pyplot as plt import numpy as np import shamrock # If we use the shamrock executable to run this script instead of the python interpreter, # we should not initialize the system as the shamrock executable needs to handle specific MPI logic if not shamrock.sys.is_initialized(): shamrock.change_loglevel(1) shamrock.sys.init("0:0") .. GENERATED FROM PYTHON SOURCE LINES 27-51 .. code-block:: Python :lineno-start: 27 kernel = "M4" C_cour = 0.3 C_force = 0.25 nx = 32 nz_thin = 12 gamma = 1.4 R_disc = 0.1 rho_disc = 10.0 rho_ambient = 1.0 P0 = 1.0 v0 = 2.0 Bx0 = 5.0 / np.sqrt(4.0 * np.pi) fact = (rho_disc / rho_ambient) ** (1.0 / 3.0) u_disc = P0 / ((gamma - 1) * rho_disc) u_ambient = P0 / ((gamma - 1) * rho_ambient) t_target = 0.15 .. GENERATED FROM PYTHON SOURCE LINES 52-77 .. code-block:: Python :lineno-start: 53 ctx = shamrock.Context() ctx.pdata_layout_new() si = shamrock.UnitSystem() codeu = shamrock.UnitSystem(unit_time=1.0, unit_length=1.0, unit_mass=1.2566370621219e-06) ucte = shamrock.Constants(codeu) model = shamrock.get_Model_SPH(context=ctx, vector_type="f64_3", sph_kernel=kernel) cfg = model.gen_default_config() cfg.set_units(codeu) mu_0 = ucte.mu_0() cfg.set_artif_viscosity_None() cfg.set_IdealMHD(sigma_mhd=1, sigma_u=1) cfg.set_boundary_periodic() cfg.set_eos_adiabatic(gamma) model.set_solver_config(cfg) cfg.print_status() crit_split = int(1e7) crit_merge = 1 model.init_scheduler(crit_split, crit_merge) .. rst-class:: sphx-glr-script-out .. code-block:: none ----- SPH Solver configuration ----- [ { "artif_viscosity": { "type": "none" }, "boundary_config": { "bc_type": "periodic" }, "cfl_config": { "cfl_cour": 0.3, "cfl_force": 0.25, "cfl_multiplier_stiffness": 2.0, "eta_sink": 0.05 }, "combined_dtdiv_divcurlv_compute": false, "debug_dump_filename": "", "do_debug_dump": false, "dust_config": { "ballabio_ts_limiter": false, "drag_mode": { "type": "none" }, "evol_mode": { "type": "none" }, "mode": { "type": "none" } }, "enable_particle_reordering": false, "eos_config": { "Tvec": "f64_3", "eos_type": "adiabatic", "gamma": 1.4 }, "epsilon_h": 1e-06, "ext_force_config": { "force_list": [] }, "gpart_mass": 0.0, "h_iter_per_subcycles": 50, "h_max_subcycles_count": 100, "htol_up_coarse_cycle": 1.1, "htol_up_fine_cycle": 1.1, "kernel_id": "M4", "mhd_config": { "alpha_AV": 1.0, "alpha_B": 1.0, "alpha_u": 1.0, "beta_AV": 1.0, "mhd_type": "ideal_mhd_constrained_hyper_para", "sigma_mhd": 1.0 }, "neigh_cache_strategy": "two_stage", "particle_killing": [], "particle_reordering_step_freq": 1000, "save_dt_to_fields": false, "scheduler_config": { "merge_load_value": 0, "split_load_value": 0 }, "self_grav_config": { "softening_length": 1e-09, "softening_mode": "plummer", "type": "none" }, "show_cfl_detail": false, "show_ghost_zone_graph": false, "show_neigh_stats": false, "smoothing_length_config": { "type": "density_based" }, "tree_reduction_level": 3, "type_id": "sycl::vec", "unit_sys": { "unit_current": 1.0, "unit_length": 1.0, "unit_lumint": 1.0, "unit_mass": 1.2566370621219e-06, "unit_qte": 1.0, "unit_temperature": 1.0, "unit_time": 1.0 } } ] ------------------------------------ .. GENERATED FROM PYTHON SOURCE LINES 78-153 .. code-block:: Python :lineno-start: 79 (xs, ys, zs) = model.get_box_dim_fcc_3d(1, nx, nx, nz_thin) dr = 1.0 / xs (xs, ys, zs) = model.get_box_dim_fcc_3d(dr, nx, nx, nz_thin) box_min = (-xs / 2, -ys / 2, -zs / 2) box_max = (xs / 2, ys / 2, zs / 2) model.resize_simulation_box(box_min, box_max) dr_disc = dr / fact disc_box_min = (-R_disc, -R_disc, -zs / 2) disc_box_max = (R_disc, R_disc, zs / 2) def in_disc(r): x, y, z = r return (x * x + y * y) < R_disc * R_disc def out_disc(r): x, y, z = r return (x * x + y * y) >= R_disc * R_disc setup = model.get_setup() gen_ambient = setup.make_generator_lattice_hcp(dr, box_min, box_max) gen_ambient_masked = setup.make_modifier_filter(parent=gen_ambient, filter=out_disc) gen_disc = setup.make_generator_lattice_hcp(dr_disc, disc_box_min, disc_box_max) gen_disc_masked = setup.make_modifier_filter(parent=gen_disc, filter=in_disc) comb = setup.make_combiner_add(gen_ambient_masked, gen_disc_masked) setup.apply_setup(comb) # Internal energy: ambient value everywhere, then overwrite inside the disc model.set_value_in_a_box("uint", "f64", u_ambient, box_min, box_max) model.set_value_in_sphere("uint", "f64", u_disc, (0, 0, 0), R_disc) def vel_func(r): x, y, z = r r2 = x * x + y * y if r2 >= R_disc * R_disc: return (0.0, 0.0, 0.0) rcyl = np.sqrt(max(r2, 1e-30)) vx = -v0 * y / rcyl vy = v0 * x / rcyl return (vx, vy, 0.0) def mag_func(r): x, y, z = r rho_local = rho_disc if (x * x + y * y) < R_disc * R_disc else rho_ambient return (Bx0 / rho_local, 0.0, 0.0) model.set_field_value_lambda_f64_3("vxyz", vel_func) model.set_field_value_lambda_f64_3("B/rho", mag_func) vol_b = xs * ys * zs vol_disc = np.pi * R_disc * R_disc * zs totmass = rho_ambient * (vol_b - vol_disc) + rho_disc * vol_disc pmass = model.total_mass_to_part_mass(totmass) model.set_particle_mass(pmass) print("Total mass :", totmass) print("Current part mass :", pmass) print("Total particle count :", model.get_total_part_count()) model.set_cfl_cour(C_cour) model.set_cfl_force(C_force) model.timestep() model.evolve_until(t_target) .. rst-class:: sphx-glr-script-out .. code-block:: none SPH setup: generating particles ... SPH setup: Nstep = 16057 ( 1.6e+04 ) Ntotal = 16057 ( 1.6e+04 rank min = 7.9e+05 max = 1.6e+04) rate = 1.605700e+04 N.s^-1 SPH setup: the generation step took : 0.023222084 s SPH setup: final particle count = 16057 beginning injection ... SPH setup: injected 16057 / 16057 => 100.0% | ranks with patchs = 1 / 1 <- global loop -> (msg count : 0) SPH setup: the injection step took : 0.004735663 s SPH setup: the setup took : 0.032173955000000004 s Total mass : 0.3517371564808811 Current part mass : 2.1905533815836154e-05 Total particle count : 16057 ---------------- 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.030624643436074166 unconverged cnt = 15757 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.033441290467993 unconverged cnt = 1412 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.03344129046799447 unconverged cnt = 612 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.03344129046799447 unconverged cnt = 24 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 = 9.361305039293543e-06 ---------------- Warning: High interface/patch volume ratio. [InterfaceGen][rank=0] This can lead to high mpi overhead, try to increase the patch split crit ---------------- t = 9.361305039293543e-06, dt = 0.0003179667331182061 ---------------- 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.00032732803815749967, dt = 0.0005068913408230215 ---------------- Warning: High interface/patch volume ratio. [InterfaceGen][rank=0] This can lead to high mpi overhead, try to increase the patch split crit Warning: the corrector tolerance are broken the step will be re rerunned [BasicGasSPH][rank=0] eps_v = 0.02140365637482218 ---------------- t = 0.0008342193789805212, dt = 0.0003088471520994329 ---------------- 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.0011430665310799542, dt = 0.0004811609334120151 ---------------- Warning: High interface/patch volume ratio. [InterfaceGen][rank=0] This can lead to high mpi overhead, try to increase the patch split crit Warning: the corrector tolerance are broken the step will be re rerunned [BasicGasSPH][rank=0] eps_v = 0.01650282944842808 ---------------- t = 0.0016242274644919694, dt = 0.00028852941834372825 ---------------- 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.0019127568828356976, dt = 0.0004513250590017535 ---------------- Warning: High interface/patch volume ratio. [InterfaceGen][rank=0] This can lead to high mpi overhead, try to increase the patch split crit Warning: the corrector tolerance are broken the step will be re rerunned [BasicGasSPH][rank=0] eps_v = 0.011857336605177899 ---------------- t = 0.002364081941837451, dt = 0.0002855860254006608 ---------------- 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.0026496679672381117, dt = 0.0004472462443549473 ---------------- 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.003096914211593059, dt = 0.0005430295859126931 ---------------- Warning: High interface/patch volume ratio. [InterfaceGen][rank=0] This can lead to high mpi overhead, try to increase the patch split crit Warning: the corrector tolerance are broken the step will be re rerunned [BasicGasSPH][rank=0] eps_v = 0.010949868689275964 ---------------- t = 0.003639943797505752, dt = 0.000299592007521199 ---------------- 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.0039395358050269514, dt = 0.00044019610690170175 ---------------- 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.004379731911928653, dt = 0.0005272079132851575 ---------------- 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.00490693982521381, dt = 0.0005800282773933735 ---------------- 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.005486968102607184, dt = 0.0006117457754147015 ---------------- 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.006098713878021886, dt = 0.0006306730880885195 ---------------- 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.006729386966110406, dt = 0.0006420545873200215 ---------------- 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.007371441553430427, dt = 0.0006674765828164316 ---------------- 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.008038918136246858, dt = 0.0006881561454382069 ---------------- 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.008727074281685065, dt = 0.0006876285612350837 ---------------- 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.009414702842920149, dt = 0.0006876819242558283 ---------------- 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.010102384767175978, dt = 0.0006884869190354837 ---------------- 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.010790871686211461, dt = 0.0006901904934907175 ---------------- 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.011481062179702179, dt = 0.0007246555609302711 ---------------- 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.012205717740632449, dt = 0.0007380571802197036 ---------------- 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.012943774920852153, dt = 0.0007432961584233591 ---------------- 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.013687071079275512, dt = 0.0007946193052820766 ---------------- 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.014481690384557589, dt = 0.0008043857606333408 ---------------- 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.01528607614519093, dt = 0.000810557241282707 ---------------- 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.016096633386473637, dt = 0.000869798410645941 ---------------- 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.01696643179711958, dt = 0.0008880191723396412 ---------------- 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.017854450969459222, dt = 0.0009223967573505114 ---------------- 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.018776847726809733, dt = 0.0009555110223614668 ---------------- 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.0197323587491712, dt = 0.0009716676481982024 ---------------- 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.0207040263973694, dt = 0.0010310594646010163 ---------------- 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.021735085861970416, dt = 0.0010887068756253 ---------------- 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.022823792737595716, dt = 0.0011425071968699525 ---------------- 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.023966299934465668, dt = 0.0011601819691180492 ---------------- 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.025126481903583717, dt = 0.0011760313594803226 ---------------- 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.02630251326306404, dt = 0.0012554323209557296 ---------------- 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.027557945584019767, dt = 0.0012969452529248983 ---------------- 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.028854890836944664, dt = 0.0013214152009816567 ---------------- 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.030176306037926322, dt = 0.0013496314261051182 ---------------- 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.03152593746403144, dt = 0.0013440212592291147 ---------------- 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.03286995872326055, dt = 0.0013384536616935679 ---------------- 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.03420841238495412, dt = 0.0013349538797200041 ---------------- 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.03554336626467412, dt = 0.001333383600951829 ---------------- 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.03687674986562595, dt = 0.001340660334262886 ---------------- 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.03821741019988884, dt = 0.00134301079756071 ---------------- 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.03956042099744955, dt = 0.0013469611038881644 ---------------- 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.040907382101337715, dt = 0.0013773134135356031 ---------------- 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.04228469551487332, dt = 0.0013902880029702594 ---------------- 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.04367498351784358, dt = 0.0013976001253743586 ---------------- 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.04507258364321794, dt = 0.001439730267489728 ---------------- 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.04651231391070766, dt = 0.0014573101904707842 ---------------- 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.04796962410117844, dt = 0.0014769751924910638 ---------------- 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.049446599293669506, dt = 0.0014977457475042247 ---------------- 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.05094434504117373, dt = 0.0016130734755654937 ---------------- 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.05255741851673923, dt = 0.001560132240544997 ---------------- 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.05411755075728423, dt = 0.0015791089188168035 ---------------- 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.05569665967610103, dt = 0.001604555836432781 ---------------- 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.05730121551253381, dt = 0.001624904882025502 ---------------- 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.058926120394559314, dt = 0.0016661173622373644 ---------------- 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.06059223775679668, dt = 0.001664261983037761 ---------------- 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.062256499739834444, dt = 0.0017138723596756653 ---------------- 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.0639703720995101, dt = 0.001688242091809867 ---------------- 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.06565861419131998, dt = 0.0017368068814297465 ---------------- 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.06739542107274972, dt = 0.0017640110309582917 ---------------- 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.06915943210370801, dt = 0.0017810701446123861 ---------------- 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.0709405022483204, dt = 0.0017802871210402858 ---------------- 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.07272078936936067, dt = 0.0018171687504317034 ---------------- 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.07453795811979237, dt = 0.0018639400204733703 ---------------- 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.07640189814026574, dt = 0.0018354687195079366 ---------------- 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.07823736685977367, dt = 0.0018934198296525436 ---------------- 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.08013078668942622, dt = 0.001888444440242755 ---------------- 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.08201923112966897, dt = 0.001962635009866464 ---------------- 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.08398186613953543, dt = 0.002024510240646751 ---------------- 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.08600637638018217, dt = 0.0019941766544564043 ---------------- 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.08800055303463858, dt = 0.0019826948515308633 ---------------- 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.08998324788616945, dt = 0.0020357798477004117 ---------------- 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.09201902773386986, dt = 0.00199297026940463 ---------------- 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.09401199800327449, dt = 0.0020549440261365464 ---------------- 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.09606694202941103, dt = 0.002134463284204309 ---------------- 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.09820140531361535, dt = 0.002201741591010624 ---------------- 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.10040314690462597, dt = 0.0022040208203038 ---------------- 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.10260716772492977, dt = 0.002239544767372461 ---------------- 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.10484671249230224, dt = 0.0022528099849507535 ---------------- 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.10709952247725299, dt = 0.002293877361630998 ---------------- 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.10939339983888399, dt = 0.0023263480987760704 ---------------- 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.11171974793766007, dt = 0.0023577993869540403 ---------------- 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.11407754732461411, dt = 0.002332722196151508 ---------------- 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.11641026952076562, dt = 0.0023578627919003925 ---------------- Warning: High interface/patch volume ratio. [InterfaceGen][rank=0] This can lead to high mpi overhead, try to increase the patch split crit Warning: the corrector tolerance are broken the step will be re rerunned [BasicGasSPH][rank=0] eps_v = 0.010095640916472443 ---------------- t = 0.11876813231266602, dt = 0.0012047836661210328 ---------------- 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.11997291597878705, dt = 0.0016243952988796132 ---------------- 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.12159731127766667, dt = 0.0019273486369679076 ---------------- 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.12352465991463457, dt = 0.001960795727008812 ---------------- Warning: High interface/patch volume ratio. [InterfaceGen][rank=0] This can lead to high mpi overhead, try to increase the patch split crit Warning: the corrector tolerance are broken the step will be re rerunned [BasicGasSPH][rank=0] eps_v = 0.010385101868818584 ---------------- t = 0.12548545564164337, dt = 0.0010247839354388978 ---------------- 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.12651023957708227, dt = 0.0014365352514479684 ---------------- 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.12794677482853023, dt = 0.00170953395210459 ---------------- 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.12965630878063483, dt = 0.0019038935488683586 ---------------- Warning: High interface/patch volume ratio. [InterfaceGen][rank=0] This can lead to high mpi overhead, try to increase the patch split crit Warning: the corrector tolerance are broken the step will be re rerunned [BasicGasSPH][rank=0] eps_v = 0.011160061032476161 ---------------- t = 0.1315602023295032, dt = 0.0010273980210009173 ---------------- 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.1325876003505041, dt = 0.0014600945047063818 ---------------- 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.13404769485521048, dt = 0.001763836488213375 ---------------- 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.13581153134342386, dt = 0.001967828367139376 ---------------- Warning: High interface/patch volume ratio. [InterfaceGen][rank=0] This can lead to high mpi overhead, try to increase the patch split crit Warning: the corrector tolerance are broken the step will be re rerunned [BasicGasSPH][rank=0] eps_v = 0.010707174971430781 ---------------- t = 0.13777935971056324, dt = 0.0010659001646381998 ---------------- 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.13884525987520144, dt = 0.0015199181104031286 ---------------- 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.14036517798560458, dt = 0.001834190117654256 ---------------- 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.14219936810325884, dt = 0.002042259790331768 ---------------- Warning: High interface/patch volume ratio. [InterfaceGen][rank=0] This can lead to high mpi overhead, try to increase the patch split crit Warning: the corrector tolerance are broken the step will be re rerunned [BasicGasSPH][rank=0] eps_v = 0.010512362264616656 ---------------- t = 0.1442416278935906, dt = 0.0011120022534457327 ---------------- 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.14535363014703634, dt = 0.001593370625881224 ---------------- 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.14694700077291756, dt = 0.0019219380479688205 ---------------- 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.14886893882088637, dt = 0.0011310611791136294 ---------------- 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=112) .. GENERATED FROM PYTHON SOURCE LINES 154-155 Diagnostics .. GENERATED FROM PYTHON SOURCE LINES 155-308 .. code-block:: Python :lineno-start: 156 kinetic_energy = shamrock.model_sph.analysisEnergyKinetic(model=model).get_kinetic_energy() total_momentum = shamrock.model_sph.analysisTotalMomentum(model=model).get_total_momentum() data = ctx.collect_data() hpart = data["hpart"] hfact = model.get_hfact() rho = pmass * (hfact / hpart) ** 3 B_on_rho = data["B/rho"] magnetic_energy = 0.5 * pmass * np.sum(rho * np.sum(B_on_rho**2, axis=1)) / mu_0 rho_max = np.max(rho) rho_min = np.min(rho) momentum_norm = np.sqrt(total_momentum[0] ** 2 + total_momentum[1] ** 2 + total_momentum[2] ** 2) # nx=256 do_plot = False dump_folder = "_to_trash" if do_plot and shamrock.sys.world_rank() == 0: os.makedirs(dump_folder, exist_ok=True) if do_plot: render_res = 300 center = (0.0, 0.0, 0.0) delta_x = (xs, 0.0, 0.0) delta_y = (0.0, ys, 0.0) def render_slice_normalized(field, field_type=None, min_normalization=1e-9): kwargs = dict(center=center, delta_x=delta_x, delta_y=delta_y, nx=render_res, ny=render_res) if field_type is None: raw = np.asarray(model.render_cartesian_slice(field, **kwargs)) else: raw = np.asarray(model.render_cartesian_slice(field, field_type, **kwargs)) unity = np.asarray(model.render_cartesian_slice("unity", "f64", **kwargs)) unity_b = unity if raw.ndim == 2 else unity[..., None] out = raw / unity_b mask = np.broadcast_to(unity_b, raw.shape) < min_normalization return np.where(mask, np.nan, out) rho_field = model.compute_field("rho", "f64") uint_field = model.compute_field("uint", "f64") def compute_pressure(size, rho, uint): return (gamma - 1.0) * rho * uint P_field = shamrock.map_fields_f64(compute_pressure, rho=rho_field, uint=uint_field) def compute_soundspeed(size, rho, P): return np.sqrt(gamma * P / rho) cs_field = shamrock.map_fields_f64(compute_soundspeed, rho=rho_field, P=P_field) rho_slice = render_slice_normalized(rho_field) P_slice = render_slice_normalized(P_field) cs_slice = render_slice_normalized(cs_field) vxyz_slice = render_slice_normalized("vxyz", "f64_3") B_on_rho_slice = render_slice_normalized("B/rho", "f64_3") vmag_slice = np.linalg.norm(vxyz_slice, axis=-1) mach_slice = vmag_slice / cs_slice # physical B = (B/rho) * rho B_slice = B_on_rho_slice * rho_slice[..., None] magnetic_pressure_slice = 0.5 * np.sum(B_slice**2, axis=-1) / mu_0 if shamrock.sys.world_rank() == 0: # Contour limits identical to Toth (2000), as quoted in the # Phantom paper's Figure 34 caption. panels = [ ("Density", rho_slice, 0.483, 12.95), ("Pressure", P_slice, 0.0202, 2.008), ("Mach number $|v|/c_s$", mach_slice, 0.0, 1.09), ("Magnetic pressure $\\frac{1}{2}B^2$", magnetic_pressure_slice, 0.0, 2.642), ] ny_render, nx_render = rho_slice.shape xg = np.linspace(box_min[0], box_max[0], nx_render) yg = np.linspace(box_min[1], box_max[1], ny_render) X, Y = np.meshgrid(xg, yg) fig, axs = plt.subplots(2, 2, figsize=(10, 9)) for ax, (title, field_slice, vmin, vmax) in zip(axs.flat, panels): levels = np.linspace(vmin, vmax, 30) contour_set = ax.contour(X, Y, field_slice, levels=levels) ax.set_title(title) ax.set_xlabel("x") ax.set_ylabel("y") ax.set_aspect("equal") fig.colorbar(contour_set, ax=ax) fig.suptitle(f"MHD rotor problem, $t={t_target}$, nx={nx}") fig.tight_layout() fig.savefig(os.path.join(dump_folder, "mhd_rotor_contours.png"), dpi=150) plt.close(fig) print(f"kinetic_energy = {kinetic_energy}") print(f"magnetic_energy = {magnetic_energy}") print(f"momentum_norm = {momentum_norm}") print(f"rho_max = {rho_max}") print(f"rho_min = {rho_min}") test_pass = True err_log = "" expect_kinetic_energy = 0.05985165376928185 expect_magnetic_energy = 0.33746158787028185 expect_rho_max = 6.3060516178274515 expect_rho_min = 0.24993176850956453 tol = 1e-3 # relative tolerance max_momentum_norm = 1e-2 def float_equal(val1, val2, prec): return abs(val1 - val2) < prec if not float_equal(kinetic_energy, expect_kinetic_energy, tol * expect_kinetic_energy): err_log += "error on the kinetic energy is outside of tolerances:\n" err_log += f" expected error = {expect_kinetic_energy} +- {tol * expect_kinetic_energy}\n" err_log += f" obtained error = {kinetic_energy}\n" test_pass = False if not float_equal(magnetic_energy, expect_magnetic_energy, tol * expect_magnetic_energy): err_log += "error on the magnetic energy is outside of tolerances:\n" err_log += f" expected error = {expect_magnetic_energy} +- {tol * expect_magnetic_energy}\n" err_log += f" obtained error = {magnetic_energy}\n" test_pass = False if momentum_norm > max_momentum_norm: err_log += "total momentum norm is outside of tolerances:\n" err_log += f" expected |p| < {max_momentum_norm}\n" err_log += f" obtained |p| = {momentum_norm}\n" test_pass = False if not float_equal(rho_max, expect_rho_max, tol * expect_rho_max): err_log += "error on rho_max is outside of tolerances:\n" err_log += f" expected error = {expect_rho_max} +- {tol * expect_rho_max}\n" err_log += f" obtained error = {rho_max}\n" test_pass = False if not float_equal(rho_min, expect_rho_min, tol * expect_rho_min): err_log += "error on rho_min is outside of tolerances:\n" err_log += f" expected error = {expect_rho_min} +- {tol * expect_rho_min}\n" err_log += f" obtained error = {rho_min}\n" test_pass = False if test_pass == False: exit("Test did not pass L2 margins : \n" + err_log) .. rst-class:: sphx-glr-script-out .. code-block:: none kinetic_energy = 0.05985165376928191 magnetic_energy = 0.33746158787028185 momentum_norm = 0.00015219398459745176 rho_max = 6.3060516178274515 rho_min = 0.24993176850956453 .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 15.989 seconds) **Estimated memory usage:** 162 MB .. _sphx_glr_download__as_gen_tests_ci_run_mhd_rotor.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: run_mhd_rotor.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: run_mhd_rotor.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: run_mhd_rotor.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_