Unexpected kernel crash with essential boundary conditions

[gthyagi]

Hi @lmoresi @julesghub @knepley ,

Issue: It is currently not possible to implement a essential boundary condition using a mesh variable that is not available in the Stokes object. In the example below, v_uw is the primary mesh variable available in the Stokes object, while v_ana is the mesh variable used just to store the top boundary velocities.

Issue lies here: https://github.com/underworldcode/underworld3/blob/d0b6f3b71ca9dea282571e7ad470c7de6b90925f/src/underworld3/cython/petsc_generic_snes_solvers.pyx#L2720C9-L2750C1

Here the block of code to replicate the essential boundary issue.

import underworld3 as uw
from underworld3.systems import Stokes
import sympy

mesh = uw.meshing.UnstructuredSimplexBox(cellSize=0.1)

# mesh variables
v_uw = uw.discretisation.MeshVariable('V_u', mesh, mesh.data.shape[1], degree=2)
p_uw = uw.discretisation.MeshVariable('P_u', mesh, 1, degree=1, continuous=True)
v_ana = uw.discretisation.MeshVariable('V_ana', mesh, mesh.data.shape[1], degree=2)

with mesh.access(v_ana):
    v_ana.data[:,0] = 1.0

# Create Stokes object
stokes = Stokes(mesh, velocityField=v_uw, pressureField=p_uw, solver_name="stokes_sph")
stokes.constitutive_model = uw.constitutive_models.ViscousFlowModel
stokes.constitutive_model.Parameters.viscosity = 1
stokes.saddle_preconditioner = 1.0/1.0

# +
# gravity
gravity_fn = -1.0 * mesh.CoordinateSystem.unit_e_1

# bodyforce term
stokes.bodyforce = 1.0*gravity_fn # 0.0 * unit_rvec
# -

# bc's
stokes.add_essential_bc(v_ana.sym, "Top") # v_ana.sym
stokes.add_essential_bc((sympy.oo,0.0), "Bottom")
stokes.add_essential_bc((0.0,sympy.oo), "Left")
stokes.add_essential_bc((0.0,sympy.oo), "Right")

stokes.solve(verbose=True, debug=False)