[WIP] Added Runge Kutta 3/8 rule particle integration scheme. #5841
Conversation
…eased error at a soon to be determined greater computational cost.
|
You will have to run |
There was a problem hiding this comment.
I think the rest is ok, though I haven't checked the formulas. It would be good to copy and adapt one of the existing tests to use this scheme, such as tests/particle_integrator_rk4.prm. You should hopefully get essentially identical output to the existing test.
bangerth
changed the title
There was a problem hiding this comment.
Thank for the contribution! Given the results of your benchmark I do not know if you want to pursue this PR further, but if you want, here are some comments.
This PR would need:
- a changelog entry about rk4_38
- a test that uses rk4_38
- a modification to the particle_integration_scheme benchmark that includes rk4_38 in the discussion
| "rk2_overhead" > $output_file | ||
|
|
||
| for refinement in 2 3 4 5; do # | ||
| # velocity_accuracy_colum=`grep "u_L2" refinement_${refinement}_higher_order_true_interpolation_bilinear_least_squares/statistics | cut -d ' ' -f 2 | head -c -2` |
There was a problem hiding this comment.
Please remove all commented code from this script (I likely gave you this script in this shape). Commented code is a hazard for future users, because it is not clear if it is old code that should be removed, future code to be implemented, or just comments explaining things.
There was a problem hiding this comment.
Did you mean to include this file in the pull request? I think it would be intended for a future PR that also measures the performance.
| @@ -0,0 +1,5 @@ | |||
| subsection Postprocess | |||
There was a problem hiding this comment.
since you add this option to the particle_integration_scheme benchmark, please also add the results and maybe a discussion to the documentation of this benchmark.
| template <int dim> | ||
| void | ||
| RK4_38<dim>::local_integrate_step(const typename ParticleHandler<dim>::particle_iterator &begin_particle, | ||
| const typename ParticleHandler<dim>::particle_iterator &end_particle, |
There was a problem hiding this comment.
You need to run make indent to fix the indentation.
| #endif | ||
| k[0] = dt * (*old_velocity); | ||
|
|
||
| Point<dim> new_location = location + 1./3. * k[0]; |
There was a problem hiding this comment.
Here is a little trick that can become important, because particle advection happens so often: It is much faster for a CPU to do multiplication than division. Because you use 1./3. and 2./3. so often and the result is always the same, it can be useful to exclude it from the loop over all particles and compute it once in line 88:
const double one_thirds = 1./3.;
const double two_thirds = 2./3.;
Then you can do here and elsewhere:
| Point<dim> new_location = location + 1./3. * k[0]; | |
| Point<dim> new_location = location + one_thirds * k[0]; |
|
Hi Rene, Thank you for your comment. Given the benchmark results I am inclined to agree however, I will run one more test before closing this pull request just to make sure. Best |
|
Based on the results of the benchmarks lets close this PR for now, the higher cost does not seem worth the marginally more accurate result. In case you want to test other integrations schemes in the future feel free to reopen this PR or open a new one. Thanks for the test. |
[WIP] Added Runga Kutta 3/8 rule particle integration scheme for decreased error at a soon to be determined greater computational cost.