-
Notifications
You must be signed in to change notification settings - Fork 21
[WIP] Consolidate sheath boundaries #382
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
base: master
Are you sure you want to change the base?
Conversation
`phi` must be calculated in field-aligned coordinates, and transformed to unaligned
Lower-y loop used wrong cell face
Possibly faster, certainly cleaner
Didn't flip noflow BCs at same time, so was not exactly symmetric
| const BoutReal _gamma_i = 2.5 + (0.5 * Mi * C_i_sq / tisheath); | ||
|
|
||
| const auto _q = | ||
| ((_gamma_i - 1 - 1 / (adiabatic - 1)) * tisheath - 0.5 * Mi * C_i_sq) | ||
| * nisheath * visheath; |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This is different from the form of gamma_e. Notably there's an explicit 2.5 here. Is this correct?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Hmmm..... I don't know of any reason why the gamma expression should be any different than in sheath_boundary_simple, apart from that sheath_boundary had it hardcoded. @bendudson needs to confirm
|
|
The |
* master: (54 commits) Silence some compiler warnings. Increase the test tolerance to work around compiler differences (seen with -O0 vs. -O2 and gcc@12+). fixed_velocity: Add description of mesh array to manual Typo in assert message. fixed_velocity: Read velocity from mesh Extend test input densities, temperatures outside clipped ranges. ADAS CX tests. Use common input ranges for all reaction regression tests. Rename test base class for ADAS izn/rec. Move AmjuelCXTest::generate_state() up to a parent class, so it can reused for ADAS. Rejig ReactionTest::sources_regression_test to check species data stored in (sub-)sections, including collision freqencies. evolve_pressure: Initialize flow_ylow if diagnose but no velocity Add H isotope CX tests. Move generate_data() for ADAS, Amjuel izn and rec tests into a separate class. Consistent naming of reaction regression data files. ADAS ionisation, recombination reaction regression tests. Correct a mistake in the docs Makefile that meant the doxygen target wasn't declared as .PHONY. Tidy up ReactionTest::generate_state(). Revert "Rename doxygen Makefile target to avoid clash with docs subdirectory name." Add missing header. ...
e94eb10 to
44e3a0e
Compare
Add to loosen tolerance on sheath test slightly as higher `Ge` test has smaller `Te,Ve` so larger relative error. Difference between old and new implementation does appear to be from build up of numerical errors rather than anything systematic. New test at intermediate `Ge` to show this.
44e3a0e to
b7b0467
Compare
Well yes, this is technically true, but |
Ah, maybe it would be wise to raise a warning or something in that case? That's a separate issue I suppose |
…s-db merge main into consolidate-sheath-boundaries
| const BoutReal Ni_im = limitFree(Ni[ip], Ni[i], density_boundary_mode); | ||
| const BoutReal Ti_im = limitFree(Ti[ip], Ti[i], temperature_boundary_mode); | ||
| const BoutReal Te_im = limitFree(Te[ip], Te[i], temperature_boundary_mode); | ||
|
|
||
| // Calculate sheath values at half-way points (cell edge) | ||
| const BoutReal nisheath = 0.5 * (Ni_im + Ni[i]); | ||
| // electron temperature | ||
| const BoutReal tesheath = floor(0.5 * (Te_im + Te[i]), 1e-5); | ||
| // ion temperature | ||
| const BoutReal tisheath = floor(0.5 * (Ti_im + Ti[i]), 1e-5); |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Does it make sense to first interpolate in log space, and then afterwards linear?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
What do you mean? Are you talking about interpolating twice?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yes.
That is what the code does currently.
Ti_im is interpolated in log-space here
and then in linear space, here, 7 lines later.
| case SheathLimitMode::exponential_free: | ||
| fp = SQ(fc) / fm; // Exponential | ||
| case SheathLimitMode::linear_free: | ||
| fp = 2.0 * fc - fm; // Linear |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Are you not missing the break statement?
Does this not always apply linear_free here?
Maybe a unit test for this would be nice.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Good spot, thanks! I thought this was a warning, but it turns out we don't enable warnings by default
| }; | ||
|
|
||
| /// Information about a point on a Y boundary | ||
| struct YBoundary { |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Would it make sense to put this somewhere separate so it can be reused in other components?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yes -- eventually I think this will form the basis of refactoring all the boundary conditions in BOUT++. Once this is in and I start on the next component, I'll see where it's sensible to move it to
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Can we please switch instead to YBoundary instead?
This allows to write sheath boundary code for FA and FCI!
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Includes #368.
e94eb10 is the commit mergingsheath_boundaryintosheath_boundary_simple.Current task is to see about mergingsheath_boundary_insulatingin as well.Merges
sheath_boundary,sheath_boundary_simple, andsheath_boundary_insulatinginto a single implementation with compile-time switches between the different conditions.Is
SheathKind::normala sensible name for the "non-simple" case? (Note that this is an internal name only, and not user-visible).TODO: