Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add maximum number of iteration steps in Propagator::AdvanceParticle #338

Merged
merged 2 commits into from
Jan 24, 2023

Conversation

Jean1995
Copy link
Member

Issue #332 shows that the propagation can stall if the PARTICLE_POSITION_RESOLUTION can not be reached due to numerical uncertainties (for very large propagation steps, complicated medium intersections, non-dense media, with multiple scattering enabled). I've been trying to find a "good" way to solve this problem, but was not successful to find a solution that made me happy.

As a "quick fix", this PR adds an (adjustable) number of maximum iteration steps in Propagator::AdvanceParticle.
After we reached the maximum number of propagation steps, we check the deviation from PARTICLE_POSITION_RESOLUTION is. If it is "still acceptable" (which is currently defined as "1 cm"), we still continue propagation. This avoids that the propagation always stops if we have the problem in issue #332. If the deviation is "not acceptable", we throw an error/exception.

Furthermore, I found a small issue with the "backscattering" (issue #267), which is solved by resampling random numbers if we run into this case twice.

I have to admit that I am still not happy with the solution, but I believe this PR provides a quickfix that doesn't complicate the algorithm even further, and it should avoid that the algorithm stalls in case someone runs into this problem. For now, I don't think that many users will run into issue #332, but if they do, they will now get a reasonable error message that we can help to debug.
At some point in the future, it is going to be necessary to revise the whole Propagator::AdvanceParticle algorithm again.

…agation if we have a discrepancy of above 1cm, otherwise log and continue.

Introduce resampling of random numbers in "backscattering" case.
@sudojan
Copy link
Contributor

sudojan commented Jan 16, 2023

Do I understand it correctly, that the introduced maximum number of iteration steps is the limit for all advancing steps in the propagation, regardless if the particle actually gets advanced or not? I suggest introducing an iteration counter for steps, where there is no progress. If e.g. the algorithm to find the correct multiple scattering does not converge (which is I think the only case where there is no progress), the number of iterations should be limited to 20, not 200.

However, I think introducing such a general iteration limit for the steps where the particle gets advanced does make sense, whether there is a bug or not. However, I would increase the default value to at least 1000, or better 10000.

@Jean1995
Copy link
Member Author

Do I understand it correctly, that the introduced maximum number of iteration steps is the limit for all advancing steps in the propagation, regardless if the particle actually gets advanced or not? I suggest introducing an iteration counter for steps, where there is no progress. If e.g. the algorithm to find the correct multiple scattering does not converge (which is I think the only case where there is no progress), the number of iterations should be limited to 20, not 200.

However, I think introducing such a general iteration limit for the steps where the particle gets advanced does make sense, whether there is a bug or not. However, I would increase the default value to at least 1000, or better 10000.

The AdvanceParticle method is called for every propagation step, and we only leave the method as soon as we found a valid propagation step (either towards the next interaction, or towards a sector transition, or towards the min energy / max distance). The counter is for the iteration process of a single propagation step. In normal cases, there are like 5 iterations. However, I had cases where the iteration correctly converged after ~100 steps. Having 20 as a default would have raised an exception. This probably indicates that the "minimalization" process needs to be improved, but this requires additional thought/work.

@sudojan
Copy link
Contributor

sudojan commented Jan 17, 2023

Do I understand it correctly, that the introduced maximum number of iteration steps is the limit for all advancing steps in the propagation, regardless if the particle actually gets advanced or not? I suggest introducing an iteration counter for steps, where there is no progress. If e.g. the algorithm to find the correct multiple scattering does not converge (which is I think the only case where there is no progress), the number of iterations should be limited to 20, not 200.
However, I think introducing such a general iteration limit for the steps where the particle gets advanced does make sense, whether there is a bug or not. However, I would increase the default value to at least 1000, or better 10000.

The AdvanceParticle method is called for every propagation step, and we only leave the method as soon as we found a valid propagation step (either towards the next interaction, or towards a sector transition, or towards the min energy / max distance). The counter is for the iteration process of a single propagation step. In normal cases, there are like 5 iterations. However, I had cases where the iteration correctly converged after ~100 steps. Having 20 as a default would have raised an exception. This probably indicates that the "minimalization" process needs to be improved, but this requires additional thought/work.

ah, so I did understand it wrong; the ADVANCE_PARTICLE_MAX_STEPS is the maximum number of iterations for a single step. Then, sorry for the misunderstanding.

But 200 iterations for a single step is still inefficient. Are the differences between the iterations steps, let's say after 10 iterations still greater than the PARTICLE_POSITION_RESOLUTION?

Maybe, if the algorithm doesn't converge, PROPOSAL should raise a warning instead of an Exception, or did we already discuss this concluding with this behaviour?

@Jean1995
Copy link
Member Author

Yes, I agree that 200 iterations is inefficient.
Under normal circumstances, and for normal propagation settings, the algorithm converges must faster, even faster than 10 iterations. However, with the very special propagation conditions in my example (air, very large steps, intersection under a large angle), I had cases where after 10 iterations, the differences were still much bigger than PARTICLE_POSITION_RESOLUTION. It might be interesting to find out statistics about the number of iterations that are needed, but they are going to vary a lot for different conditions.

In this PR, if the algorithm does not converge, we check the remaining difference: If it is smaller than 1 cm (this value is currently more or less arbitrary), we ignore this difference and just issue a logging info. Only if the difference is bigger, we throw the exception. We might introduce a different behavior here, like a warning, and just continue the propagation.

@sudojan
Copy link
Contributor

sudojan commented Jan 18, 2023

Only if the difference is bigger, we throw the exception. We might introduce a different behavior here, like a warning, and just continue the propagation.

I think when reaching the 200 iterations, PROPOSAL should always give a warning with the remaining uncertainty but should continue to propagate.

@Jean1995
Copy link
Member Author

Only if the difference is bigger, we throw the exception. We might introduce a different behavior here, like a warning, and just continue the propagation.

I think when reaching the 200 iterations, PROPOSAL should always give a warning with the remaining uncertainty but should continue to propagate.

Ok. I think a valid idea would be to just continue propagating (and giving a warning) no matter the remaining uncertainty. But in addition, if we have reached the iteration limit, we will recalculate energy and grammage with the current (and therefore accepted) actual distance_to_border. This way, the only "bias" that we may introduce is in the multiple scattering direction (which might be under/overestimated), but everything else will be fine.

…cle to the border with current scattering angles.
@Jean1995 Jean1995 merged commit 1bb15a0 into master Jan 24, 2023
@Jean1995 Jean1995 deleted the fix_advance_particles branch January 24, 2023 19:44
@Jean1995 Jean1995 mentioned this pull request Jan 25, 2023
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants