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.

Sign up for GitHub

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

Start writing automated tests for MMS; renormalise manufactured distribution functions #64

Draft
wants to merge 30 commits into
base: radial-vperp-standard-DKE-with-neutrals
Choose a base branch
from
Draft

Start writing automated tests for MMS; renormalise manufactured distribution functions #64

Changes from 1 commit
Commits
Show all changes
30 commits
Select commit Hold shift + click to select a range
f10f8fd
Refactor loading of coordinates data
johnomotani May 26, 2022
ac40519
Start working on automated test for MMS
johnomotani May 25, 2022
4de5828
Use absolute not relative errors
johnomotani May 30, 2022
895a33b
Fix missing mk_float argument to ensure Chebyshev weights are mk_float
johnomotani May 31, 2022
03b20ea
Use cos instead of sin so manufactured solution has perturbation at t=0
johnomotani Jun 1, 2022
2efe60a
Functions to build manufactured RHS expression/function
johnomotani Jun 1, 2022
1803b44
Bugfix - use begin_s_r_z_region() in ssp_rk!()
johnomotani Jun 1, 2022
9c49fd5
Hacky function for evaluating df/dt
johnomotani Jun 1, 2022
1e8a10c
Move advance_info definition to input_structs.jl
johnomotani Jun 4, 2022
a4ef258
Add vpa variation in manufactured solution for dfni
johnomotani Jun 4, 2022
3ea98bd
Test case evaluating RHS with manufactured f
johnomotani Jun 2, 2022
4d3d192
Parallelise loop initialising manufactured solutions
johnomotani Jun 6, 2022
1ea92f2
Fix so "zero" bc works for vpa coordinate
johnomotani Jun 6, 2022
1fd2e9d
MMS test: use "zero" bc for vpa, non-zero upar, correct f normalization
johnomotani Jun 6, 2022
2e3e252
Remove points affected by bc in MMS test
johnomotani Jun 7, 2022
27c7d42
Support neutrals when calculating rhs in manufactured_solns.jl
johnomotani Jun 20, 2022
5693c4e
Support neutrals in test/manufactured_ddt_tests.jl
johnomotani Jun 20, 2022
859c3c4
Make initial neutral density perturbation non-zero
johnomotani Jun 20, 2022
49bfb9a
Fix dfnn normalisation, add odd components, rationalise argument orders
johnomotani Jun 20, 2022
b4fa510
Fix cartesian_dfni_sym(), gyroaveraged_dfnn_sym() to include norm
johnomotani Jun 20, 2022
8d166a8
(Re-)add charge exchange term in manufactured rhs for ion equation
johnomotani Jun 20, 2022
db2cf67
Fix indexing typo in charge_exchange_collisions_3V!()
johnomotani Jun 21, 2022
ee465d7
Simplify declaration of netcdf_info struct
johnomotani Jun 22, 2022
f404095
Update load_coordinate_data() for refactored coordinates setup
johnomotani Jun 22, 2022
7023930
Update sound wave tests for refactored coordinates setup
johnomotani Jun 22, 2022
2fae0b8
Fix some accidental allocations
johnomotani Jun 23, 2022
0917b14
Run cases without neutrals in manufactured_ddt_tests.jl
johnomotani Jun 23, 2022
d3b3d67
Use slightly lower v-space resolution for MMS tests
johnomotani Jun 23, 2022
15a3b70
Parallelise manufactured_solutions_as_arrays(), manufactured_rhs_as_a…
johnomotani Jun 23, 2022
3f63c75
More convenient choice of maximum n_element in manufactured_ddt_tests.jl
johnomotani Jun 23, 2022
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
27 changes: 13 additions & 14 deletions src/manufactured_solns.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -102,12 +102,13 @@ using ..type_definitions
end
return dfnn
end
function gyroaveraged_dfnn_sym(Lr,Lz,r_bc,z_bc,geometry,composition)
densn = densn_sym(Lr,Lz,r_bc,z_bc,geometry,composition)
#if (r_bc == "periodic" && z_bc == "periodic")
dfnn = densn * exp( - vpa^2 - vperp^2 )
#end
return dfnn
function gyroaveraged_dfnn_sym(Lr,Lz,Lvzeta,Lvr,Lvz,r_bc,z_bc,geometry,composition)
dfnn = dfnn_sym(Lr,Lz,Lvzeta,Lvr,Lvz,r_bc,z_bc,geometry,composition)
# Because of particular form of odd component in dfnn_sym, if we set vzeta=-vr
# then we eliminate the odd component in the perpendicular velocities,
# effectively gyro-averaging
gav_dfnn = substitute(dfnn, Dict(vzeta=>vperp/sqrt(two), vr=>-vperp/sqrt(two), vz=>vpa))
Copy link
Collaborator

Choose a reason for hiding this comment

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

I am worried that this operation is only correct when the field line pitch is inifinite (or zero), so that z and b are codirectional. I think that you are wanting vzeta and vr to be the velocities in the plane perpendicular to the magnetic field line, which they are not. In my original 'gyroaveraged' functions code, I was exploiting that vperp^2 + vpa^2 = vtot^2 = vz^2 + vr^2 + vzeta^2, so the gyroaverage was trivial. Am I missing the point here?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

I made a mistake here - was thinking z and b are in the same direction. Was trying to be clever and give non-zero flow, but this doesn't work 😞

return gav_dfnn
end

# Define some constants in mk_float. Avoids loss of precision due to implicit
Expand Down Expand Up @@ -212,12 +213,10 @@ using ..type_definitions
end
return dfni
end
function cartesian_dfni_sym(Lr,Lz,r_bc,z_bc)
densi = densi_sym(Lr,Lz,r_bc,z_bc)
#if (r_bc == "periodic" && z_bc == "periodic") || (r_bc == "Dirichlet" && z_bc == "periodic")
dfni = densi * exp( - vz^2 - vr^2 - vzeta^2)
#end
return dfni
function cartesian_dfni_sym(Lr,Lz,Lvperp,Lvpa,r_bc,z_bc,geometry,nr)
dfni = dfni_sym(Lr,Lz,Lvperp,Lvpa,r_bc,z_bc,geometry,nr)
vrvzvzeta_dfni = substitute(dfni, Dict(vperp=>sqrt(vzeta^2+vr^2), vpa=>vz))
return vrvzvzeta_dfni
end

function manufactured_solutions(Lr,Lz,Lvperp,Lvpa,Lvzeta,Lvr,Lvz,r_bc,z_bc,geometry,composition,nr)
Expand Down Expand Up @@ -259,12 +258,12 @@ using ..type_definitions
# ion manufactured solutions
densi = densi_sym(Lr,Lz,r_bc,z_bc)
dfni = dfni_sym(Lr,Lz,Lvperp,Lvpa,r_bc,z_bc,geometry,nr)
vrvzvzeta_dfni = cartesian_dfni_sym(Lr,Lz,r_bc,z_bc) #dfni in vr vz vzeta coordinates
vrvzvzeta_dfni = cartesian_dfni_sym(Lr,Lz,Lvperp,Lvpa,r_bc,z_bc,geometry,nr) #dfni in vr vz vzeta coordinates

# neutral manufactured solutions
densn = densn_sym(Lr,Lz,r_bc,z_bc,geometry,composition)
dfnn = dfnn_sym(Lr,Lz,Lvzeta,Lvr,Lvz,r_bc,z_bc,geometry,composition)
gav_dfnn = gyroaveraged_dfnn_sym(Lr,Lz,r_bc,z_bc,geometry,composition) # gyroaverage < dfnn > in vpa vperp coordinates
gav_dfnn = gyroaveraged_dfnn_sym(Lr,Lz,Lvzeta,Lvr,Lvz,r_bc,z_bc,geometry,composition) # gyroaverage < dfnn > in vpa vperp coordinates

# define derivative operators
Dr = Differential(r)
Expand Down