I have a problem where I am considering multiple domains. I am attempting to formulate the problem in one of the subdomains by conveniently pulling mesh, subdomain and boundary data from .xml files derived from my .msh. I define the submesh and subdomain without a problem, but when trying to define the boundaries of my submesh, I get the following error:
IndexError: boolean index did not match indexed array along dimension 0; dimension is 27115 but corresponding boolean dimension is 29582
I believe that this comes from the inconsistency in the codimensions between facets in the overall mesh, and what would be edges in the subdomain. I can simply define the boundaries geometrically, although this would only be a temporary solution that would need to be fixed later on. The easiest approach may be to define separate meshes, seeing as I’ll need to interpolate regardless. Although, if there is an easy solution to this then I’d rather keep the single mesh.
Any thoughts would be appreciated!
Update: I realize that the two boundary arrays are of different length. Is there any straightforward approach to extracting the relevant information (the boundary data from the larger array with only considering one subdomain)? Thank you.
from fenics import *
from mshr import *
import numpy as np
from msh2xdmf import import_mesh
#from fenicstools import interpolate_nonmatching_mesh
#Create mesh, domains and boundaries
mesh, bound, subd, association_table= import_mesh(
File(“Boundaries.xml”) << bound
fix_subdomains = MeshFunction(“size_t”, mesh, mesh.topology().dim())
fix_subdomains.array()[subd.array()==11] = 1
submesh = SubMesh(mesh, fix_subdomains, 1)
boundaries = MeshFunction(“size_t”, mesh, “Boundaries.xml”)
submesh_boundaries = MeshFunction(‘size_t’, submesh, mesh.topology().dim()-1)
submesh_boundaries.array()[boundaries.array()==6] = 1