dolfinX : How to calculate a gradient/derivative and get it as a numpy.array?

Hello everyone,

Dolfinx version is 0.7.3 and have been installed with conda.
Other version I need in the minimal code (at the bottom of the post) are :

numpy : 1.26.4
mpi4py : 3.1.6
ufl : 2.23.2.0

Context : I’m solving an acoustic equation with Neumann condition with an iteratived domain decomposition method (no overlapping). I have 2 subdomains named \Omega _1 and \Omega_2. The common boundary (interface) between the two subdomains is named \Gamma .

\begin{equation} \mathcal{L}_{1}[\mathbf{u}_{1}^{m}] = 0 \; \; \; \text{on} \; \; \Omega _{1} \\ \mathcal{L}_{2}[\mathbf{u}_{2}^{m}] = 0 \; \; \; \text{on} \; \; \Omega _{2} \\ \mathbf{t}_{1}^{m} + \alpha \, \mathbf{u}_1^{m} = -\mathbf{t}_{2}^{m-1} + \alpha \, \mathbf{u}_2^{m-1} \; \; \; \text{ sur } \Gamma \\ \mathbf{t}_{2}^{m} + \alpha \, \mathbf{u}_2^{m} = -\mathbf{t}_{1}^{m-1} + \alpha \, \mathbf{u}_1^{m-1} \; \; \; \text{ sur } \Gamma \\ \text{ } \\ \text{where } \; \mathbf{t}_{i}^{m} = \nabla [\mathbf{u}_{i}^{m}]\cdot \mathbf{n}_{i}^{ } \end{equation} \; \text{ and } \; \mathbf{n}_{i}^{ } \; \text{ a normal vector to } \; \Gamma .

For example, the weak formulation is the domain 2 is :

\mathcal{A}_{2}^{m} (\mathbf{u}_{2}^{m},\mathbf{v}_{2}^{}) = \displaystyle {\int _ {\Gamma } } \left( - \mathbf{t} _1 ^{m-1} + \alpha \mathbf{u} _1 ^{m-1} \right) \cdot \mathbf{v}_2

Problem :
I must calculate the gradient at each iteration and get it with a numpy.array type to initialise the next calculation.

What i find to this purpose :

# Define a problem
problem = LinearProblem(a_form, L_form, bcs = [ ])

# Solution
uh = problem.solve()

It is possible to use .x.array method on the solution uh of a probleme :

array = uh.x.array

Then, the gradient and the traction are defined by :

def sigma(u):
    return ufl.nabla_grad(u)

traction = ufl.dot(sigma(uh),ufl.as_vector([1.0,0.0,0.0]))

Type of these three objets are :

 type(uh) : <class 'dolfinx.fem.function.Function'>
 type(sigma(uh)) : <class 'ufl.differentiation.NablaGrad'>
 type(traction) : <class 'ufl.tensoralgebra.Dot'>

This is to say i must convert these two last types in dolfinx.fem.function.Function to use the .x.array method.

So, my idea was to create a trivial system to get the gradient such that :

\displaystyle{\int _\Omega} \mathbf{t} \cdot \mathbf{v} = \displaystyle{\int _\Omega} \mathbf{traction} \cdot \mathbf{v} \;

where \mathbf{t} is a trial function, \mathbf{v} is a test function and \mathbf{traction} the second member.

Then \mathbf{t} will have the <dolfinx.fem.function.Function> type and \mathbf{t}.x.array will initialise the next calculation.

Question Are there other ways to switch a type to a dolfinx.fem.function.Function or to get an array associated to it ? Especially because this approach double my computation time …

Here is the minimal code : where I need the traction to be accessible in a numpy.array

## Library

import numpy as np                               # 1.26.4
from dolfinx import fem, la, mesh, plot          # 0.7.3
from mpi4py import MPI                           # 3.1.6
import ufl                                       # 2.23.2.0
from dolfinx.fem.petsc import LinearProblem      # 0.7.3

## Domain parameter

length, height = 10, 3
Nx, Ny = 80, 60
extent = [[0.0, 0.0], [length, height]]
domain = mesh.create_rectangle(MPI.COMM_WORLD, extent, [Nx, Ny], mesh.CellType.quadrilateral)

## Finite Element Space/functions

V = fem.FunctionSpace(domain,("P",2,(3,)))

## Define a Fucntion

def Solution(x):
    X = np.zeros_like(x)
    X[0] = x[0]
    X[1] = np.cos(x[0]+x[2]**2)
    X[2] = 0.0
    return X

sol = fem.Function(V)
sol.interpolate(Solution)

## Traction of sol 

def sigma(u):
    return ufl.nabla_grad(u)

traction = ufl.dot(sigma(uh),ufl.as_vector([1.0,0.0,0.0]))

Yes,
You can use:

1 Like

Dear Dokken,

Thank you for your quick answer.

I have followed the first proposition cause the second looks a bit harder.

Now, it works very well and divide my old computation time by 2. For example (same computation) :

Time before (s) Time now (s)
Solve the problem : 106 106
Get the traction/gradient : 106 0.6

Thanks again.