-
Notifications
You must be signed in to change notification settings - Fork 50
Add hyperresistivity term to electric field in fieldsolver #1157
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: dev
Are you sure you want to change the base?
Conversation
…th, add switch to select whether to accelerate with hyper or not
ursg
left a comment
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.
Nice! Pretty straightforward implementation, too.
If only we had the fsgrid refactor already in place, so this wouldn't be such a huge amount of boilerplate code when passing the new fsgrid struct everywhere. But such is life.
| std::array<Real, fsgrids::bfield::N_BFIELD> * centPerB = perBGrid.get(i,j,k); | ||
| std::array<Real, fsgrids::bgbfield::N_BGB> * centBgB = BgBGrid.get(i,j,k); | ||
| std::array<Real, fsgrids::moments::N_MOMENTS> * centMom = momentsGrid.get(i,j,k); | ||
|
|
||
| Real Bmag = sqrt((centBgB->at(fsgrids::bgbfield::BGBX)+centPerB->at(fsgrids::bfield::PERBX))* | ||
| (centBgB->at(fsgrids::bgbfield::BGBX)+centPerB->at(fsgrids::bfield::PERBX)) + | ||
| (centBgB->at(fsgrids::bgbfield::BGBY)+centPerB->at(fsgrids::bfield::PERBY))* | ||
| (centBgB->at(fsgrids::bgbfield::BGBY)+centPerB->at(fsgrids::bfield::PERBY)) + | ||
| (centBgB->at(fsgrids::bgbfield::BGBZ)+centPerB->at(fsgrids::bfield::PERBZ))* | ||
| (centBgB->at(fsgrids::bgbfield::BGBZ)+centPerB->at(fsgrids::bfield::PERBZ)) | ||
| ); | ||
| rhoq = centMom->at(fsgrids::moments::RHOQ); | ||
| limitedRhoq = (rhoq <= Parameters::hallMinimumRhoq ) ? Parameters::hallMinimumRhoq : rhoq ; |
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.
Note that these are now simply taken from their locations on the yee-lattice and used "as is", with no balsara interpolation or volumetric averages.
Since we are deaeling with wiggling on the grid-scale here, should we do better than that? Or can we afford to not care, as long as we get the damping order-of-magnitude correct so those wiggles never even grow?
Anyway, it might at least warrent a comment, that the "center" thing there is not the cell centre.
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.
My logic for doing it on the Yee-lattice without any sort of interpolation or averaging was that this is essentially an artificial term, so as long as it does what we need it to, does it need to be physically accurate?
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.
Thing is, this code computes a Bmag value and calculates something from that and the charge density - and charge density is a spatial-cell-average value. Thus, to compare apples with apples, it should also calculate Bmag from (Balsara-reconstructed) cell-averages IMO.
(See if the ordering means those values have already been calculated?)
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.
Getting edge-rhoq is trivial and likely already done elsewhere (for Hall) hopefully. Getting the Balsara-reconstructed, correct-order (2nd spatial order to be consistent with the remaining things) would probably require Hall-term-grade arithmetics which was famously handled with wxMaxima in its time and made into code with a careful set of sed commands...
parameters.cpp
Outdated
| "The minimum CFL limit for field propagation. Used to set timestep if dynamic_timestep is true.", 0.4); | ||
| RP::add( | ||
| "fieldsolver.ohmHyperFactor", | ||
| "Factor to multiply hyperresistivity term by. Default: 1 (damping becomes strong at grid scale)", |
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.
For the future uninitiated, it might be helpful here to add a comment of the sort
"more means larger scales are being damped".
fieldsolver/ldz_hyper.cpp
Outdated
| int computeTimerId {phiprof::initializeTimer("Ehyper compute cells")}; | ||
|
|
||
| phiprof::Timer mpiTimer {"Hyper field update ghosts MPI", {"MPI"}}; | ||
| dPerBGrid.updateGhostCells(); |
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 now leads to a duplicate ghost update of dPerBGrid if both Ohm and Hyperresistivity are active
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 catch, I had that in at some point but apparently removed it.
Now it's back in
|
|
||
| default: | ||
| cerr << __FILE__ << ":" << __LINE__ << "Derivative component should be 0 (xx), 1 (yy), or 2 (zz)." << endl; | ||
| break; |
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.
I would suggest aborting here
|
|
||
| default: | ||
| cerr << __FILE__ << ":" << __LINE__ << "Curl component should be 0 (x), 1 (y), or 2 (z)." << endl; | ||
| break; |
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.
Same
Here's the draft PR for the hyperresistivity term, with latest dev merged in and extra unnecessary stuff removed.
I've tested that it compiles (on carrington).