How to return "ufl.algebra" instead of "numpy.ndarray"?

Hi Jorgen,

in general I am trying to define material properties (as a function of solution) using data points instead of polynomials in dolfinx, see here. The solver works fine with the “polynom” function (which has an output type: ufl.algebra). However if I use the “interpolation” function to define properties (output type: numpy.ndarray) the solver gives the following error: “AttributeError: ‘numpy.ndarray’ object has no attribute ‘subdomain_data’”. Therefore I though that if the “interpolation” function returns exactly the same “ufl.algebra” data type it could work.

MWE:

import dolfinx
import numpy as np
import ufl
from mpi4py import MPI
import scipy.interpolate as sp

mesh = dolfinx.RectangleMesh(MPI.COMM_WORLD, [np.array([0, 0, 0]), np.array([2, 2, 0])], [1, 1])
V = dolfinx.FunctionSpace(mesh, ("CG", 1))
uh = dolfinx.Function(V)
u_n = dolfinx.Function(V)

# def polynom(u):
#     return 0.00004*u + 0.02

def interpolation(var):
    data_points = sp.interp1d([0,1000],[0.02,0.06])
    return data_points(var.vector.array)

v = ufl.TestFunction(V)
f = dolfinx.Constant(mesh, 0)

# F = polynom(uh)*ufl.dot(ufl.grad(uh), ufl.grad(v))*ufl.dx
F = interpolation(uh)*ufl.dot(ufl.grad(uh), ufl.grad(v))*ufl.dx

problem = dolfinx.fem.NonlinearProblem(F, uh)
solver = dolfinx.NewtonSolver(MPI.COMM_WORLD, problem)
n, converged = solver.solve(uh)

Can I somehow adapt the “interpolation” function, such that it could be used for the existing nonlinear solver?

I could add:

q_ = dolfinx.Function(V)
q_.vector.array[:] = interpolation(uh)
F = q_*ufl.dot(ufl.grad(uh), ufl.grad(v))*ufl.dx

which doesn’t arise any error, however does not update material properties after each iteration as I showed here