Sc/mhd treemesh 3d #1323
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,102 @@ | ||
| using OrdinaryDiffEq | ||
| using Trixi | ||
|
|
||
| ############################################################################### | ||
| # semidiscretization of the visco-resistive compressible MHD equations | ||
|
|
||
| prandtl_number() = 0.72 | ||
| mu_const = 2e-2 | ||
| eta_const = 2e-2 | ||
|
|
||
| equations = IdealGlmMhdEquations3D(5 / 3) | ||
| equations_parabolic = ViscoResistiveMhd3D(equations, mu = mu_const, | ||
| Prandtl = prandtl_number(), | ||
| eta = eta_const, | ||
| gradient_variables = GradientVariablesPrimitive()) | ||
|
|
||
| volume_flux = (flux_hindenlang_gassner, flux_nonconservative_powell) | ||
| solver = DGSEM(polydeg = 3, | ||
| surface_flux = (flux_lax_friedrichs, flux_nonconservative_powell), | ||
| volume_integral = VolumeIntegralFluxDifferencing(volume_flux)) | ||
|
|
||
| coordinates_min = (-1.0, -1.0, -1.0) # minimum coordinates (min(x), min(y), min(z)) | ||
| coordinates_max = (1.0, 1.0, 1.0) # maximum coordinates (max(x), max(y), max(z)) | ||
|
|
||
| # Create a uniformly refined mesh | ||
| mesh = TreeMesh(coordinates_min, coordinates_max, | ||
| initial_refinement_level = 2, | ||
| n_cells_max = 200_000) # set maximum capacity of tree data structure | ||
|
|
||
| function initial_condition_constant_alfven(x, t, equations) | ||
SimonCan marked this conversation as resolved.
Show resolved
Hide resolved
|
||
| # Alfvén wave in three space dimensions modified by a periodic density variation. | ||
| # For the system without the density variations see: Altmann thesis http://dx.doi.org/10.18419/opus-3895. | ||
| # Domain must be set to [-1, 1]^3, γ = 5/3. | ||
| omega = 2.0 * pi # may be multiplied by frequency | ||
| # r = length-variable = length of computational domain | ||
| r = 2.0 | ||
| # e = epsilon | ||
| e = 0.02 | ||
| nx = 1 / sqrt(r^2 + 1.0) | ||
| ny = r / sqrt(r^2 + 1.0) | ||
| sqr = 1.0 | ||
| Va = omega / (ny * sqr) | ||
| phi_alv = omega / ny * (nx * (x[1] - 0.5 * r) + ny * (x[2] - 0.5 * r)) - Va * t | ||
|
|
||
| rho = 1.0 + e*cos(phi_alv + 1.0) | ||
| v1 = -e * ny * cos(phi_alv) / rho | ||
| v2 = e * nx * cos(phi_alv) / rho | ||
| v3 = e * sin(phi_alv) / rho | ||
| p = 1.0 | ||
| B1 = nx - rho * v1 * sqr | ||
| B2 = ny - rho * v2 * sqr | ||
| B3 = -rho * v3 * sqr | ||
| psi = 0.0 | ||
|
|
||
| return prim2cons(SVector(rho, v1, v2, v3, p, B1, B2, B3, psi), equations) | ||
| end | ||
|
|
||
| initial_condition = initial_condition_constant_alfven | ||
|
|
||
| semi = SemidiscretizationHyperbolicParabolic(mesh, (equations, equations_parabolic), | ||
| initial_condition, solver) | ||
|
|
||
|
|
||
| ############################################################################### | ||
| # ODE solvers, callbacks etc. | ||
|
|
||
| # Create ODE problem with time span `tspan` | ||
| tspan = (0.0, 0.1) | ||
SimonCan marked this conversation as resolved.
Show resolved
Hide resolved
|
||
| ode = semidiscretize(semi, tspan) | ||
|
|
||
| # At the beginning of the main loop, the SummaryCallback prints a summary of the simulation setup | ||
| # and resets the timers | ||
| summary_callback = SummaryCallback() | ||
|
|
||
| # The AnalysisCallback allows to analyse the solution in regular intervals and prints the results | ||
| analysis_callback = AnalysisCallback(semi, interval = 100) | ||
|
|
||
| # The SaveSolutionCallback allows to save the solution to a file in regular intervals | ||
| save_solution = SaveSolutionCallback(interval = 100, | ||
| solution_variables = cons2prim) | ||
|
|
||
| # The StepsizeCallback handles the re-calculation of the maximum Δt after each time step | ||
| cfl = 0.5 | ||
| stepsize_callback = StepsizeCallback(cfl = cfl) | ||
|
|
||
| glm_speed_callback = GlmSpeedCallback(glm_scale = 0.5, cfl = cfl) | ||
|
|
||
| callbacks = CallbackSet(summary_callback, | ||
| analysis_callback, | ||
| save_solution, | ||
| stepsize_callback, | ||
| glm_speed_callback) | ||
|
|
||
| ############################################################################### | ||
| # run the simulation | ||
|
|
||
| sol = solve(ode, CarpenterKennedy2N54(williamson_condition = false), | ||
| dt = 1e-5, # solve needs some value here but it will be overwritten by the stepsize_callback | ||
| save_everystep = false, callback = callbacks); | ||
|
|
||
| # Print the timer summary. | ||
| summary_callback() | ||
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,148 @@ | ||
| using OrdinaryDiffEq | ||
| using Trixi | ||
|
|
||
| ############################################################################### | ||
| # semidiscretization of the visco-resistive compressible MHD equations | ||
|
|
||
| prandtl_number() = 0.72 | ||
| mu_const = 0e-2 | ||
| eta_const = 0e-2 | ||
|
|
||
| equations = IdealGlmMhdEquations3D(5 / 3) | ||
| equations_parabolic = ViscoResistiveMhd3D(equations, mu = mu_const, | ||
| Prandtl = prandtl_number(), | ||
| eta = eta_const, | ||
| gradient_variables = GradientVariablesPrimitive()) | ||
|
|
||
| volume_flux = (flux_hindenlang_gassner, flux_nonconservative_powell) | ||
| solver = DGSEM(polydeg = 3, | ||
| surface_flux = (flux_lax_friedrichs, flux_nonconservative_powell), | ||
| volume_integral = VolumeIntegralFluxDifferencing(volume_flux)) | ||
|
|
||
| coordinates_min = (-1.0, -1.0, -1.0) # minimum coordinates (min(x), min(y), min(z)) | ||
| coordinates_max = (1.0, 1.0, 1.0) # maximum coordinates (max(x), max(y), max(z)) | ||
|
|
||
| # Create a uniformly refined mesh | ||
| mesh = TreeMesh(coordinates_min, coordinates_max, | ||
| initial_refinement_level = 2, | ||
| n_cells_max = 200_000) # set maximum capacity of tree data structure | ||
|
|
||
| function initial_condition_constant_alfven(x, t, equations) | ||
| # Alfvén wave in three space dimensions modified by a periodic density variation. | ||
| # For the system without the density variations see: Altmann thesis http://dx.doi.org/10.18419/opus-3895. | ||
| # Domain must be set to [-1, 1]^3, γ = 5/3. | ||
| omega = 2.0 * pi # may be multiplied by frequency | ||
| # r = length-variable = length of computational domain | ||
| r = 2.0 | ||
| # e = epsilon | ||
| e = 0.02 | ||
| nx = 1 / sqrt(r^2 + 1.0) | ||
| ny = r / sqrt(r^2 + 1.0) | ||
| sqr = 1.0 | ||
| Va = omega / (ny * sqr) | ||
| phi_alv = omega / ny * (nx * (x[1] - 0.5 * r) + ny * (x[2] - 0.5 * r)) - Va * t | ||
|
|
||
| rho = 1.0 + e*cos(phi_alv + 1.0) | ||
| v1 = -e * ny * cos(phi_alv) / rho | ||
| v2 = e * nx * cos(phi_alv) / rho | ||
| v3 = e * sin(phi_alv) / rho | ||
| p = 1.0 | ||
| B1 = nx - rho * v1 * sqr | ||
| B2 = ny - rho * v2 * sqr | ||
| B3 = -rho * v3 * sqr | ||
| psi = 0.0 | ||
|
There was a problem hiding this comment. This looks like a fairly complicated manufactured solution that would be the ideal MHD Alfvén wave (so no source term required if viscosity and resistivity are set to 0). If I may propose an alternative, here is an MMS from my old (2D) Fortran code: alpha = 2.0_RP*PI*(x+y)-4.0_RP*t
Resu(1) = SIN(alpha)+4.0_RP
Resu(2) = SIN(alpha)+4.0_RP
Resu(3) = SIN(alpha)+4.0_RP
Resu(4) = 0.0_RP
Resu(5) = 2.0_RP*(SIN(alpha)+4.0_RP)**2
Resu(6) = SIN(alpha)+4.0_RP
Resu(7) = -SIN(alpha)-4.0_RP
Resu(8) = 0.0_RP
Resu(9) = 0.0_RPThe associated source term is then (fairly) simple and takes the form: alpha = 2.0_RP*PI*(x+y)-4.0_RP*t
phi = SIN(alpha)+4.0_RP
phi_t = -4.0_RP*COS(alpha)
phi_x = 2.0_RP*PI*COS(alpha)
phi_xx = -4.0_RP*PI**2*SIN(alpha)
source(1) = phi_t+2.0_RP*phi_x
source(2) = phi_t+4.0_RP*phi*phi_x+phi_x
source(3) = source(2)
source(4) = 0.0_RP
source(5) = 4.0_RP*phi*phi_t+16.0_RP*phi*phi_x-2.0_RP*phi_x-4.0_RP*eta*(phi_x**2+phi*phi_xx)-4.0_RP*mu/Prandtl*phi_xx
source(6) = source(1)-2.0_RP*eta*phi_xx
source(7) = -source(6)
source(8) = 0.0_RP
source(9) = 0.0_RPMaybe a somewhat simpler MMS for the 3D case could help in this instance. There was a problem hiding this comment. What are domain and final time for this case? Domain probably There was a problem hiding this comment. Yes, sorry, the domain in [0,1]^2 and the final time I typically used 1.5. The There was a problem hiding this comment. Thanks for that. I had a go with this and tried some modifications, but the convergence rate is close to 0, even without viscosity or magnetic resistivity. I am running to time 0.1, though, as it takes rather long. |
||
|
|
||
| # k = 2*pi | ||
| # rho = 1.0 | ||
| # v1 = 0 | ||
| # v2 = -e*sin(k*x[1])*sqrt(rho) | ||
| # v3 = 0 | ||
| # B1 = 1 | ||
| # B2 = e*sin(k*x[1]) | ||
| # B3 = 0 | ||
| # p = 1 | ||
| # psi = 0 | ||
|
|
||
| return prim2cons(SVector(rho, v1, v2, v3, p, B1, B2, B3, psi), equations) | ||
| end | ||
|
|
||
|
|
||
| @inline function source_terms_mhd_convergence_test(u, x, t, equations) | ||
|
There was a problem hiding this comment. From the looks of it there source term contribution were determined with computer algebra software (maybe Maple?) and then copied into here. It is hard to parse what exactly is going on. From the initial condition routine it seems that you are using the periodic manufactured solution (MMS) from Altmann's thesis. The energy source term is especially ugly. Would it be possible that you pick a simpler MMS ansatz and then compute the source term again (even by hand if necessary)? This is what I did for the compressible Navier-Stokes MMS test. There was a problem hiding this comment. Yes, that result is from using sympy. I apply the differential operator from the PDE (visco-resistive MHD) on the manufactured solution. There are two test functions. One is a simple Alfven wave in the x-direction and the other is indeed from Altmann's thesis, which is also an Alfven wave. The simpler ansatz is the Alfven wave in x, which fails in the energy term when mu and eta > 0. |
||
| r_1 = 0.02*sqrt(5)*pi*sin(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) | ||
|
|
||
| r_2 = sqrt(5)*pi^2*mu_const*(-0.04*(0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)^2*cos(pi*(sqrt(5)*t - x[1] - 2*x[2] +3.0)) + (0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)*(0.0012*cos(-2*sqrt(5)*pi*t + 2*pi*x[1] + 4*pi*x[2] + 1) - 0.0004*cos(1)) + 3.2e-5*sin(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1)^2*cos(pi*(sqrt(5)*t - x[1] - 2*x[2] + 3.0)))/(0.02*cos(-sqrt(5)*pi*t + pi*(x[1]+ 2*x[2] - 3.0) + 1) + 1)^3 | ||
|
|
||
| r_3 = -sqrt(5)*pi^2*mu_const*(-0.02*(0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)^2*cos(pi*(sqrt(5)*t - x[1] - 2*x[2] + 3.0)) + (0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)*(0.0006*cos(-2*sqrt(5)*pi*t + 2*pi*x[1] + 4*pi*x[2] + 1) - 0.0002*cos(1)) + 1.6e-5*sin(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1)^2*cos(pi*(sqrt(5)*t - x[1] - 2*x[2] + 3.0)))/(0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)^3 | ||
|
|
||
| r_4 = -pi^2*mu_const*((0.003*sin(-2*sqrt(5)*pi*t + 2*pi*x[1] + 4*pi*x[2] + 1) + 0.001*sin(1))*(0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1) + 0.1*(0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)^2*sin(pi*(-sqrt(5)*t + x[1] + 2*x[2])) - 8.0e-5*sin(pi*(-sqrt(5)*t + x[1] + 2*x[2]))*sin(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1)^2)/(0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)^3 | ||
|
|
||
| r_5 = 0.2*pi*(-0.00138888888888889*pi*mu_const^2*sin(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1)^2*cos(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1) + 0.0694444444444444*pi*mu_const^2*sin(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1)^2 - 0.000694444444444444*pi*mu_const^2*cos(-sqrt(5)*pi*t+ pi*x[1] + 2*pi*x[2] + 1)^3 + 0.0694444444444444*pi*mu_const^2*cos(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1)^2 - 1.73611111111111*pi*mu_const^2*cos(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1) - 1.2e-5*pi*mu_const*sin(pi*(-sqrt(5)*t + x[1] + 2*x[2]))^2*sin(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1)^2 - 4.0e-6*pi*mu_const*sin(pi*(-sqrt(5)*t + x[1] + 2*x[2]))^2*cos(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1)^2 + 0.0002*pi*mu_const*sin(pi*(-sqrt(5)*t + x[1] + 2*x[2]))^2*cos(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1) - 1.2e-5*pi*mu_const*sin(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1)^2*cos(pi*(sqrt(5)*t - x[1] - 2*x[2] + 3.0))^2 - 4.0e-6*pi*mu_const*cos(pi*(sqrt(5)*t - x[1] - 2*x[2]+ 3.0))^2*cos(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1)^2 + 0.0002*pi*mu_const*cos(pi*(sqrt(5)*t - x[1] - 2*x[2] + 3.0))^2*cos(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1) + 1.0e-7*sqrt(5)*(-sin(-4*sqrt(5)*pi*t + 4*pi*x[1] + 8*pi*x[2] + 2) + 2*sin(-2*sqrt(5)*pi*t + 2*pi*x[1] + 4*pi*x[2] + 2) - sin(2)) + 1.0e-7*sqrt(5)*(sin(-4*sqrt(5)*pi*t + 4*pi*x[1] + 8*pi*x[2] + 2) + 2*sin(-2*sqrt(5)*pi*t + 2*pi*x[1] + 4*pi*x[2] + 2) + sin(2)) - 8.0e-9*sqrt(5)*sin(pi*(-sqrt(5)*t + x[1] + 2*x[2]))^2*sin(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1)*cos(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1)^2 - 2.0e-5*sqrt(5)*sin(pi*(-sqrt(5)*t + x[1] + 2*x[2]))^2*sin(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1) - 8.0e-9*sqrt(5)*sin(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1)*cos(pi*(sqrt(5)*t - x[1] - 2*x[2] + 3.0))^2*cos(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1)^2 - 2.0e-5*sqrt(5)*sin(-sqrt(5)*pi*t + pi*x[1] + 2*pi*x[2] + 1)*cos(pi*(sqrt(5)*t - x[1] - 2*x[2] + 3.0))^2)/(-0.04*sin(pi*x[2])*sin(-sqrt(5)*pi*t + pi*x[1] + pi*x[2] + 1) + 0.02*cos(-sqrt(5)*pi*t + pi*x[1] + 1) - 1.0)^4 | ||
|
|
||
| r_6 = pi*(-0.04*(sqrt(5)*pi*eta_const*cos(pi*(-sqrt(5)*t + x[1] + 2*x[2])) + sin(pi*(-sqrt(5)*t + x[1] + 2*x[2])))*(0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)^2 + 0.04*(0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)*sin(pi*(-sqrt(5)*t + x[1] + 2*x[2])) +0.0004*sin(-2*sqrt(5)*pi*t + 2*pi*x[1] + 4*pi*x[2] + 1) + 0.0004*sin(1))/(0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)^2 | ||
|
|
||
| r_7 = pi*(0.02*(sqrt(5)*pi*eta_const*cos(pi*(-sqrt(5)*t + x[1] + 2*x[2])) + sin(pi*(-sqrt(5)*t + x[1] + 2*x[2])))*(0.02*cos(-sqrt(5)*pi*t +pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)^2 - 0.02*(0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)*sin(pi*(-sqrt(5)*t + x[1] + 2*x[2])) - 0.0002*sin(-2*sqrt(5)*pi*t + 2*pi*x[1] + 4*pi*x[2] + 1) - 0.0002*sin(1))/(0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)^2 | ||
|
|
||
| r_8 = pi*((0.1*pi*eta_const*sin(pi*(-sqrt(5)*t + x[1] + 2*x[2])) + 0.02*sqrt(5)*cos(pi*(sqrt(5)*t - x[1] - 2*x[2] + 3.0)))*(0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)^2 - 0.02*sqrt(5)*(0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)*cos(pi*(sqrt(5)*t - x[1] - 2*x[2] + 3.0)) + 0.0002*sqrt(5)*(cos(-2*sqrt(5)*pi*t + 2*pi*x[1] + 4*pi*x[2] + 1) - cos(1)))/(0.02*cos(-sqrt(5)*pi*t + pi*(x[1] + 2*x[2] - 3.0) + 1) + 1)^2 | ||
|
|
||
| r_9 = 0.0 | ||
|
|
||
| # r_1 = 0.0 | ||
| # r_2 = 0.0004*pi*sin(4*pi*x[1]) | ||
| # r_3 = pi*(0.08*pi*mu_const*sin(2*pi*x[1]) - 0.04*cos(2*pi*x[1])) | ||
| # r_4 = 0 | ||
| # r_5 = pi*(-0.0016*pi*eta_const*sin(2*pi*x[1])^2 + 0.0016*pi*eta_const*cos(2*pi*x[1])^2 + pi*mu_const*(0.0016 - 0.0111111111111111*mu_const)*cos(4*pi*x[1]) + 0.0008*sin(4*pi*x[1])) | ||
| # r_6 = 0 | ||
| # r_7 = pi*(-0.08*pi*eta_const*sin(2*pi*x[1]) + 0.04*cos(2*pi*x[1])) | ||
| # r_8 = 0 | ||
| # r_9 = 0.0 | ||
|
|
||
| return SVector(r_1, r_2, r_3, r_4, r_5, r_6, r_7, r_8, r_9) | ||
| end | ||
|
|
||
|
|
||
| initial_condition = initial_condition_constant_alfven | ||
|
|
||
| semi = SemidiscretizationHyperbolicParabolic(mesh, (equations, equations_parabolic), | ||
| initial_condition, solver; | ||
| source_terms = source_terms_mhd_convergence_test) | ||
|
|
||
|
|
||
| ############################################################################### | ||
| # ODE solvers, callbacks etc. | ||
|
|
||
| # Create ODE problem with time span `tspan` | ||
| tspan = (0.0, 0.1) | ||
| ode = semidiscretize(semi, tspan) | ||
|
|
||
| # At the beginning of the main loop, the SummaryCallback prints a summary of the simulation setup | ||
| # and resets the timers | ||
| summary_callback = SummaryCallback() | ||
|
|
||
| # The AnalysisCallback allows to analyse the solution in regular intervals and prints the results | ||
| analysis_callback = AnalysisCallback(semi, interval = 100) | ||
|
|
||
| # The SaveSolutionCallback allows to save the solution to a file in regular intervals | ||
| save_solution = SaveSolutionCallback(interval = 100, | ||
| solution_variables = cons2prim) | ||
|
|
||
| # The StepsizeCallback handles the re-calculation of the maximum Δt after each time step | ||
| cfl = 0.5 | ||
| stepsize_callback = StepsizeCallback(cfl = cfl) | ||
|
|
||
| glm_speed_callback = GlmSpeedCallback(glm_scale = 0.5, cfl = cfl) | ||
|
|
||
| callbacks = CallbackSet(summary_callback, | ||
| analysis_callback, | ||
| save_solution, | ||
| stepsize_callback, | ||
| glm_speed_callback) | ||
|
|
||
| ############################################################################### | ||
| # run the simulation | ||
|
|
||
| sol = solve(ode, CarpenterKennedy2N54(williamson_condition = false), | ||
| dt = 1e-5, # solve needs some value here but it will be overwritten by the stepsize_callback | ||
| save_everystep = false, callback = callbacks); | ||
|
|
||
| # Print the timer summary. | ||
| summary_callback() | ||
Oops, something went wrong.
Oops, something went wrong.
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.