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

WIP: Some benchmarks for the documentation #123

Draft
wants to merge 69 commits into
base: main
Choose a base branch
from
Draft

Conversation

SKopecz
Copy link
Owner

@SKopecz SKopecz commented Sep 6, 2024

This will add some benchmarks for #7.

  • NPZD model
  • Robertson
  • Stratospheric reaction

@codecov-commenter
Copy link

codecov-commenter commented Sep 6, 2024

⚠️ Please install the 'codecov app svg image' to ensure uploads and comments are reliably processed by Codecov.

Codecov Report

Attention: Patch coverage is 95.60440% with 4 lines in your changes missing coverage. Please review.

Project coverage is 98.05%. Comparing base (6ca36f5) to head (da96e89).

Files with missing lines Patch % Lines
src/utilities.jl 95.60% 4 Missing ⚠️

❗ Your organization needs to install the Codecov GitHub app to enable full functionality.

Additional details and impacted files
@@            Coverage Diff             @@
##             main     #123      +/-   ##
==========================================
- Coverage   98.26%   98.05%   -0.21%     
==========================================
  Files           7        7              
  Lines        1556     1647      +91     
==========================================
+ Hits         1529     1615      +86     
- Misses         27       32       +5     

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

@SKopecz SKopecz marked this pull request as draft September 6, 2024 10:59
@SKopecz SKopecz changed the title Some benchmarks for the documentation WIP: Some benchmarks for the documentation Sep 6, 2024
#```@example eoc
using PrettyTables

# collect data and create headers
Copy link
Owner Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Coding this almost drove me crazy!

Originally I used something like data = dts outside the loop and data = [data err[i] [NaN; eoc[i]]] within the loop. This was perfectly fine on my local machine, but caused an error like "data not defined" in the build process.

@ranocha: Any idea what the problem was/is?

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Did you built the docs locally or copy and paste the code to the REPL?

Copy link
Owner Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Oh. I didn't build the docs locally, I copied & pasted everything from my vscode repl to the markdown files. But I ensured to use the docs environment.

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Then it might be an issue of soft global scope, see https://docs.julialang.org/en/v1/manual/variables-and-scoping/#on-soft-scope

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You may need to add a keyword global

@SKopecz
Copy link
Owner Author

SKopecz commented Sep 9, 2024

@ranocha: This is still WIP, but do you have any thoughts on current state of these benchmarks?

@ranocha: I played around with the WorkPrecisionSet from DiffEqDevTools. I think its not flexible enough. So I'm gonna redo all work-precision diagrams in a way similar to what you did here: https://github.com/JoshuaLampert/2024_fokker_planck

I'm also going to add WP diagrams for constant time step sizes.

Copy link
Collaborator

@ranocha ranocha left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Interesting, thanks! The MPRK methods do not look very good, it seems...

Comment on lines 63 to 64
formatters = (ft_printf("%5.4e", [1, 2, 4, 6, 8]),
ft_printf("%5.4f", [3, 5, 7, 9])))
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think we don't need all of the digits - this would also keep the tables smaller so that we don't have to scroll as much. And maybe we can just show "EOC" right next to the value so that the headers in the first row are not duplicated.

# select problem
prob = prob_pds_npzd

# compute reference solution (standard tolerances are too low)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do you mean not strict enough? The question is if you consider 1e-5 to be lower (smaller) than 1e-13

docs/src/npzd_model_benchmark.md Outdated Show resolved Hide resolved
Comment on lines 61 to 68
labels = ["MPRK22(0.5)"
"MPPRK22(2/3)"
"MPRK22(1.0)"
"SSPMPRK22(0.5,1.0)"
"MPRK43I(1.0,0.5)"
"MPRK43I(0.5,0.75)"
"MPRK43II(0.5)"
"MPRK43II(2.0/3.0)"]
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This pattern seems to be present several times. Maybe we shall overload show to automate this? But we would likely lose some information so I am not fully sure...

Although the SSPMPRK22 solution seems to be more accurate at first glance, the `l∞`-error of the SSPMPRK22 scheme is 0.506359, whereas the `l∞`-error of the MPRK43 scheme is 0.413915. Both errors occurs at approximately $t=2$, where there is a sharp kink in the first component.


Next we compare `SSPMPRK22(0.5, 1.0)` and `MPRK43I(1.0, 0.5)` with some second and third order methods from [OrdinaryDiffEq.jl](https://docs.sciml.ai/OrdinaryDiffEq/stable/). To guarantee positive solutions with these methods, we must select the solver option `isoutofdomain = isnegative`.
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is isoutofdomain = isnegative to get good work-precision diagrams? How would they look like without it?

# set tolerances and error
abstols = 1.0 ./ 10.0 .^ (2:0.5:8)
reltols = 1.0 ./ 10.0 .^ (1:0.5:7)
err_est = :l∞
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do other error choices lead to significantly different results?

plot!(p1, sol1; tspan = (1e-6, 1e11), xaxis = :log, denseplot = false, markers = :circle, ylims = (-0.2, 1.2), idxs = [(0, 1), ((x, y) -> (x, 1e4 .* y), 0, 2), (0, 3)], title = "SSPMPRK22(1.0, 0.5)", xticks =10.0 .^ (-6:4:10))
```

With `abstol=1e-6` and `reltol = 1e-5` the `MPRK22(0.5)` schemes needs over 800.000 steps to integrate the Robertson problem. In comparison, `MPRK22(0.5)` needs less than 1.000.
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
With `abstol=1e-6` and `reltol = 1e-5` the `MPRK22(0.5)` schemes needs over 800.000 steps to integrate the Robertson problem. In comparison, `MPRK22(0.5)` needs less than 1.000.
With `abstol=1e-6` and `reltol = 1e-5` the `MPRK22(0.5)` schemes needs over 800.000 steps to integrate the Robertson problem.
In comparison, `MPRK22(1.0)` needs less than 1.000.

Is this related to solutions close to zero and the issues we discovered with Davide and Philipp? Or due to maybe a bad error estimator?

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Relative errors would make sense here as well

Copy link
Owner Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, especially for the stratospheric reaction problem.

@ranocha
Copy link
Collaborator

ranocha commented Sep 9, 2024

Maybe it would also make sense to check the controllers Thomas optimized in https://arxiv.org/abs/2312.01796

@SKopecz
Copy link
Owner Author

SKopecz commented Nov 4, 2024

In the Robertson benchmark tutorial we see

plot_3_web

But on my local machine the situation is different. With abstol=1e-6 and reltol = 1e-5 the MPRK22(0.5) scheme needs over 800.000 steps to integrate the Robertson problem and the corresponding error is large. See the following picture.

plot_3

When I originally created the benchmark tutorial I excluded MPRK22(0.5) from the comparison, but now it's suddenly the best scheme. Of course, this is some strange exception.

In addition, the question arises if schemes that deliver good results in the benchmark tutorials perform significantly worse on other machines?

I suggest to keep MPRK22(0.5) in the benchmark and to mention that this superior performance is not visible in the other benchmarks.

@ranocha: Any thoughts on this?

@SKopecz
Copy link
Owner Author

SKopecz commented Nov 4, 2024

@ranocha The CI tests fail because there are issues with Aqua.jl and ExplicitImports.jl. Do you know what is going wrong there?

@ranocha
Copy link
Collaborator

ranocha commented Nov 5, 2024

I suggest to keep MPRK22(0.5) in the benchmark and to mention that this superior performance is not visible in the other benchmarks.

Sounds reasonable to me 👍

test/runtests.jl Outdated
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Since you're using the median below, you should change the line

using Statistics: mean

in line 4 to

using Statistics: mean, median

I can't comment on the line directly

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.

3 participants