Use `VectorOfArray` in `wrap_array` for `DGMulti` solvers #2150

jlchan · 2024-11-07T15:29:42Z

@huiyuxie FYI this PR goes towards addressing #1789.

d

github-actions · 2024-11-07T15:30:00Z

Review checklist

This checklist is meant to assist creators of PRs (to let them know what reviewers will typically look for) and reviewers (to guide them in a structured review process). Items do not need to be checked explicitly for a PR to be eligible for merging.

Purpose and scope

The PR has a single goal that is clear from the PR title and/or description.
All code changes represent a single set of modifications that logically belong together.
No more than 500 lines of code are changed or there is no obvious way to split the PR into multiple PRs.

Code quality

The code can be understood easily.
Newly introduced names for variables etc. are self-descriptive and consistent with existing naming conventions.
There are no redundancies that can be removed by simple modularization/refactoring.
There are no leftover debug statements or commented code sections.
The code adheres to our conventions and style guide, and to the Julia guidelines.

Documentation

New functions and types are documented with a docstring or top-level comment.
Relevant publications are referenced in docstrings (see example for formatting).
Inline comments are used to document longer or unusual code sections.
Comments describe intent ("why?") and not just functionality ("what?").
If the PR introduces a significant change or new feature, it is documented in NEWS.md with its PR number.

Testing

The PR passes all tests.
New or modified lines of code are covered by tests.
New or modified tests run in less then 10 seconds.

Performance

There are no type instabilities or memory allocations in performance-critical parts.
If the PR intent is to improve performance, before/after time measurements are posted in the PR.

Verification

The correctness of the code was verified using appropriate tests.
If new equations/methods are added, a convergence test has been run and the results
are posted in the PR.

Created with ❤️ by the Trixi.jl community.

codecov · 2024-11-07T16:12:17Z

Codecov Report

All modified and coverable lines are covered by tests ✅

Project coverage is 95.44%. Comparing base (91eaaf6) to head (9b21e12).
Report is 2 commits behind head on main.

Additional details and impacted files

@@            Coverage Diff             @@
##             main    #2150      +/-   ##
==========================================
- Coverage   96.36%   95.44%   -0.92%     
==========================================
  Files         480      480              
  Lines       38028    38025       -3     
==========================================
- Hits        36645    36292     -353     
- Misses       1383     1733     +350

Flag	Coverage Δ
unittests	`95.44% <ø> (-0.92%)`	⬇️

Flags with carried forward coverage won't be shown. Click here to find out more.

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

Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>

jlchan · 2024-11-07T23:09:17Z

@DanielDoehring if I remember correctly, you implemented @allocated tests, right?

I am running into some @allocated failures (for example https://github.com/trixi-framework/Trixi.jl/actions/runs/11732466598/job/32684751742?pr=2150)

elixir_euler_weakform.jl (SBP, EC): Test Failed at /home/runner/work/Trixi.jl/Trixi.jl/test/test_threaded.jl:388
  Expression: #= /home/runner/work/Trixi.jl/Trixi.jl/test/test_threaded.jl:388 =# @allocated(Trixi.rhs!(du_ode, u_ode, semi, t)) < 5000
   Evaluated: 118704 < 5000

However, at least locally, this appears to only happen the first time Trixi.rhs! is run. I tried running Trixi.rhs! to avoid these allocations in CI, but it doesn't seem to work.

Did you run into this issue before?

Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>

Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com> Update test/test_dgmulti_2d.jl Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com> Apply suggestions from code review Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>

…ework/Trixi.jl into jc/wrap_VectorOfArray

DanielDoehring · 2024-11-08T07:36:13Z

However, at least locally, this appears to only happen the first time Trixi.rhs! is run. I tried running Trixi.rhs! to avoid these allocations in CI, but it doesn't seem to work.

Did you run into this issue before?

Hm, I have not encountered this. Can you maybe track down the allocations using the SummaryCallback ? I guess there must be some type instability. But I have no idea why the behaviour after executing rhs! once is different depending on the machine.

jlchan · 2024-11-08T21:07:48Z

Hm, I have not encountered this. Can you maybe track down the allocations using the SummaryCallback ? I guess there must be some type instability. But I have no idea why the behaviour after executing rhs! once is different depending on the machine.

There doesn't seem to be a type instability. There's a weird GC call stack that @profview_allocs shows, I have no idea what's going on. I'll bring it up at the next Trixi meeting.

DanielDoehring · 2024-11-09T12:15:00Z

There doesn't seem to be a type instability. There's a weird GC call stack that @profview_allocs shows, I have no idea what's going on. I'll bring it up at the next Trixi meeting.

Oh wow, that sounds strange. So the allocations are only present in the test "environment"?

jlchan · 2024-11-10T16:05:22Z

There doesn't seem to be a type instability. There's a weird GC call stack that @profview_allocs shows, I have no idea what's going on. I'll bring it up at the next Trixi meeting.

Oh wow, that sounds strange. So the allocations are only present in the test "environment"?

Not exactly - on my machine, the allocations disappear after running rhs! once. In the test environment (or on CI, I can't tell) this isn't the case.

Here's the output from Profile.Allocs.@profile sample_rate=0.001 Trixi.rhs!(du, ode.u0, semi, 0.0) printed to a txt file: profile_allocs.txt.

Project.toml

src/auxiliary/precompile.jl

test/test_dgmulti_1d.jl

huiyuxie · 2024-11-11T03:12:23Z

https://github.com/JuliaLang/julia/blob/3318941e585db632423366b8b703ea55a6ba8421/base/timing.jl#L479-L489
If the allocation is less in the second call compared to the first during a single session, the GC might have been triggered between these two gc_bytes calls. But I also suspect it could be caused by type instability introduced by RecursiveArrayTools.jl. Why do you think it doesn’t seem to be a type instability issue?

huiyuxie · 2024-11-11T03:12:48Z

Here's the output from Profile.Allocs.@Profile sample_rate=0.001 Trixi.rhs!(du, ode.u0, semi, 0.0) printed to a txt file: profile_allocs.txt.

The text file only shows half of the information on each line. Could you provide it in a different way?

jlchan · 2024-11-11T03:37:27Z

Why do you think it doesn’t seem to be a type instability issue?

I checked @code_warntype on the first call; nothing is red and I don't see an Any instance anywhere. If you have any suggestions for where else to check, please let me know.

The text file only shows half of the information on each line. Could you provide it in a different way?

Thanks - how is this?
profile_allocs.txt

huiyuxie · 2024-11-11T03:41:38Z

Same issue - how about telling me how to reproduce it? Use your latest commit and then run Profile.Allocs.@profile sample_rate=0.001 Trixi.rhs!(du, ode.u0, semi, 0.0)?

huiyuxie · 2024-11-15T03:59:36Z

But if you said that u is the same type for both, I highly suspect there is a type instability problem with this package...

jlchan · 2024-11-15T04:15:42Z

But if you said that u is the same type for both, I highly suspect there is a type instability problem with this package...

No, I think it's a weird subtlety with how OrdinaryDiffEq.jl initializes the du container for specifically the SSPRK43 time-stepper. Here's a MWE:

using Trixi, OrdinaryDiffEq

dg = DGMulti(polydeg = 3, element_type = Line(),
             surface_integral = SurfaceIntegralWeakForm(FluxLaxFriedrichs()),
             volume_integral = VolumeIntegralWeakForm())

equations = CompressibleEulerEquations1D(1.4)
initial_condition = initial_condition_constant

mesh = DGMultiMesh(dg, (2, ), periodicity = true)
semi = SemidiscretizationHyperbolic(mesh, equations, initial_condition, dg)
ode = semidiscretize(semi, (0.0, .01))

integrator = init(ode, CarpenterKennedy2N54(williamson_condition = false), dt = 0.01) # works
integrator = init(ode, Tsit5()) # works
integrator = init(ode, RDPK3SpFSAL49()) # works
integrator = init(ode, SSPRK43(), dt = 0.01) # works

integrator = init(ode, SSPRK43()) # fails?!?

The error for the last line isn't very helpful, but if I run @show typeof(du) in Trixi.rhs!, I get typeof(du) = Matrix{SVector{3, Float64}}. It should be StructArray{SVector{3, Float64}, ...} instead, so I think OrdinaryDiffEq.jl is initializing the integrator.k field (that gets passed into Trixi.rhs! as du) incorrectly.

@ranocha - you know the SSPRK time-steppers in OrdinaryDiffEq.jl better than I do. Any idea what might be causing this?

ranocha · 2024-11-20T19:19:30Z

Which versions of the packages are you using?

julia> using Pkg; Pkg.activate(temp = true); Pkg.add(["OrdinaryDiffEq", "StaticArrays", "StructArrays"])
[...]
⌃ [1dea7af3] + OrdinaryDiffEq v6.89.0
  [90137ffa] + StaticArrays v1.9.8
⌃ [09ab397b] + StructArrays v0.6.18
[...]

julia> using OrdinaryDiffEq, StructArrays, StaticArrays

julia> u = StructArray([SVector(1.0, 1.0), SVector(2.0, 2.0)])
2-element StructArray(::Vector{Float64}, ::Vector{Float64}) with eltype SVector{2, Float64}:
 [1.0, 1.0]
 [2.0, 2.0]

julia> ode = ODEProblem(u, (0.0, 1.0)) do du, u, p, t
           @show typeof(du)
           du .= u
           return nothing
       end

julia> solve(ode, Tsit5())
typeof(du) = StructVector{SVector{2, Float64}, Tuple{Vector{Float64}, Vector{Float64}}, Int64}
typeof(du) = StructVector{SVector{2, Float64}, Tuple{Vector{Float64}, Vector{Float64}}, Int64}
typeof(du) = StructVector{SVector{2, Float64}, Tuple{Vector{Float64}, Vector{Float64}}, Int64}
[...]

julia> solve(ode, SSPRK43())
typeof(du) = StructVector{SVector{2, Float64}, Tuple{Vector{Float64}, Vector{Float64}}, Int64}
typeof(du) = StructVector{SVector{2, Float64}, Tuple{Vector{Float64}, Vector{Float64}}, Int64}
typeof(du) = StructVector{SVector{2, Float64}, Tuple{Vector{Float64}, Vector{Float64}}, Int64}
[...]

I can't see the difference here that you describe 🤷

jlchan · 2024-11-20T20:47:20Z

Which versions of the packages are you using?

julia> using Pkg; Pkg.activate(temp = true); Pkg.add(["OrdinaryDiffEq", "StaticArrays", "StructArrays"])
[...]
⌃ [1dea7af3] + OrdinaryDiffEq v6.89.0
  [90137ffa] + StaticArrays v1.9.8
⌃ [09ab397b] + StructArrays v0.6.18
[...]

julia> using OrdinaryDiffEq, StructArrays, StaticArrays

julia> u = StructArray([SVector(1.0, 1.0), SVector(2.0, 2.0)])
2-element StructArray(::Vector{Float64}, ::Vector{Float64}) with eltype SVector{2, Float64}:
 [1.0, 1.0]
 [2.0, 2.0]

julia> ode = ODEProblem(u, (0.0, 1.0)) do du, u, p, t
           @show typeof(du)
           du .= u
           return nothing
       end

julia> solve(ode, Tsit5())
typeof(du) = StructVector{SVector{2, Float64}, Tuple{Vector{Float64}, Vector{Float64}}, Int64}
typeof(du) = StructVector{SVector{2, Float64}, Tuple{Vector{Float64}, Vector{Float64}}, Int64}
typeof(du) = StructVector{SVector{2, Float64}, Tuple{Vector{Float64}, Vector{Float64}}, Int64}
[...]

julia> solve(ode, SSPRK43())
typeof(du) = StructVector{SVector{2, Float64}, Tuple{Vector{Float64}, Vector{Float64}}, Int64}
typeof(du) = StructVector{SVector{2, Float64}, Tuple{Vector{Float64}, Vector{Float64}}, Int64}
typeof(du) = StructVector{SVector{2, Float64}, Tuple{Vector{Float64}, Vector{Float64}}, Int64}
[...]

I can't see the difference here that you describe 🤷

Odd! I'll check tomorrow (Thursdays are my Trixi PR day).

jlchan · 2024-11-20T20:48:19Z

@ranocha I misread your MWE; I couldn't reproduce this with a simplified MWE either. It only occurs in Trixi.jl.

I'll try to isolate it into a MWE tomorrow.

huiyuxie · 2024-11-22T00:46:42Z

Project.toml

@@ -94,7 +94,7 @@ PrecompileTools = "1.1"
 Preferences = "1.3"
 Printf = "1"
 RecipesBase = "1.1"
-RecursiveArrayTools = "2.38.10"
+RecursiveArrayTools = "3.27.1"


You are setting this compatibility version so high that SciMLBase will also adopt a higher version, forcing you to use the new version of the DiscreteCallback struct.

That's why even if you're using Julia 1.10 but it still failed here https://github.com/trixi-framework/Trixi.jl/actions/runs/11871869411/job/33084910314?pr=2150#step:7:24

Either change this to a lower version or adopt the new DiscreteCallback struct in precompile.

JoshuaLampert · 2024-12-03T16:08:14Z

Is there any update on this? I also need Trixi.jl being compatible with newer versions of RecursiveArrayTools.jl and OrdinaryDiffEq.jl.

jlchan · 2024-12-03T19:13:22Z

Is there any update on this? I also need Trixi.jl being compatible with newer versions of RecursiveArrayTools.jl and OrdinaryDiffEq.jl.

I've got a conference this week so I won't be able to return to this until next Monday. However, there aren't any major numerical issues to address IMO. The main failing tests are related to allocation (e.g., #2150 (comment)) and precompilation (#2150 (comment)). I believe @huiyuxie has identified some of these issues already.

If you have any time to look at those before I get back to it next week I'd welcome the help.

huiyuxie · 2024-12-04T00:05:36Z

I will open a PR to address this but I'm not familiar with new DiscreteCallback struct I may need help from @ranocha

ranocha · 2024-12-04T05:44:46Z

Please let me know if you have questions

github-actions · 2025-01-28T12:42:05Z

test/test_parabolic_2d.jl

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/test/test_parabolic_2d.jl

Lines 137 to 148 in ac33292

#= FIXME: Values from above (PR) or values from below (main)?

0.0015355076812512347,

0.0033843168272690953,

0.0036531858107444093,

0.009948436427520873

],

linf=[

0.005522560467186466,

0.013425258500736614,

0.013962115643462503,

0.0274831021205042

=#

github-actions · 2025-01-28T12:42:06Z

test/test_parabolic_2d.jl

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/test/test_parabolic_2d.jl

Lines 175 to 186 in ac33292

#= FIXME: Values from above (PR) or values from below (main)?

0.004255101916146402,

0.011118488923215394,

0.011281831283462784,

0.035736564473886

],

linf=[

0.015071710669707805,

0.04103132025860057,

0.03990424085748012,

0.13094017185987106

=#

github-actions · 2025-01-28T12:42:06Z

test/test_parabolic_3d.jl

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/test/test_parabolic_3d.jl

Lines 32 to 45 in ac33292

#= FIXME: Values from above (PR) or values from below (main)?

0.000553284711585015,

0.0006592634909666629,

0.0007776436127373607,

0.0006592634909664286,

0.0038073628897819524

],

linf=[

0.0017039861523628907,

0.002628561703550747,

0.0035310574250866367,

0.002628561703585053,

0.015587829540340437

=#

github-actions · 2025-01-28T12:42:07Z

test/test_parabolic_3d.jl

[JuliaFormatter] _{reported by reviewdog 🐶}

Trixi.jl/test/test_parabolic_3d.jl

Lines 74 to 87 in ac33292

#= FIXME: Values from above (PR) or values from below (main)?

0.001402722725120774,

0.0021322235533272546,

0.002787374144745514,

0.002458747307062751,

0.009978368180191861

],

linf=[

0.0063417504028402405,

0.010306014252245865,

0.015207402509253565,

0.010968264045485343,

0.04745438983155026

=#

huiyuxie · 2025-01-30T07:24:54Z

The helper PR is closed please change as the bot suggested @jlchan

…p_VectorOfArray"

mark all 1D allocated tests as broken for now

…rap_VectorOfArray

jlchan · 2025-01-31T21:44:41Z

This test failure appears to be an issue in OrdinaryDiffEq.jl / DiffEqBase.jl.

jlchan · 2025-01-31T22:11:36Z

This test failure appears to be an issue in OrdinaryDiffEq.jl / DiffEqBase.jl.

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Use `VectorOfArray` in `wrap_array` for `DGMulti` solvers #2150

Use `VectorOfArray` in `wrap_array` for `DGMulti` solvers #2150

jlchan commented Nov 7, 2024 •

edited

Loading

github-actions bot commented Nov 7, 2024

codecov bot commented Nov 7, 2024 •

edited

Loading

jlchan commented Nov 7, 2024 •

edited

Loading

DanielDoehring commented Nov 8, 2024

jlchan commented Nov 8, 2024

DanielDoehring commented Nov 9, 2024

jlchan commented Nov 10, 2024 •

edited

Loading

huiyuxie commented Nov 11, 2024

huiyuxie commented Nov 11, 2024

jlchan commented Nov 11, 2024

huiyuxie commented Nov 11, 2024

huiyuxie commented Nov 15, 2024 •

edited

Loading

jlchan commented Nov 15, 2024

ranocha commented Nov 20, 2024

jlchan commented Nov 20, 2024

jlchan commented Nov 20, 2024

huiyuxie Nov 22, 2024

JoshuaLampert commented Dec 3, 2024

jlchan commented Dec 3, 2024

huiyuxie commented Dec 4, 2024

ranocha commented Dec 4, 2024

github-actions bot Jan 28, 2025

github-actions bot Jan 28, 2025

github-actions bot Jan 28, 2025

github-actions bot Jan 28, 2025

huiyuxie commented Jan 30, 2025

jlchan commented Jan 31, 2025

jlchan commented Jan 31, 2025

	#= FIXME: Values from above (PR) or values from below (main)?
	0.0015355076812512347,
	0.0033843168272690953,
	0.0036531858107444093,
	0.009948436427520873
	],
	linf=[
	0.005522560467186466,
	0.013425258500736614,
	0.013962115643462503,
	0.0274831021205042
	=#

	#= FIXME: Values from above (PR) or values from below (main)?
	0.004255101916146402,
	0.011118488923215394,
	0.011281831283462784,
	0.035736564473886
	],
	linf=[
	0.015071710669707805,
	0.04103132025860057,
	0.03990424085748012,
	0.13094017185987106
	=#

	#= FIXME: Values from above (PR) or values from below (main)?
	0.000553284711585015,
	0.0006592634909666629,
	0.0007776436127373607,
	0.0006592634909664286,
	0.0038073628897819524
	],
	linf=[
	0.0017039861523628907,
	0.002628561703550747,
	0.0035310574250866367,
	0.002628561703585053,
	0.015587829540340437
	=#

	#= FIXME: Values from above (PR) or values from below (main)?
	0.001402722725120774,
	0.0021322235533272546,
	0.002787374144745514,
	0.002458747307062751,
	0.009978368180191861
	],
	linf=[
	0.0063417504028402405,
	0.010306014252245865,
	0.015207402509253565,
	0.010968264045485343,
	0.04745438983155026
	=#

Use VectorOfArray in wrap_array for DGMulti solvers #2150

Are you sure you want to change the base?

Use VectorOfArray in wrap_array for DGMulti solvers #2150

Conversation

jlchan commented Nov 7, 2024 • edited Loading

github-actions bot commented Nov 7, 2024

Review checklist

Purpose and scope

Code quality

Documentation

Testing

Performance

Verification

codecov bot commented Nov 7, 2024 • edited Loading

Codecov Report

jlchan commented Nov 7, 2024 • edited Loading

DanielDoehring commented Nov 8, 2024

jlchan commented Nov 8, 2024

DanielDoehring commented Nov 9, 2024

jlchan commented Nov 10, 2024 • edited Loading

huiyuxie commented Nov 11, 2024

huiyuxie commented Nov 11, 2024

jlchan commented Nov 11, 2024

huiyuxie commented Nov 11, 2024

huiyuxie commented Nov 15, 2024 • edited Loading

jlchan commented Nov 15, 2024

ranocha commented Nov 20, 2024

jlchan commented Nov 20, 2024

jlchan commented Nov 20, 2024

huiyuxie Nov 22, 2024

Choose a reason for hiding this comment

JoshuaLampert commented Dec 3, 2024

jlchan commented Dec 3, 2024

huiyuxie commented Dec 4, 2024

ranocha commented Dec 4, 2024

github-actions bot Jan 28, 2025

Choose a reason for hiding this comment

github-actions bot Jan 28, 2025

Choose a reason for hiding this comment

github-actions bot Jan 28, 2025

Choose a reason for hiding this comment

github-actions bot Jan 28, 2025

Choose a reason for hiding this comment

huiyuxie commented Jan 30, 2025

jlchan commented Jan 31, 2025

jlchan commented Jan 31, 2025

Use `VectorOfArray` in `wrap_array` for `DGMulti` solvers #2150

Use `VectorOfArray` in `wrap_array` for `DGMulti` solvers #2150

jlchan commented Nov 7, 2024 •

edited

Loading

codecov bot commented Nov 7, 2024 •

edited

Loading

jlchan commented Nov 7, 2024 •

edited

Loading

jlchan commented Nov 10, 2024 •

edited

Loading

huiyuxie commented Nov 15, 2024 •

edited

Loading