Dear all,

I’m trying to solve Navier Stokes equations in a 2D toy problem with time-dependent Dirichlet inflow boundary conditions using a monolithic solver. The code I created does not produce any error but when it comes to the velocity solution (at `t_end`

for instance) it seems like the inflow boundary condition is not applied correctly to all the nodes, as you can see in the following picture:

The same phenomenon happens also when refining the mesh on what seem to be random nodes on the boundary. I ran the code in serial, if this may be of any interest.

Here you can find the minimal working example that produces this result:

```
import dolfinx as df
from mpi4py import MPI
from petsc4py import PETSc
import ufl
import numpy as np
from dolfinx.mesh import create_unit_square, locate_entities, meshtags
from dolfinx.fem import (Constant, Function, FunctionSpace, assemble_scalar, dirichletbc, form, locate_dofs_geometrical,
locate_dofs_topological)
from dolfinx.fem.petsc import assemble_matrix_block, assemble_vector_block, LinearProblem
from dolfinx.fem.petsc import assemble_matrix, assemble_vector, apply_lifting, create_vector, set_bc
from ufl import (FacetNormal, FiniteElement, MixedElement, Identity, TestFunction, TrialFunction, VectorElement,
TrialFunctions, TestFunctions, Measure, div, grad, dot, ds, dx, inner, lhs, nabla_grad, rhs, sym)
msh = create_unit_square(MPI.COMM_WORLD, 16, 16)
# set fluid parameters
mu = 0.035
rho = 1.0
Re = 1200.
# spaces
VE = VectorElement('P', msh.ufl_cell(), 2)
PE = FiniteElement('P', msh.ufl_cell(), 1)
elem = MixedElement([VE, PE])
W = FunctionSpace(msh, elem)
# time settings
t_end = 0.2
ts = 20
Dt = t_end/ts
t = 0. # initial time
# boundary conditions
def inflow(x):
return np.isclose(x[1], 0.)
class time_dependent_inflow:
def __init__(self):
self.t = 0.0
def eval(self, x):
return np.stack((np.zeros(x.shape[1]), np.full(x.shape[1], self.t*20)))
V0, _ = W.sub(0).collapse()
inflow_expr = time_dependent_inflow()
inflow_expr.t = 0.
bc_fun, _ = Function(W).split()
bc_fun.interpolate(inflow_expr.eval)
inflow_dofs = locate_dofs_geometrical((W.sub(0),V0), inflow)
bcu = dirichletbc(bc_fun, inflow_dofs, W.sub(0))
# Set noslip BC
def walls(x):
return np.logical_or(np.isclose(x[0], 0.), np.isclose(x[0], 1.))
u_bc = Function(V0)
facets = df.mesh.locate_entities_boundary(msh, 1, walls)
bc_noslip = dirichletbc(u_bc, df.fem.locate_dofs_topological((W.sub(0), V0), 1, facets), W.sub(0))
bcs = [bcu, bc_noslip]
bnd = [(1, lambda x: np.isclose(x[0], 1.)),
(2, lambda x: np.isclose(x[0], 0.)),
(3, lambda x: np.isclose(x[1], 1.)),
(4, lambda x: np.isclose(x[1], 0.))]
facet_indices, facet_markers = [], []
fdim = msh.topology.dim - 1
for (marker, locator) in bnd:
facets = locate_entities(msh, fdim, locator)
facet_indices.append(facets)
facet_markers.append(np.full_like(facets, marker))
facet_indices = np.hstack(facet_indices).astype(np.int32)
facet_markers = np.hstack(facet_markers).astype(np.int32)
sorted_facets = np.argsort(facet_indices)
facet_tag = meshtags(msh, fdim, facet_indices[sorted_facets], facet_markers[sorted_facets])
ds = Measure("ds", domain=msh, subdomain_data=facet_tag)
# variational form
w_n = Function(W)
u_n, p_n = w_n.split()
u, p = TrialFunctions(W)
v, q = TestFunctions(W)
n = FacetNormal(msh)
h = ufl.CellDiameter(msh)
problem_type = 'NS'
if problem_type == 'NS':
a = rho*inner( u,v )*dx - Dt*( - rho*inner( grad(u)*u_n, v )*dx - mu*inner( grad(u),grad(v) )*dx + p*div(v)*dx - q*div(u)*dx)
L = rho*inner( u_n,v )*dx
# add backflow stabilization
boundary_id = [3]
for el in boundary_id:
a += 0.25*rho*Dt*( abs( inner( u_n,n ) - abs(inner( u_n,n )) )*inner( u,v )*ds(el) )
# time loop
while t < t_end:
t += Dt
inflow_expr.t = t
bc_fun.interpolate(inflow_expr.eval)
problem = LinearProblem(a, L, bcs, petsc_options={"ksp_type": "preonly",
"pc_type": "lu",
"pc_factor_mat_solver_type": "mumps"})
w_n = problem.solve()
u_n, p_n = w_n.split()
u_n.name = 'u'
# write function for visualization
filepath = './u_n.xdmf'
with df.io.XDMFFile(msh.comm, filepath, 'w') as xdmf:
u_n.x.scatter_forward()
xdmf.write_mesh(msh)
xdmf.write_function(u_n)
```

Any advice on how to fix it would be great.

Thanks in advance!