Vector/Tensor Material Properties

Hi,

I am currently setting material tensors and vectors like here:

dx = ufl.Measure("dx", domain=mesh, subdomain_data=subdomains)
ds = ufl.Measure("ds", domain=mesh, subdomain_data=boundaries)

m11 = dolfinx.Constant(mesh, 1)
m12 = dolfinx.Constant(mesh, 3)
m21 = dolfinx.Constant(mesh, 2)
m22 = dolfinx.Constant(mesh, 4)

M1 = ufl.as_matrix([[m11, m12], [m21, m22]])

# similar for M2

a = inner(M1 * grad(u), grad(v)) * dx(1) \
  + inner(M2 * grad(u), grad(v)) * dx(2)

However, this method gets long and confusing with more parameters and subdomains.

How can I define an equivalent vector- or tensor-valued function that is constant on each subdomain?

So far, I tried TensorFunctionSpace with (“DG”, 0) and setValueLocal, but I didn’t work

1 Like

Could you share the non-working code for tensorfunctionspace, so that it is possible to investigate what’s wrong with it?

I was trying along these lines:

W = dolfinx.VectorFunctionSpace(mesh, ("DG", 0), dim=3)

M = dolfinx.Function(W)
dof = W.tabulate_dof_coordinates()

for i in range(dof.shape[0]):
   if subdomains.values[i] == 1:
      M.vector.setValueLocal(i, ufl.as_vector([1, 2, 3]))
   else:
      M.vector.setValueLocal(i, ufl.as_vector([4, 5, 6]))

Consider the following example:

import dolfinx
import dolfinx.cpp as cpp
from mpi4py import MPI
import numpy as np
mesh = dolfinx.UnitSquareMesh(MPI.COMM_WORLD, 2, 2)

W = dolfinx.VectorFunctionSpace(mesh, ("DG", 0), dim=3)
Wsubs = [W.sub(i) for i in range(3)]
M = dolfinx.Function(W)

topology = mesh.topology
map_c = topology.index_map(topology.dim)
num_cells = map_c.size_local + map_c.num_ghosts

entities = np.array([0,2,4], dtype=np.int32)
values = np.array([1,2,1], dtype=np.int32)
subdomains = dolfinx.MeshTags(mesh, topology.dim, entities, values)

case1 = [1,2,3]
case2 = [4,5,6]

for i in range(num_cells):
    for j in range(3):
        wj_dofs = Wsubs[j].dofmap.cell_dofs(i)
        for dof in wj_dofs:
            # Check if i is in the MeshTag
            if i in subdomains.indices:
                # Find where in the MeshTag the cell is
                loc = np.flatnonzero(i==subdomains.indices)[0]
                if subdomains.values[loc] == 1:
                    M.vector.setValueLocal(dof, case1[j])
                    continue
            M.vector.setValueLocal(dof, case2[j])
print(M.vector.array)

Here, I have used the sub dof-maps per dimension, to set the values for cell 0 and 4 to 1,2,3 for each subspace, while all others are set to 4,5,6.