Closed mrhardman closed 9 months ago
I've modified the calculus.jl setup back to the version in master
. I'm working on merging the current master
- just chasing down some small bug that I seem to have added doing that.
I noticed there are some new dependencies being added which don't actually seem to be used (any more): ArbNumerics, Cubature, Elliptic, IterativeSolvers. @mrhardman is it OK to remove these?
I've modified the calculus.jl setup back to the version in master. I'm working on merging the current master - just chasing down some small bug that I seem to have added doing that.
This will be in a PR to this PR?
I noticed there are some new dependencies being added which don't actually seem to be used (any more): ArbNumerics, Cubature, Elliptic, IterativeSolvers. @mrhardman is it OK to remove these?
Yes, I forgot to remove these. Please remove them.
This will be in a PR to this PR?
What do you prefer? I'm happy to do it as a PR, or can just push the changes directly.
What do you prefer? I'm happy to do it as a PR, or can just push the changes directly.
Let's go with a PR so that I can see what has to change.
We should also sit down to discuss what is needed for https://github.com/mabarnes/moment_kinetics/issues/160.
@mrhardman could you try the test/fokker_planck_tests.jl
in parallel? When I run on 2 cores, I get test failures from the self collisions case without numerical conserving terms
Fokker Planck tests
test_self_operator=true test_numerical_conserving_terms=false test_parallelism = true test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic_solv
e_for_d2Gdvperp2=false: Test Failed at /home/john/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:367
Expression: isapprox(C_M_L2, rtol_L2; atol = atol_L2)
Evaluated: isapprox(1.4224224829260572e-5, 7.0e-6; atol = 7.0e-6)
Stacktrace:
[1] macro expansion
@ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:478 [inlined]
[2] macro expansion
@ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:367 [inlined]
[3] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/src/looping.jl:463 [inlined]
[4] macro expansion
@ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:357 [inlined]
[5] macro expansion
@ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1586 [inlined]
[6] macro expansion
@ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:303 [inlined]
[7] macro expansion
@ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined]
[8] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:293 [inlined]
[9] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined]
[10] runtests() @ Main.FokkerPlanckTests ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:72
test_self_operator=true test_numerical_conserving_terms=false test_parallelism = true test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic_solve_for_d2Gdvperp2=false: Test Failed at /home/john/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:387
Expression: isapprox(delta_n, rtol; atol = atol) Evaluated: isapprox(3.869261784452312e-6, 0.0; atol = 1.0e-12)
Stacktrace:
[1] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:478 [inlined]
[2] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:387 [inlined]
[3] macro expansion
@ ~/physics/moment_kinetics-merge_fkpl_collisions/src/looping.jl:463 [inlined]
[4] macro expansion
@ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:357 [inlined]
[5] macro expansion
@ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1586 [inlined]
[6] macro expansion
@ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:303 [inlined]
[7] macro expansion
@ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined]
[8] macro expansion
@ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:293 [inlined]
[9] macro expansion
@ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined]
[10] runtests()
@ Main.FokkerPlanckTests ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:72
test_self_operator=true test_numerical_conserving_terms=false test_parallelism = true test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic_solv
e_for_d2Gdvperp2=false: Test Failed at /home/john/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:389
Expression: isapprox(delta_upar, rtol; atol = atol)
Evaluated: isapprox(2.1175349685924036e-8, 0.0; atol = 1.0e-9)
Stacktrace:
[1] macro expansion
@ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:478 [inlined]
[2] macro expansion
@ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:389 [inlined]
[3] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/src/looping.jl:463 [inlined]
[4] macro expansion
@ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:357 [inlined]
[5] macro expansion
@ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1586 [inlined]
[6] macro expansion
@ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:303 [inlined]
[7] macro expansion
@ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined]
[8] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:293 [inlined]
[9] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined]
[10] runtests() @ Main.FokkerPlanckTests ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:72
test_self_operator=true test_numerical_conserving_terms=false test_parallelism = true test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic_solve_for_d2Gdvperp2=false: Test Failed at /home/john/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:395
Expression: isapprox(delta_pressure, rtol; atol = atol) Evaluated: isapprox(1.3462986951440539e-5, 0.0; atol = 1.0e-8)
Stacktrace:
[1] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:478 [inlined]
[2] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:395 [inlined]
[3] macro expansion
@ ~/physics/moment_kinetics-merge_fkpl_collisions/src/looping.jl:463 [inlined]
[4] macro expansion
@ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:357 [inlined]
[5] macro expansion
@ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1586 [inlined]
[6] macro expansion
@ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:303 [inlined]
[7] macro expansion
@ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined]
[8] macro expansion
@ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:293 [inlined]
[9] macro expansion
@ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined]
[10] runtests()
@ Main.FokkerPlanckTests ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:72
Test Summary:
| Pass Fail Total Time
Fokker Planck tests
| 71 4 75 25.5s
- test weak-form 2D differentiation
| 5 5 7.7s
- test weak-form Rosenbluth potential calculation
| 24 24 7.5s
- test weak-form collision operator calculation
| 42 4 46 10.2s
test_self_operator=true test_numerical_conserving_terms=false test_parallelism = false test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic
_solve_for_d2Gdvperp2=false | 6 6 1.0s test_self_operator=false test_numerical_conserving_terms=false test_parallelism = false test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebrai
c_solve_for_d2Gdvperp2=false | 4 4 0.1s
test_self_operator=true test_numerical_conserving_terms=true test_parallelism = false test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic_
solve_for_d2Gdvperp2=false | 6 6 0.2s
test_self_operator=true test_numerical_conserving_terms=false test_parallelism = true test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic_
solve_for_d2Gdvperp2=false | 2 4 6 1.8s
test_self_operator=true test_numerical_conserving_terms=false test_parallelism = false test_dense_construction=true use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic_
solve_for_d2Gdvperp2=false | 6 6 0.1s
test_self_operator=true test_numerical_conserving_terms=false test_parallelism = false test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=true use_Maxwellian_field_particle_distribution=false algebraic_
solve_for_d2Gdvperp2=false | 6 6 0.1s test_self_operator=true test_numerical_conserving_terms=false test_parallelism = false test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=true algebraic_
solve_for_d2Gdvperp2=false | 6 6 0.1s test_self_operator=true test_numerical_conserving_terms=false test_parallelism = false test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic
_solve_for_d2Gdvperp2=true | 6 6 0.1s ERROR: LoadError: Some tests did not pass: 71 passed, 4 failed, 0 errored, 0 broken.
in expression starting at ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:445
Regarding this point
Having restored the functionality of "fkpl_functional_test.jl" for the strong form operator, I get the following run-time error
ERROR: UndefVarError: boundary_distributions
not defined
This is because euler_time_advance no longer has any structs used for imposing boundary conditions. These structs were necessary for the conserving terms in my original implementation.
@johnomotani We should choose now whether or not to upgrade this functionality (it should be possible using the same method as in the weak-form version, but this would require some time to test) or remove the code altogether from the time advance.
@mrhardman could you try the
test/fokker_planck_tests.jl
in parallel? When I run on 2 cores, I get test failures from the self collisions case without numerical conserving termsFokker Planck tests test_self_operator=true test_numerical_conserving_terms=false test_parallelism = true test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic_solv e_for_d2Gdvperp2=false: Test Failed at /home/john/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:367 Expression: isapprox(C_M_L2, rtol_L2; atol = atol_L2) Evaluated: isapprox(1.4224224829260572e-5, 7.0e-6; atol = 7.0e-6) Stacktrace: [1] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:478 [inlined] [2] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:367 [inlined] [3] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/src/looping.jl:463 [inlined] [4] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:357 [inlined] [5] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1586 [inlined] [6] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:303 [inlined] [7] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined] [8] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:293 [inlined] [9] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined] [10] runtests() @ Main.FokkerPlanckTests ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:72 test_self_operator=true test_numerical_conserving_terms=false test_parallelism = true test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic_solve_for_d2Gdvperp2=false: Test Failed at /home/john/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:387 Expression: isapprox(delta_n, rtol; atol = atol) Evaluated: isapprox(3.869261784452312e-6, 0.0; atol = 1.0e-12) Stacktrace: [1] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:478 [inlined] [2] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:387 [inlined] [3] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/src/looping.jl:463 [inlined] [4] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:357 [inlined] [5] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1586 [inlined] [6] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:303 [inlined] [7] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined] [8] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:293 [inlined] [9] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined] [10] runtests() @ Main.FokkerPlanckTests ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:72 test_self_operator=true test_numerical_conserving_terms=false test_parallelism = true test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic_solv e_for_d2Gdvperp2=false: Test Failed at /home/john/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:389 Expression: isapprox(delta_upar, rtol; atol = atol) Evaluated: isapprox(2.1175349685924036e-8, 0.0; atol = 1.0e-9) Stacktrace: [1] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:478 [inlined] [2] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:389 [inlined] [3] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/src/looping.jl:463 [inlined] [4] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:357 [inlined] [5] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1586 [inlined] [6] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:303 [inlined] [7] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined] [8] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:293 [inlined] [9] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined] [10] runtests() @ Main.FokkerPlanckTests ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:72 test_self_operator=true test_numerical_conserving_terms=false test_parallelism = true test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic_solve_for_d2Gdvperp2=false: Test Failed at /home/john/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:395 Expression: isapprox(delta_pressure, rtol; atol = atol) Evaluated: isapprox(1.3462986951440539e-5, 0.0; atol = 1.0e-8) Stacktrace: [1] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:478 [inlined] [2] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:395 [inlined] [3] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/src/looping.jl:463 [inlined] [4] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:357 [inlined] [5] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1586 [inlined] [6] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:303 [inlined] [7] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined] [8] macro expansion @ ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:293 [inlined] [9] macro expansion @ ~/pkg/julia-1.9.3/share/julia/stdlib/v1.9/Test/src/Test.jl:1498 [inlined] [10] runtests() @ Main.FokkerPlanckTests ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:72 Test Summary: | Pass Fail Total Time Fokker Planck tests | 71 4 75 25.5s - test weak-form 2D differentiation | 5 5 7.7s - test weak-form Rosenbluth potential calculation | 24 24 7.5s - test weak-form collision operator calculation | 42 4 46 10.2s test_self_operator=true test_numerical_conserving_terms=false test_parallelism = false test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic _solve_for_d2Gdvperp2=false | 6 6 1.0s test_self_operator=false test_numerical_conserving_terms=false test_parallelism = false test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebrai c_solve_for_d2Gdvperp2=false | 4 4 0.1s test_self_operator=true test_numerical_conserving_terms=true test_parallelism = false test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic_ solve_for_d2Gdvperp2=false | 6 6 0.2s test_self_operator=true test_numerical_conserving_terms=false test_parallelism = true test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic_ solve_for_d2Gdvperp2=false | 2 4 6 1.8s test_self_operator=true test_numerical_conserving_terms=false test_parallelism = false test_dense_construction=true use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic_ solve_for_d2Gdvperp2=false | 6 6 0.1s test_self_operator=true test_numerical_conserving_terms=false test_parallelism = false test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=true use_Maxwellian_field_particle_distribution=false algebraic_ solve_for_d2Gdvperp2=false | 6 6 0.1s test_self_operator=true test_numerical_conserving_terms=false test_parallelism = false test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=true algebraic_ solve_for_d2Gdvperp2=false | 6 6 0.1s test_self_operator=true test_numerical_conserving_terms=false test_parallelism = false test_dense_construction=false use_Maxwellian_Rosenbluth_coefficients=false use_Maxwellian_field_particle_distribution=false algebraic _solve_for_d2Gdvperp2=true | 6 6 0.1s ERROR: LoadError: Some tests did not pass: 71 passed, 4 failed, 0 errored, 0 broken. in expression starting at ~/physics/moment_kinetics-merge_fkpl_collisions/test/fokker_planck_tests.jl:445
On commit e1023ce6536bb1c0badc6c7b54fbe57734a6ebaa my tests pass in parallel with two cores $ mpirun -n 2 julia --project -e "using Pkg; Pkg.test()" Test Summary: | Pass Broken Total Time moment_kinetics tests | 117871 2 117873 13m14.4s Testing moment_kinetics tests passed
Bug seems to have been fixed by the commits between 6f5193d and e1023ce. Don't know what the problem was, but when I merged your latest few commits in, it went away!
Bug seems to have been fixed
...or maybe not. I thought I tested and it went away, but now running test/fokker_planck_tests.jl
on 2 cores at 6045a81 gives me huge errors for the same test.
...or maybe not. I thought I tested and it went away, but now running
test/fokker_planck_tests.jl
on 2 cores at 6045a81 gives me huge errors for the same test.
Do you mean that the test fails for you in the branch of PR #149 or in the local merged version of PR #162 ? Recall that I did have a shared-memory bug which should have been fixed with https://github.com/mabarnes/moment_kinetics/pull/149/commits/c4ef1e6d242e1062d517befeea732a7ea7e23ee5.
On my local machine the tests still pass in parallel in the latest commit https://github.com/mabarnes/moment_kinetics/pull/149/commits/6045a81c648a06b674faba11bb77a5da4bb95ba2 to this PR. (#149)
The failure happened on both #149 and #162. It was a bug - fix coming shortly, I'm just testing some fixes to the debugging functions before I push.
The failure happened on both #149 and #162. It was a bug - fix coming shortly, I'm just testing some fixes to the debugging functions before I push.
I haven't seen any further bugs in #149 on the two machines that I am testing, even up to 64 cores. I am surprised!
The bug only affects the use_Maxwellian_Rosenbluth_coefficients=true
or use_Maxwellian_field_particle_distribution=true
branches, see c6b081b, so probably wouldn't affect most of your testing?
The bug only affects the
use_Maxwellian_Rosenbluth_coefficients=true
oruse_Maxwellian_field_particle_distribution=true
branches, see c6b081b, so probably wouldn't affect most of your testing?
Do you see the bug in the automatic tests? If yes, then I should have seen the same bug. (EDIT: Have now read the code and agreed with your changes!)
time_evolving_new.zip time_evolving_1D2V_new.zip Uploads including example input files and output .gif/.pdfs for a 1D2V sheath problem with C[F,F] and wall boundaries, and the corresponding Maxwellian relaxation problem (0D2V). These simulations suggest that commit c3b6261115ba671cf338e88bbb762de121f59c82 behaves acceptably for 1) small enough time steps and 2) simulations where the stagnation point z=0 are eliminated with the choice of grids.
Determine how to usefully store interactive tests scripts "2D_FEM_assembly_test.jl" "fkpl_functional_test.jl" and "GaussLobattoLegendre_test.jl".
Maybe we could move these and similar scripts to a test_scripts
subdirectory?
- Clean up the earlier "strong-form" implementation of the collision operator.
@johnomotani We should choose now whether or not to upgrade this functionality (it should be possible using the same method as in the weak-form version, but this would require some time to test) or remove the code altogether from the time advance.
I'd lean towards removing it rather than invest more time in something we aren't necessarily going to use?
I'd lean towards removing it rather than invest more time in something we aren't necessarily going to use?
Sounds good to me. I think the only feature that we need to save is the ability to use the integration weights to find the potentials everywhere in vpa,vperp. That should be fairly simple to do, and we can make sure the functionality is tested even if it isn't used in the main time advance. Another benefit would be that we can remove 4 or 5 pdf sized buffer arrays from scratch_dummy, which should help with memory usage. I'll take a look at this ASAP.
Determine how to usefully store interactive tests scripts "2D_FEM_assembly_test.jl" "fkpl_functional_test.jl" and "GaussLobattoLegendre_test.jl".
Maybe we could move these and similar scripts to a
test_scripts
subdirectory?
That sounds like a good solution. "fkpl_functional_test.jl" will become redundant if we remove the older version of the collision operator.
I found the bugs in the MMS test of the Krook collision operator. Fixed in 59c2a8f.
I just made the feature and test for the direct integration method for calculating the Rosenbluth potentials. In doing so I accidentally ressurected the branch "merge_fkpl_collisions-merge-master". My next step will be to review the comments in fokker_planck.jl whilst removing the "strong form" collision operator and associated flags.
I also note that a nice feature would be one where nu_ii is computed using reference parameters like for the Krook operator. That should probably go in a separate PR, but I just note the idea here.
The only remaining outstanding point is to test the manufactured solutions test. I am working on this now. Since @johnomotani found the bug in the Krook operator, I hope that this will not take long.
I have tested
runs/2D-wall_MMS_nel_r_8_z_8_vpa_16_vperp_1_diss.toml runs/1D-wall_MMS_new_nel_r_1_z_16_vpa_16_vperp_1_krook.toml
with good results showing small errors.
Testing
runs/2D-sound-wave_cheb_ion_only_nel_r_8_z_8_vpa_8_vperp_8.toml
I find that ppar and vth have O(1) errors.
julia> analyze_and_plot_data("runs/2D-sound-wave_cheb_ion_only_nel_r_8_z_8_vpa_8_vperp_8/") time/ (Lref/cref): [0.0, 0.32000000000000023] test: phi: [1.737592677941611e-8, 5.801355406684373e-9] test: Er: [3.575165147234937e-10, 1.0235171027365706e-8] test: Ez: [3.575165147234937e-10, 4.2536024303175266e-8] test: dens: ion [2.6120856384952635e-8, 8.793044046291286e-9] test: upar: ion [6.448856152718709e-18, 8.168700163528552e-9] test: ppar: ion [0.7516393356062538, 0.7516393734155177] test: pperp: ion [1.3440885806176344e-8, 4.154681744568655e-9] test: vthi: ion [0.18350344427438467, 0.1835034489276695] test: pdf: ion [7.760007061812827e-19, 9.390632736607075e-8]
Testing
runs/2D-sound-wave_cheb_ion_only_nel_r_8_z_8_vpa_8_vperp_8.toml
I find that ppar and vth have O(1) errors.
This is also true on master
. Factor of 2 from normalisation somewhere maybe?
This is also true on
master
. Factor of 2 from normalisation somewhere maybe?
I think this is the most likely reason. Will take a quick look.
It's better than a bug:
It's a feature that I didn't support yet!
Using the input .toml here https://github.com/mabarnes/moment_kinetics/blob/master/runs/2D-sound-wave_cheb_ion_only_nel_r_4_z_4_vpa_4_vperp_4.toml, with 16 cores and z_nelement_local = 1 and r_nelement_local = 1, and changing vperp_ngrid => 1 and vperp_nelement => 1, I receive a segmentation fault. EDIT: this fault doesn't seem to be reproducible...
I cannot reproduce the segmentation fault with distributed memory-mpi and periodic wall boundary conditions. If this is a bug in our source, it shouldn't have any connection to the collision operator. I would now be happy for the merge to go ahead. Please leave reviews @johnomotani @mabarnes.
@johnomotani In your next review it would be useful to consider whether and how we plan to address issue https://github.com/mabarnes/moment_kinetics/issues/140. The choice of parallelisation also affects fokker_planck.jl, the testing modules, and much of fokker_planck_calculus.jl (apart from the functions associated with initialisation of the integration weights). My brief investigations suggest that for small velocity grids parallelisation can give at best a factor of 2 speedup, but this factor may be larger for larger problem sizes. Removing the ability to parallelise over vpa and vperp entirely does not seem that attractive from the point of view of wishing to keep the ability to profile the operator.
The last commit https://github.com/mabarnes/moment_kinetics/pull/149/commits/0fda1f6902726a09578ac75c4e3fffd9f761cce1 partially addresses issue https://github.com/mabarnes/moment_kinetics/issues/157. All that is required is to go and change the appropriate if statements in the boundary conditions function.
I think I'm done reviewing the code now - once we sort out the discussion on setup_global_weak_form_matrix!()
, I think it's ready to merge.
Do we have any ideas on why the MacOS checks are failing consistently?
I have noticed that in initial_conditions.jl we are now using the centred version of reconcile_element_boundaries_MPI! rather than the upwind version. This might explain the bugs that I am seeing with spikes at element boundaries on distriubted memory -- commits to follow. We should make sure distributed memory is tested to prevent bugs like this creeping in..
https://github.com/mabarnes/moment_kinetics/blob/ee12e1a669375adf512f3e591697bd2566aeb606/src/initial_conditions.jl#L1197 https://github.com/mabarnes/moment_kinetics/blob/ee12e1a669375adf512f3e591697bd2566aeb606/src/initial_conditions.jl#L1247 https://github.com/mabarnes/moment_kinetics/blob/ee12e1a669375adf512f3e591697bd2566aeb606/src/initial_conditions.jl#L1379 https://github.com/mabarnes/moment_kinetics/blob/ee12e1a669375adf512f3e591697bd2566aeb606/src/initial_conditions.jl#L1430
Do we have any ideas on why the MacOS checks are failing consistently?
Don't really know. Sometimes seems to just get timed out for being slow, but sometimes they segfault or seem to hang. My current hope is that when we move the postprocessing to separate packages, getting rid of some dependencies (especially Python) might help.
I have noticed that in initial_conditions.jl we are now using the centred version of reconcile_element_boundaries_MPI! rather than the upwind version. This might explain the bugs that I am seeing with spikes at element boundaries on distriubted memory -- commits to follow. We should make sure distributed memory is tested to prevent bugs like this creeping in..
What's changed? I can't see that we ever used a different version of reconcile_element_boundaries_MPI!()
in these functions that enforce boundary conditions. I also wouldn't expect to - the element boundary points at MPI boundaries should just have the same value on both processes. There is no 'advection speed' passed to these functions to use for upwinding.
@johnomotani It would be handy for you to take a look at the last commit for your comments. It doesn't fix my current bug issue #169 but I think it is an important correction.
I have noticed that in initial_conditions.jl we are now using the centred version of reconcile_element_boundaries_MPI! rather than the upwind version. This might explain the bugs that I am seeing with spikes at element boundaries on distriubted memory -- commits to follow. We should make sure distributed memory is tested to prevent bugs like this creeping in..
What's changed? I can't see that we ever used a different version of
reconcile_element_boundaries_MPI!()
in these functions that enforce boundary conditions. I also wouldn't expect to - the element boundary points at MPI boundaries should just have the same value on both processes. There is no 'advection speed' passed to these functions to use for upwinding.
Perhaps I am getting confused. In any case, when we reconcile the boundaries here, I would have thought that we would want to take as the "correct" value of the pdf to be the upwind one. Do you disagree?
Perhaps I am getting confused. In any case, when we reconcile the boundaries here, I would have thought that we would want to take as the "correct" value of the pdf to be the upwind one. Do you disagree?
I think it shouldn't make any difference. Communication was already taken care of in the derivative routines, so the two values should be identical. The communication being done in the boundary conditions functions is, I think, only to be doubly sure that some rounding error difference does not cause the values to diverge.
If it did make any difference, I think your update would not be quite correct, as the boundary conditions are applied at the end of a time-step stage and after the RK update is done, so the distribution function, etc. are the updated values, but the adv_fac
was calculated with the old values, so it doesn't seem consistent to use it here.
So I think you'll find that 816e058 makes not difference at all, and we should revert it. If it does make any difference then I am wrong and need to look much harder at the communication in those boundary condition functions.
Perhaps I am getting confused. In any case, when we reconcile the boundaries here, I would have thought that we would want to take as the "correct" value of the pdf to be the upwind one. Do you disagree?
I think it shouldn't make any difference. Communication was already taken care of in the derivative routines, so the two values should be identical. The communication being done in the boundary conditions functions is, I think, only to be doubly sure that some rounding error difference does not cause the values to diverge. If it did make any difference, I think your update would not be quite correct, as the boundary conditions are applied at the end of a time-step stage and after the RK update is done, so the distribution function, etc. are the updated values, but the
adv_fac
was calculated with the old values, so it doesn't seem consistent to use it here.So I think you'll find that 816e058 makes not difference at all, and we should revert it. If it does make any difference then I am wrong and need to look much harder at the communication in those boundary condition functions.
OK. I am open to doing this. There was a minor bug that I did catch there where one of the advect structs was used incorrectly. I should have fixed it separately but it was tempting to commit it all at once.
Ah, sorry - l.1016-1017 in gauss_legendre.jl https://github.com/mabarnes/moment_kinetics/blob/816e058a125b773df1af649a4a7145570d737998/src/gauss_legendre.jl#L1016-L1017
Those boundary terms end up being something like the difference of the first derivative evaluated in the elements on either side of the boundary, once you add the contributions from the two elements that share the boundary point, right? So for smooth input data they converge to zero as the resolution increases? That probably makes them a bit hard to test. Maybe you could cook something up like calculate the second derivative of abs(sin(x))
? I don't know if it'd be possible to calculate analytically the expected result for that though...
Ah, sorry - l.1016-1017 in gauss_legendre.jl
I thought I had tested these terms in https://github.com/mabarnes/moment_kinetics/blob/merge_fkpl_collisions/test_scripts/GaussLobattoLegendre_test.jl. The signs look correct to me. What we want is
int phi d^2 f / dx^2 d x = [ phi d f / d x ]^{xmax}{x_min} - int (dphi /d x) (d f / d x ) d x.
Does that make sense?
Those boundary terms end up being something like the difference of the first derivative evaluated in the elements on either side of the boundary, once you add the contributions from the two elements that share the boundary point, right? So for smooth input data they converge to zero as the resolution increases? That probably makes them a bit hard to test. Maybe you could cook something up like calculate the second derivative of
abs(sin(x))
? I don't know if it'd be possible to calculate analytically the expected result for that though...
Initially I would remove the internal boundary points and check that the result was the same.
I agree, flipping those signs seems to be wrong - I tried calculating a second derivative of abs(sin(2*pi*vpa.grid/vpa.L))
, and the version with flipped signs gave large errors everywhere.
It does look like a problem at element boundaries though - in your simulation that 'goes wrong', one symptom of the problem seems to be that there is a peak in f at vpa=0, but the vpa-diffusion evaluates the second derivative as positive at vpa=0 and so makes the peak grow, which seems clearly wrong. I don't understand why it does that though!
Using the updated test script in the REPL, we can confirm the correctness for a single derivative.
julia> include("test_scripts/GaussLobattoLegendre_test.jl")
julia> gausslegendre_test(ngrid=17,nelement=16, L_in=6.0)
vpa test
F_err: 2.220446049250313e-16 F_exact: 1.7724538509055159 F_num: 1.772453850905516 max(df_err) (interpolation): 8.770761894538737e-15 max(d2f_err) (double first derivative by interpolation): 6.433603649824704e-13 max(df_err) (weak form): 6.140921104957897e-14 max(d2f_err) (weak form): 2.6008240006092365e-10
vperp test
F_err: 3.3306690738754696e-16 F_exact: 1.0 F_num: 0.9999999999999997 max(df_err) (interpolation): 1.2101430968414206e-14 max(d2f_err) (double first derivative by interpolation): 1.3371082019375535e-11 max(divg_err) (weak form): 8.908651594197181e-12 max(d2f_err) (weak form): 2.0351857976663723e-9 max(laph_err) (weak form): 1.0388290228036112e-9 max(divg_err) (interpolation): 2.1036505870597466e-12 max(laph_err) (interpolation): 4.0971048775872987e-10
However, I did turn off the explicit boundary terms in the collision operator already because they caused bad behaviour in the elliptic solvers. It might be that the terms should never be used because of a numerical instability.
I think turning off the terms at interior element boundaries might be the way to go. Eliminating them at (non-periodic) global boundaries would give large errors, but eliminating them at internal boundaries I think is equivalent to assuming that the first derivative is continuous between elements. For smooth functions this makes no noticeable difference, but if there is a cusp on an element boundary, the version without 'explicit boundary terms' at interior boundaries evaluates in a way that would decrease the size of the cusp (although it also introduces large 'errors' at nearby points when the cusp is large, presumably due to Gibbs ringing) - I checked this with the abs(sin(x)) test. I think that would give better numerical stability, because if there is a kink at an element boundary, the version without 'explicit bonudary terms' at internal element boundaries will smooth it.
I think turning off the terms at interior element boundaries might be the way to go. Eliminating them at (non-periodic) global boundaries would give large errors, but eliminating them at internal boundaries I think is equivalent to assuming that the first derivative is continuous between elements. For smooth functions this makes no noticeable difference, but if there is a cusp on an element boundary, the version without 'explicit boundary terms' at interior boundaries evaluates in a way that would decrease the size of the cusp (although it also introduces large 'errors' at nearby points when the cusp is large, presumably due to Gibbs ringing) - I checked this with the abs(sin(x)) test. I think that would give better numerical stability, because if there is a kink at an element boundary, the version without 'explicit bonudary terms' at internal element boundaries will smooth it.
I have carried this out successfully (according to checks so far) in https://github.com/mabarnes/moment_kinetics/pull/149/commits/d306695e713bbd35c5dadc5c3fd7ecc194de7c66. There is an example input file in https://github.com/mabarnes/moment_kinetics/pull/149/commits/6dfd0d9e310ab9de8fc9e385d9e3dc325aeda934 that could form the basis for a check in test. The z and r dissipations could be tested in separate 1D tests to give coverage to these features.
See update to https://github.com/mabarnes/moment_kinetics/issues/169.
Perhaps I am getting confused. In any case, when we reconcile the boundaries here, I would have thought that we would want to take as the "correct" value of the pdf to be the upwind one. Do you disagree?
I think it shouldn't make any difference. Communication was already taken care of in the derivative routines, so the two values should be identical. The communication being done in the boundary conditions functions is, I think, only to be doubly sure that some rounding error difference does not cause the values to diverge. If it did make any difference, I think your update would not be quite correct, as the boundary conditions are applied at the end of a time-step stage and after the RK update is done, so the distribution function, etc. are the updated values, but the
adv_fac
was calculated with the old values, so it doesn't seem consistent to use it here.So I think you'll find that 816e058 makes not difference at all, and we should revert it. If it does make any difference then I am wrong and need to look much harder at the communication in those boundary condition functions.
Reverting my commit affecting the initial_conditions.jl module with https://github.com/mabarnes/moment_kinetics/pull/149/commits/ca15efd7e2e1bc74999f7e1dcc6947cfe0dd808c does not fix issue #169. Another change which affects distributed memory MPI must be responsible.
Pull request containing the weak-form Fokker-Planck self-collision operator for ions and associated features (check-in tests, underlying ability to support cross-ion-species collisions). Creating this request now to gather feedback before the merge. The major additions are the fokker_planck*.jl modules and the gauss_legendre.jl module, which provides the weak-form matrices necessary for the implementation.
To do list (before merging):
To do list (after merging):