DiscreteOperators.build_gradient in dolfinx

dolfinx
1

Hello, I am a new user of dolfinx and I am going to convert the code in this link, https://github.com/ryanfreckleton/dolfin/blob/f8cd6c0b2722bef82f1d429a9e838a4491b30154/python/demo/undocumented/curl-curl/demo_curl-curl.py# , from dolfin to dolfinx version. However, I have trouble in line 121, which is “G = DiscreteOperators.build_gradient(V, P1)”. I cannot find the equivalent command in dolfinx. The code in dolfinx is provided below.

from mpi4py import MPI
from dolfinx import mesh, fem, plot
from dolfinx.fem import FunctionSpace, Function, Constant, Expression, dirichletbc
import numpy as np
from petsc4py import PETSc
from petsc4py.PETSc import ScalarType
import ufl
from ufl import TrialFunction, TestFunction
import pyvista
import dolfinx

domain = mesh.create_box(MPI.COMM_WORLD, [np.array([0,0,0]), np.array([10, 10, 1])],[20,20,4])

V = FunctionSpace(domain, ("Nedelec 1st kind H(curl)", 1))

v = TestFunction(V)
u = TrialFunction(V)

dbdt = Constant(domain, (0.0, 0.0, 1.0))
zero = Function(V)

T = Function(V)

tdim = domain.topology.dim
fdim = tdim - 1
domain.topology.create_connectivity(fdim, tdim)
boundary_facets = mesh.exterior_facet_indices(domain.topology)
boundary_dofs = fem.locate_dofs_topological(V, fdim, boundary_facets)
bc = dirichletbc(zero, boundary_dofs)

a = ufl.inner(ufl.curl(v), ufl.curl(u))*ufl.dx
L = -ufl.inner(v, dbdt)*ufl.dx

a_assembled = fem.petsc.assemble_matrix(fem.form(a))
a_assembled.assemble()
L_assembled = fem.petsc.assemble_vector(fem.form(L))
L_assembled.assemble()

ksp = PETSc.KSP()
ksp.create(PETSc.COMM_WORLD)

ksp.setType("cg")
ksp.setTolerances(rtol=1.0e-8, atol=1.0e-12, divtol=1.0e10, max_it=300)

pc = ksp.getPC()
pc.setType("hypre")
pc.setHYPREType("ams")

P1 = FunctionSpace(domain, ("Lagrange", 1))
G = dolfinx.cpp.fem.petsc.discrete_gradient(V, P1)

pc.setHYPREDiscreteGradient(as_backend_type(G).mat())

constants = [fem.Function(V) for i in range(3)]
for i, c in enumerate(constants):
direction = Constant(domain, [1.0 if i == j else 0.0 for j in range(3)])
c_expr = Expression(direction , V.element.interpolation_points())
c.interpolate(c_expr)

cvecs = [constant.vector for constant in constants]

pc.setHYPRESetEdgeConstantVectors(cvecs[0], cvecs[1], cvecs[2])

pc.setHYPRESetBetaPoissonMatrix(None)

ksp.setOperators(a_assembled)

ksp.setOptionsPrefix("eddy_")

opts = PETSc.Options()

opts.setValue("-eddy_ksp_monitor_true_residual", None)

ksp.setFromOptions()

ksp.solve(L_assembled, T.vector)

I receive the following error.

TypeError: discrete_gradient(): incompatible function arguments. The following argument types are supported:
    1. (V0: dolfinx::fem::FunctionSpace, V1: dolfinx::fem::FunctionSpace) -> mat

Invoked with: FunctionSpace(Mesh(VectorElement(Basix element (P, tetrahedron, 1, gll_warped, unset, False), 3), 2), Basix element (N1E, tetrahedron, 1, legendre, unset, False)), FunctionSpace(Mesh(VectorElement(Basix element (P, tetrahedron, 1, gll_warped, unset, False), 3), 2), Basix element (P, tetrahedron, 1, gll_warped, unset, False))

I would be grateful if you could please help me in this regard.

2

Please format your code with 3x`,
Ie.

```python
import dolfinx

# add code here
```

and i would suggest looking at:

ie change

to

G = dolfinx.cpp.fem.petsc.discrete_gradient(P1._cpp_object, V._cpp_object)
3

Thanks for the prompt response. The solution you provided worked for the minimal example. However, in the problem I am solving I have used mixed elements for defining function space as shown below:

Element_Mec = ufl.VectorElement("Lagrange", domain.ufl_cell(), 1, dim=3)
Element_EM_vec = ufl.FiniteElement("Nedelec 1st kind H(curl)", domain.ufl_cell(), 1)
Element_EM_sca = ufl.FiniteElement("Lagrange", domain.ufl_cell(), 1)

MixedElement = ufl.MixedElement(Element_Mec, Element_EM_vec, Element_EM_sca)

V = fem.FunctionSpace(domain, MixedElement)

num_subs = V.num_sub_spaces
spaces = []
maps = []
for i in range(num_subs):
    space_i, map_i = V.sub(i).collapse()
    spaces.append(space_i)
    maps.append(map_i)

Based on your suggestion, I defined discrete gradient as:

G = dolfinx.cpp.fem.petsc.discrete_gradient(spaces[0]._cpp_object, spaces[1]._cpp_object)

But I receive the following error:

RuntimeError: Block size is greater than 1 for V0.

What should I do for this problem?
I would like to thank you in advance for your guidance.

4

The first input space to the discrete gradient should be a scalar Lagrange space, not a vector space

1 Like