Plotting the boundary function of a 2d mesh

mesh
1

Below is a simple example that creates a mesh, and marks different parts of the boundary with a boundary function. I would like to be able to plot this boundary function to check whether the markings are correct, however, when I try to call plot(boundary_function) (which is how I would plot a mesh function) I get the warning,

*** Warning: Matplotlib plotting backend does not support mesh function of dim 1. Continuing without plotting...

I am running the latest stable version of FEniCS available from Docker on a mac. Does anyone know a good way of plotting a boundary function so that the marked parts of the boundary show up as different colours?

from dolfin import *
import matplotlib.pyplot as plt

# create mesh
mesh = UnitSquareMesh(15,15)

# mark mesh with boundary function
class left(SubDomain):
    def inside(self, x, on_boundary):
        return on_boundary and abs(x[0])<DOLFIN_EPS

left = left()
boundary_function = MeshFunction("size_t", mesh, mesh.geometric_dimension()-1)
boundary_function.set_all(0)
left.mark(boundary_function, 1)

plot(boundary_function)
plt.savefig("tmp.png")
1 Like
2

You can write the function to a file and then inspect it in paraview for example.

boundary_file = df.File("boundaries.pvd")
boundary_file << boundary_function
1 Like
3

Hi Lukas:
I tried your solution. It runs without error and writes the .pvd files (boundaries.pvd and boundaries000000.vtu).
When I open the file boundaries.pvd in Paraview I get the following error:
Algorithm vtkXMLUnstructuredGridReader returned failure for request: vtkInformation
and it does not show anything.

4

Can you please provide a minimal example?

5

Yes, this is the minimal working example based on demo_stokes-taylor-hood.py

import matplotlib.pyplot as plt
from dolfin import *
# Load mesh and subdomains
mesh = Mesh("dolfin_fine.xml.gz")
sub_domains = MeshFunction("size_t", mesh, "dolfin_fine_subdomains.xml.gz")
# this runs but does not work when reading in paraview: Algorithm vtkXMLUnstructuredGridReader returned failure for request: vtkInformation
boundary_file = File("boundaries.pvd")
boundary_file << sub_domains
6

Does the plotting from the demo work for you?
I get an error if I try to use plotting because it is not supported for the subdomains.
I get the same error as you if I try to save as pvd.
I managed to save it as XDMF file but if I look at it there is no information contained about the subdomains…

boundary_file = XDMFFile(“boundaries.xdmf”)
boundary_file.write(sub_domains)

7

No, it doesn’t. That’s the reason I tried your solution. It gives the same error as alastair at the beginning of the thread:

*** Warning: Matplotlib plotting backend does not support mesh function of dim 1. Continuing without plotting...

8

It might be because sub_domains.dim() is 1.

9

yes
Your solution prior to editing seemed promising but importing vtkplotter.dolfin gives me error.