The overall procedure is generating Mesh(), MeshValueCollection objects first, and reading information into them. Then generating MeshFunctionSizet accordingly to define dx, ds... for subdomains.
Step 0. Open corresponding XDMF files; e.g.,
with XDMFFile(MPI.comm_world, "mesh.xdmf") as xdmf_infile:Step 1. Create Mesh object, using
mesh = dolfin.cpp.mesh.Mesh()Step 2. Read Mesh from XDMF file, e.g.,
xdmf_infile.read(mesh)Step 3. Create MeshValueCollection object as container of subdomain information;
mvc_subdomain = dolfin.MeshValueCollection("size_t", mesh, mesh.topology().dim())Step 4. Read subdomain information,
xdmf_infile.read(mvc_subdomain, "name-to-read")The "name-to-read" here is the name given in XDMF file generated by meshio.
Step 5. Using step 3, 4 to read boundaries (open another file if necessary):
mvc_boundaries = dolfin.MeshValueCollection("size_t", mesh, mesh.topology().dim()-1)
xdmf_infile.read(mvc_boundaries, "name-to-read")Step 6. Generate MeshFunctionSizet object, this is crucial for the latter solving process:
domains = dolfin.cpp.mesh.MeshFunctionSizet(mesh, mvc_subdomain)
boundaries = dolfin.cpp.mesh.MeshFunctionSizet(mesh, mvc_boundaries)Step 7. Define dx, ds for subdomains and boundaries:
dx = Measure('dx', domain=mesh, subdomain_data=domains)
ds = Measure('ds', domain=mesh, subdomain_data=boundaries)Then we could use dx(1), dx(2)... to define problem in specific sub-domains. The index of sub-domains must be consistent with domains.array(), and same for ds. One other thing notable, the default is to have a unique subdomain with value 0. For example, if domains.array() gives [5,5,5,6,6,6], an error would be raised
*** -------------------------------------------------------------------------
*** DOLFIN encountered an error. If you are not able to resolve this issue
*** using the information listed below, you can ask for help at
***
*** fenics-support@googlegroups.com
***
*** Remember to include the error message listed below and, if possible,
*** include a *minimal* running example to reproduce the error.
***
*** -------------------------------------------------------------------------
*** Error: Unable to set given (local) rows to identity matrix.
*** Reason: some diagonal elements not preallocated (try assembler option keep_diagonal).
*** Where: This error was encountered inside PETScMatrix.cpp.
*** Process: 0
***
*** DOLFIN version: 2019.2.0.dev0
*** Git changeset: unknown
*** -------------------------------------------------------------------------
domains.set_values([0,0,0,1,1,1]).
Step 8. Give the variational problem for each sub-domain, e.g., for poisson problem:
F = inner(a0 * grad(u), grad(v)) * dx(0) + \
inner(a1 * grad(u), grad(v)) * dx(1) + \
inner(g_L, v) * ds(1) - \
inner(g_R, v) * ds(3) - \
inner(f, v) * dx(0) - \
inner(f, v) * dx(1)
Posts
