Hello FEniCS Team, I try to solve a convection-diffusion-reaction NonlinearProblem which involves a network consisting of many 1D line elements in 3D space. Because there are so many line elements in 3D space, I need to compute the directional derivative rather than ‘u.dx(0)’. But there are some problems related to ‘ufl.sqrt’. Any comments/suggestions would be greatly appreciated. I use “DOLFINx version: 0.6.0.0 based on GIT commit: 72c70a59b8bed7cc5caa4c68a84182c9c4be87cf of GitHub - FEniCS/dolfinx: Next generation FEniCS problem solving environment”. Thanks a lot!

Here is a Minimal Reproducible Example:

####################################

```
import numpy as np
import dolfinx
from ufl import (TrialFunction, TestFunction, dx, ds, dS, dot, inner, jump, div, grad,
lhs, rhs, as_vector, FacetNormal, Measure, derivative)
from ufl.operators import sqrt
from dolfinx.fem import *
from dolfinx.mesh import create_interval, exterior_facet_indices, locate_entities, meshtags
from mpi4py import MPI
from petsc4py.PETSc import ScalarType
from dolfinx.fem.petsc import NonlinearProblem, assemble_vector, assemble_matrix, create_vector, apply_lifting, set_bc
from dolfinx.nls.petsc import NewtonSolver
# Compute directional derivative
def dudx(u):
directional_derivative_u = sqrt(dot(grad(u), grad(u)))
return directional_derivative_u
# Create mesh
L = 1.0; dt = 0.1
domain = create_interval(MPI.COMM_WORLD, 10, (0.0, L)) # In fact, I have many 1D lines in 3D space!
V = FunctionSpace(domain, ('CG', 1))
du = TrialFunction(V)
v = TestFunction(V)
u = Function(V)
u0 = Function(V)
# BCs
dofs_L = locate_dofs_geometrical(V, lambda x: np.isclose(x[0], 0.0))
dofs_R = locate_dofs_geometrical(V, lambda x: np.isclose(x[0], L))
uL, uR = Function(V), Function(V)
uL.vector.array = 0.0
uR.vector.array = 0.0
bcL, bcR = dirichletbc(uL, dofs_L), dirichletbc(uR, dofs_R)
bcs = [bcL, bcR]
F = v*(u-u0)/dt*dx + 0.01*v*dudx(u)*dx + 0.2*dudx(u)*dudx(v)*dx + 0.001*u/(1.0+u)*v*dx
a, L = lhs(F), rhs(F)
Jac = derivative(F, u, du)
# initial conditions
formFb = form(rhs(F))
init_cond = create_vector(formFb)
with init_cond.localForm() as loc:
loc.set(ScalarType(50.0))
u0.interpolate(init_cond)
problem = NonlinearProblem(F, u, bcs=bcs, J=Jac)
solver = NewtonSolver(domain.comm, problem)
solver.convergence_criterion = "incremental"
solver.rtol = 1e-6
solver.report = True
num_steps = 10
for n in range(num_steps):
num_its, converged = solver.solve(u)
print(f"--> Time step {n}, Number of iterations {num_its}")
```

####################################

Here is part of the error message:

…

Found Argument in , this is an invalid expression.

ERROR:UFL:Found Argument in , this is an invalid expression.