BOUT-dev icon indicating copy to clipboard operation
BOUT-dev copied to clipboard

FCI: taking derivative without communicating first is resulting in confusing error

Open dschwoerer opened this issue 4 years ago • 4 comments

        Option grid = data/gridv4.128.32.fci.tmp.nc (data/BOUT.inp)
terminate called after throwing an instance of 'BoutException'
  what():  FCI method cannot transform into field aligned grid

With backtrace:

Catchpoint 1 (exception thrown), 0x00002aaabaa51492 in __cxa_throw () from /lib64/libstdc++.so.6
(gdb) bt
#0  0x00002aaabaa51492 in __cxa_throw () from /lib64/libstdc++.so.6
#1  0x00002aaab8b443d6 in FCITransform::toFieldAligned (this=<optimized out>, UNUSED_f=..., UNUSED_region=...) at fci.hxx:111
#2  0x00002aaab8cbc172 in toFieldAligned<Field3D> (region="RGN_NOX", f=...) at ../../include/field.hxx:337
#3  bout::derivatives::index::DDY<Field3D> (f=..., outloc=outloc@entry=CELL_LOC::deflt, method="DEFAULT", region="RGN_NOBNDRY") at ../../include/bout/index_derivs_interface.hxx:209
#4  0x00002aaab8cada2a in DDY (f=..., outloc=outloc@entry=CELL_LOC::deflt, method="DEFAULT", region="RGN_NOBNDRY") at derivs.cxx:75
#5  0x00002aaab8c86703 in Grad_perp (f=..., outloc=outloc@entry=CELL_LOC::deflt, method="DEFAULT") at vecops.cxx:109
#6  0x00002aaab89d14e7 in __pyx_pf_8boutcore_36Grad_perp (__pyx_v_field=0x2aacc4b9bd30, __pyx_v_outloc=<optimized out>, __pyx_v_method=0x2aacc56fceb0, __pyx_v_region=<optimized out>, __pyx_self=<optimized out>) at boutcore.cpp:42517

Would it make sense to add a try: statement to catch this issue and write a nicer error?

dschwoerer avatar Feb 16 '21 10:02 dschwoerer

Maybe simpler to just add to the error message that it might be due to this?

johnomotani avatar Feb 16 '21 10:02 johnomotani

Maybe simpler to just add to the error message that it might be due to this?

Unfortunately this is not just an issue with this particular operator. It could be helpful to add to the "cannot transform" error that a communication could be helpful, but this might not always be the best solution.

bshanahan avatar Feb 16 '21 13:02 bshanahan

The correct way to communicate with FCI is:

var.applyBoundary("neumann");
mesh->communicate(var);
var..applyParallelBoundary("parallel_neumann");

This can be tricky to get right. Especially communicating again results in invalid field data, as field lines leaving the domain are invalidated.

Any ideas how to make this less error prone, besides documenting this?

dschwoerer avatar Sep 19 '22 14:09 dschwoerer

Perhaps this could somehow be put in the FCI ParallelTransform class, so that any calls to mesh->communicate invoke these boundary condition setters? The invalidation is tricky though: Is neumann always the right boundary condition to set?

bendudson avatar Sep 19 '22 14:09 bendudson