Domain with multiple different orthotropic materials in dolphinx

mesh
3

With the help of the following to post assemble matrix of subdomain I tried to assemble the stiffness matrix over each subdomain:

from dolfinx.fem import VectorFunctionSpace
import ufl
import dolfinx

V = VectorFunctionSpace(mesh, ('CG', 2)) # CG-Lagrange element, 2-degree of freedom (2nd order element)

# Define variational problem
u = ufl.TrialFunction(V)
v = ufl.TestFunction(V)

material_tags = np.unique(cell_markers.values)
stiffnesses = []
for tag in material_tags:
    cells = cell_markers.find(tag)
    if tag < len(stack):
        material = stack[tag]
        dx_domain = ufl.Measure('dx', domain=mesh, subdomain_data=cell_markers.find(tag), subdomain_id=tag)
        stiffness = ufl.as_matrix(material["stiff_tens"])
        a = ufl.inner(sigma(u, stiffness), epsilon(v))* dx_domain
        stiffnesses.append(a)

This worked (kind of) until I wanted to solve my problem:

from dolfinx.fem.petsc import LinearProblem


mat = stiffnesses
result = mat[0] + mat[1] + mat[2] # trying to add the different domains back together

problem = LinearProblem(result, L, bcs=[bc], petsc_options={"ksp_type": "preonly", "pc_type": "lu"})
uh = problem.solve()

results in:

---------------------------------------------------------------------------
AssertionError                            Traceback (most recent call last)
Cell In[18], line 5
      3 mat = stiffnesses
      4 result = mat[0] + mat[1] + mat[2]
----> 5 problem = LinearProblem(result, L, bcs=[bc], petsc_options={"ksp_type": "preonly", "pc_type": "lu"})
      6 uh = problem.solve()

File /usr/local/dolfinx-complex/lib/python3.10/dist-packages/dolfinx/fem/petsc.py:566, in LinearProblem.__init__(self, a, L, bcs, u, petsc_options, form_compiler_options, jit_options)
    538 def __init__(self, a: ufl.Form, L: ufl.Form, bcs: typing.List[DirichletBCMetaClass] = [],
    539              u: typing.Optional[_Function] = None, petsc_options={}, form_compiler_options={}, jit_options={}):
    540     """Initialize solver for a linear variational problem.
    541 
    542     Args:
   (...)
    564                                                "pc_factor_mat_solver_type": "mumps"})
    565     """
--> 566     self._a = _create_form(a, form_compiler_options=form_compiler_options, jit_options=jit_options)
    567     self._A = create_matrix(self._a)
    569     self._L = _create_form(L, form_compiler_options=form_compiler_options, jit_options=jit_options)

File /usr/local/dolfinx-complex/lib/python3.10/dist-packages/dolfinx/fem/forms.py:176, in form(form, dtype, form_compiler_options, jit_options)
    173         return list(map(lambda sub_form: _create_form(sub_form), form))
    174     return form
--> 176 return _create_form(form)
...
--> 138     assert all([id(d) == id(data[0]) for d in data])
    139     subdomains[integral_type] = data[0]
    141 mesh = domain.ufl_cargo()

AssertionError:

Can someone, PLEASE, tell me how to merge the matrixies of the different domains back together? OR is there a better way?