ammarhakim / gkyl

This is the main source repo for the Gkeyll 2.0 code. Please see gkeyll.rtfd.io for details.
https://gkeyll.readthedocs.io/en/latest/
54 stars 14 forks source link

pre-g0: MomentSpecies:updateInDirection raises nondescriptive error #174

Open skylord-a52 opened 3 months ago

skylord-a52 commented 3 months ago
function MomentSpecies:updateInDirection(dir, tCurr, dt, fIn, fOut, tryInv)
   local status, dtSuggested = true, GKYL_MAX_DOUBLE
   local tryInv_next = false
   if self.evolve then
      self:applyBc(tCurr, fIn, dir)
      assert(self:checkInv(fIn))

Line 495 of MomentSpecies.lua raises an error in... some cases, but because there's no descriptive text on the error and no documentation on what fIn or checkInv, it's impossible to tell what exactly causes the error. I'm guessing there's some issue with either a non-invertible matrix somewhere, or with the definition of equationInv in the species definition.

Heres my code that's causing the error. (Please note -- this is a bit of a mishmash of different examples at the moment, and because I'm still learning, some of the initial conditions may be nonphysical.)

https://pastebin.pl/view/42736129

When run, the console prints

*** LOAD ERROR ***
 ...e/andrea/gkylsoft/gkyl/bin/App/Species/MomentSpecies.lua:495: assertion failed!

with no further information.

JunoRavin commented 3 months ago

If you look at the function being assigned to checkInv, it is the invariant domain of the equation system being solved for. So for example, the invariant domain of Euler's equations is that the density and pressure are both positive definite.

This assertion, while perhaps not as descriptive as it could be, is just breaking out of the simulation flow when the invariant domain is not satisfied and there is no recovering from this violation (either because you initialized a simulation which violates the invariant domain, or because the simulation could not be prevented from driving outside of the invariant domain, for example if you try to do a fluid simulation with a true vacuum state on one side and the Riemann solve has difficulty resolving the transition from finite density to genuinely 0 density).

skylord-a52 commented 3 months ago

Sounds like it's an issue with physicality then. Thanks for your continued help! I'll adjust my initial conditions and see if I can get something that works a bit better.

skylord-a52 commented 3 months ago

It would be very helpful if the error not only said what the problem was (that some invariant was violated) but specifically which particular parameter was the issue. In the meantime I can use the debugger to figure things out, but a message to console saying "Density is negative [at x,y=...]" would make things much faster.