Test difops - #1315
Open
johnomotani wants to merge 73 commits into
Open
Test difops#1315johnomotani wants to merge 73 commits into
johnomotani wants to merge 73 commits into
Conversation
When iterating over flags from h5cc, the HDF5_{CPPFLAGS,LDFLAGS,LIBS}
variables were assembled in reverse order. This causes problems when
linking against static libraries.
Fixes #961
Fix order of HDF5 flags
If the simulation is restarted (`restart`), and the old dump files have already been squashed, dump on restart is enabled by default. This patch checks whether the two times are the same, and if so, does not duplicate that time slice.
For some reason the time index is called `x2` - and thus this condition is not matched. This fixes this, as the newly created output should be more reliable in this regard.
`glob` was previously not correctly imported, if append mode was disabled, but delete was enabled.
* master: AX_LIB_HDF5: fix order of flags
Advective and flux derivatives, which take v and f as inputs, had cases using yup/ydown fields for one of v and f, but not the other. These do not make sense as multiplication of a field in field-aligned coordinates with another in non-field-aligned coordinates is incorrect. Delete last of these cases and fall back to converting both v and f to field-aligned if either does not have yup/ydown fields. Also add some comments where mixed cases were removed before.
In several y-derivatives where we have to convert the input to field-aligned coordinates, the result was not transformed back to the original coordinates. Also remove uses of global 'mesh' in Vpar_Grad_par_LCtoC()
Was previously a typo that resulted in using non-field-aligned f in VDDY, which is incorrect.
The location of the result will in future need to be set before shifting it from field-aligned coordinates, if this is necessary. So it is safer to set the location to outloc as soon as 'Field3D result' is declared. Also in a couple of places just 'return apply*diff(...)' instead of creating an unneeded intermediate variable 'result'.
* next: Move handling of outloc==CELL_DEFAULT above location checks Use ASSERT1 instead of ASSERT0 for location checking in derivatives Remove unnecessary conditionals in applyXdiff() Tidy up division in D2DYDZ Simplify location checking in index_derivs.cxx Remove mixed yup/ydown-fieldaligned cases in advective/flux derivatives Move checks for no guard cells to ASSERT1 at top of functions Add ASSERT1s to check that location of Coordinates is outloc Add missing outloc argument in indexFDDY(Field2D,Field2D) Check 'vloc != diffloc' in indexFDDX Use RGN_NOY for y-derivative in D2DXDY(Field2D) Fix location of Bxy in Div_par_flux Fix location checking in VDD* and FDD* Fix staggering in SPLIT branch of FDD* Fix staggering in Div_par_K_Grad_par Don't use global mesh in bracket_location Correct dy in D2DYDZ Fix sign in Div_par_flux Fix Bxy location in Div_par
…-update Bugfixes from shiftedmetric-update
Also rename timer wrappers with timed_ prefix
Someone nice who knows IDL, Matlab/Octave and/or Mathematica could maybe actually fix these routines and submit pull requests :)
Test squashoutput shell wrapper
d7919
reviewed
Oct 11, 2018
| // from 2pi->0 | ||
| bool need_comms = !mesh->sourceHasYGuards(); | ||
|
|
||
| if (localmesh->get(dx, "dx", 0., need_comms)) { |
Member
There was a problem hiding this comment.
I'm not sure, but it looks like need_comms could be a member of source/griddatasource or at least mesh as it doesn't seem to need anything in coordinates?
Contributor
Author
There was a problem hiding this comment.
think i'm missing your point... need_comms is a local variable, the functionality is in sourceHasYGuards() which is a member of Mesh.
Member
There was a problem hiding this comment.
If sourceHasGuards is all that is required then couldn't we put the call to this inside mesh::get so that mesh::get does the right thing?
Contributor
Author
There was a problem hiding this comment.
I was originally going to load zShift in Coordinates (instead of only in ShiftedMetric), and pass 'false' for that, but I think there might be a better solution, and in any case is better to add the extra argument only when it's actually needed, so I'll put the call inside Mesh::get as you suggest.
Member
There was a problem hiding this comment.
This should all only be done once (or a few times) during initialisation so even if we pay a few extra communications that aren't strictly needed (as we overwrite guards etc.) then I think that is ok if it helps make the code cleaner (e.g. not having to repeat the default argument everywhere).
d7919
reviewed
Oct 11, 2018
| DIFF_METHOD UNUSED(method)) { | ||
| DIFF_METHOD method) { | ||
| TRACE("Coordinates::Grad_par( Field2D )"); | ||
| ASSERT1(location == outloc || outloc == CELL_DEFAULT); |
Member
There was a problem hiding this comment.
You didn't change this here but I'm not sure this assert is correct -- I think it should be (location == outloc || (outloc == CELL_DEFAULT && location == var.getLocation())`?
d7919
reviewed
Oct 11, 2018
| CELL_LOC outloc, | ||
| DIFF_METHOD UNUSED(method)) { | ||
| DIFF_METHOD method) { | ||
| ASSERT1(location == outloc || outloc == CELL_DEFAULT); |
d7919
reviewed
Oct 11, 2018
| const Field2D Coordinates::Grad2_par2(const Field2D &f, CELL_LOC outloc) { | ||
| const Field2D Coordinates::Grad2_par2(const Field2D &f, CELL_LOC outloc, DIFF_METHOD method) { | ||
| TRACE("Coordinates::Grad2_par2( Field2D )"); | ||
| ASSERT1(location == outloc || outloc == CELL_DEFAULT); |
Makes Mesh.coordinates() a regular method instead of a property and adds an optional location argument.
Backticks are special characters for bash, so it's safer to avoid them. Use single quotes instead.
Use a bash function to generate repeated pattern of 'isinstance(a, Field3D)...'. The pattern is needed to check types in Python functions in order to call the right version of overloaded C++ functions.
This will be useful for tests written using boutcore that depend on the CHECK level in BOUT++, e.g. to know when exceptions should be expected.
bracket operator with method BRACKET_STD is converging faster than 2nd order for 64->128. It approaches closer to 2nd order for higher resolutions, but that uses too much memory. Instead, fudge the expected order of the operator so that the test passes.
Used by MMS/difops/runtest.py
Older versions of numpy.logspace(), such as the one on Travis, don't support the dtype argument.
The integrated shear may be discontinuous when poloidal angle goes from 2pi->0. In this case communication will give incorrect values in the y-boundary guard cells for quantities that depend on the integrated shear. The incorrect values cause bugs with interpolation to staggered grids and y-derivatives used to calculate, e.g., G1/G2/G3. Previously all 2D fields loaded from a grid were communicated. Now only communicate if the grid file does not contain y-boundary guard cells. This maintains backward compatibility although it is not 'correct'. Replace 'communicate' argument to Mesh::get() with 'allow_communicate' for Field3D. Use source->hasYGuards() inside Mesh::get(). Field3D is communicated by default only if source->hasYGuards() returns false, but is never communicated if allow_communicate is set to false.
Previously used Neumann boundary conditions for metric components (and a couple of other geometric quantities) to set the guard cells when interpolating to staggered grids. This caused things that depend on derivatives of the metric, like G1/G2/G3 to have very different, and wrong, values at the boundaries. Now extrapolate instead, with the same stencil as free_o3 boundary conditions, but extrapolating only the points past the boundary. An extra interpolation has to be added after calling interp_to to get the interpolated field at the upper/outer boundary.
When metrics are calculated on staggered grids, they are extrapolated into the boundary cells. Extrapolation requires 3 grid points. If there are fewer than 3, just copy the nearest value into the guard cells, since probably derivatives aren't being taken in that direction anyway. If only 1 guard cell and 1 grid point per processor in a certain direction are being used, then only have 2 points to extrapolate from on a processor, so cannot extrapolate to guard cells with 3rd order accuracy. In this case throw an exception because the user needs to change MXG/MYG or NXPE.
We have 'result = interp_to(f, location, RGN_NOBNDRY)' but then communicate result. Need to ensure result's data is unique. Otherwise result might be a duplicate of f (if no interpolation is needed, e.g. if interpolation is in the z-direction); then f would be communicated. Since this function is used on geometrical quantities that might not be periodic in y even on closed field lines (due to dependence on integrated shear), we don't want to communicate f. We will sort out result's boundary guard cells, but not f's so we don't want to change f.
Previously, fields at CELL_ZLOW had nothing done by any boundary condition operator. This commit passes them through to the unstaggered case, which is correct since there are no z-boundaries so no special handling is needed.
Requires adding a Mesh* pointer and getDataMesh() method to FieldData, so that FieldData::setBoundary() can use the local mesh. getDataMesh() method is needed in case the FieldData is constructed before the global mesh is created (i.e. when Field3D, etc. are declared in global scope). It must not be called getMesh() because that would clash with the method of Field.
johnomotani
force-pushed
the
test-difops
branch
from
October 12, 2018 18:42
2a83dd0 to
72d293e
Compare
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
4 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
There are quite a few commits here, but it does add a test and a bugfix, so worth considering for v4.2?
The test is an MMS test using boutcore for all the differential operators and derivatives, with all combinations of inloc/outloc, and vloc if necessary. It checks that all the operators either converge at 2nd order (a couple use U1, so first order for those) or throw an exception (if CHECK>0). The full test is pretty long (oddly, testing for Exceptions via cython is slower than most of the successful tests), so there's a short version for the test-suite that only tests the difops (not the derivatives), skips the checks for exceptions and skips operators that require 3d grid refinement. Arguably, the short version is still too slow, but if so we need to decide what to take out (maybe operators that few people use, or some unlikely locations like ZLOW/ZLOW).
Summary of changes:
boutcore.bout_CHECKvariable to access the value ofCHECK; some extra keywords to allow differential operators to raise C++ exceptions as Python exceptions.