Small drift in MHD blast time integration #2
Description
I found a drift in the divergence values of the magnetic field when time_integrating a mhd blast with the following initial state and configuration with num_cells=512. The divergence goes over the threshold of e-8 set in the checkify check and by looking at the difference of the divergence values this happens around the last 10 timesteps approximately.
def get_blast_setup(num_cells):
# spatial domain
box_size = 1.0
# setup simulation config
config = SimulationConfig(
progress_bar=True,
mhd=True,
dimensionality=2,
box_size=box_size,
num_cells=num_cells,
differentiation_mode=BACKWARDS,
limiter=MINMOD,
riemann_solver=HLL,
exact_end_time=True,
runtime_debugging=True,
)
helper_data = get_helper_data(config)
params = SimulationParams(t_end=0.2, C_cfl=0.1)
registered_variables = get_registered_variables(config)
# Grid size and configuration
num_cells = config.num_cells
# --- Initial Conditions ---
grid_spacing = config.box_size / config.num_cells
x = jnp.linspace(
grid_spacing / 2, config.box_size - grid_spacing / 2, config.num_cells
)
y = jnp.linspace(
grid_spacing / 2, config.box_size - grid_spacing / 2, config.num_cells
)
X, Y = jnp.meshgrid(x, y, indexing="ij")
r = helper_data.r
# Initialize state
rho = jnp.ones_like(X)
P = jnp.ones_like(X) * 0.1
r_inj = 0.1 * box_size
p_inj = 10.0
P = jnp.where(r**2 < r_inj**2, p_inj, P)
V_x = jnp.zeros_like(X)
V_y = jnp.zeros_like(X)
B_0 = 1 / jnp.sqrt(2)
B_x = B_0 * jnp.ones_like(X)
B_y = jnp.zeros_like(X)
B_z = jnp.zeros_like(X)
initial_state = construct_primitive_state(
config=config,
registered_variables=registered_variables,
density=rho,
velocity_x=V_x,
velocity_y=V_y,
magnetic_field_x=B_x,
magnetic_field_y=B_y,
magnetic_field_z=B_z,
gas_pressure=P,
)
print("inital state built")
config = finalize_config(config, initial_state.shape)
return initial_state, config, params, helper_data, registered_variables