SNES solver with neutral_full_velocity model

[malamast]

Hi everyone,

I try to test the neutral_full_velocity component with the SNES solver but it gives me an PETSC ERROR "Argument out of range".
The same model / input file works fine with CVODE.
The model is here: https://github.com/boutproject/hermes-3/blob/master/src/neutral_full_velocity.cxx#L76C3-L76C9

The output file says:

            SNES steady state solver
            3d fields = 4, 2d fields = 5 neq=28800, local_N=2880
    Creating vector
    Create SNES
    Setting Jacobian matrix sizes
            Option solver:stencil:square = 0 (default)
            Option solver:stencil:cross = 0 (default)
            Option solver:stencil:taxi = 2 (default)
    Pre-allocating Jacobian
    Marking non-zero Jacobian entries
    Assembling Jacobian matrix

the error file says:

      [9]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------
      [9]PETSC ERROR: Argument out of range
      [9]PETSC ERROR: New nonzero at (5,89) caused a malloc. Use MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE) to turn off this check

With the neutral_mixed component we solve only 3d fields = 7. Now, with the neutral_full_velocity component we have 3d fields = 4, 2d fields = 5. I use the latest "next" bout++ branch. Anyone have a clue what’s going down?

[bendudson]

Thanks @malamast ! This might be a bug in how SNES constructs the Jacobian for 2D fields: I don't think we've used evolving 2D fields with the SNES solver previously, so would not have seen this. The full velocity model assumes all neutral gas fields are axisymmetric, hence 5 2d fields (n, p, vx, vy, vz).