the codes will break down after several time steps.
I don’t know reasons. I guess this syntax add exact one 4th order polynomial in the interior of tetrahedral element.

However, I found

velocity.SetOrder(TRIG,3)

is stable, can I understand it as adding a bubble function for each face of terahedral elements? Since each face just is a 2D domain.

Hello, it appears that the “Set” method of functions in 3D Mini element space is yielding unexpected “nan” values. Please refer to the example below, where the output is:

2031 nan in u_mini
0 nan in u_p1

The version of ngsolve used is NGSolve-6.2.2305.

from ngsolve import *
from netgen.occ import OCCGeometry, Box
import numpy as np
geo = OCCGeometry(Box((0,0,0), (1, 1, 1)))
mesh = Mesh(geo.GenerateMesh(maxh=0.1))
V1 = VectorH1(mesh, order=1)
V2 = VectorH1(mesh, order=1)
V2.SetOrder(TET, 4)
V2.Update()
u = CF((x, y, z))
u_p1 = GridFunction(V1)
u_mini = GridFunction(V2)
u_p1.Set(u)
u_mini.Set(u)
print(f"{len(np.where(np.isnan(u_mini.vec))[0]):4d} nan in u_mini")
print(f"{len(np.where(np.isnan(u_p1.vec))[0]):4d} nan in u_p1")