DOLFIN

Read /Extract Gmsh subdomain

Asked by Paul Scott

I used Gmsh to create a small sphere inside a big sphere. I labeled volumes (inner sphere, and outer sphere) in Gmsh. I can convert .msh to .xml with no problem. But when using this file in FEniCS (I use python), how to I get this formation out and use it to define subdomain (which is just the inner sphere).

More information: I am trying to solve Poisson equation and I am modify this demo code:
http://fenicsproject.org/documentation/dolfin/dev/python/demo/pde/subdomains-poisson/python/documentation.html

Question information

Language:
English Edit question
Status:
Solved
For:
DOLFIN Edit question
Assignee:
No assignee Edit question
Solved by:
M Waste
Solved:
Last query:
Last reply:
Revision history for this message
M Waste (michael-waste) said : 		#1

I'm not sure if I get you right.

but with the following it was no problem to define boundary conditions on the corresponding gmsh facets.

In the "geo file" you have to rename manually the subdomains with 0 and 1. Also the volumes need integer numbers beginning from 0.

then you can reference to the corresponding boundaries and subdomains in the python file.

mesh = Mesh("yourmeshfile.xml")
subdomains = MeshFunction("size_t", mesh, "yourmeshfile_physical_region.xml")
boundaries = MeshFunction("size_t", mesh, "yourmeshfile_facet_region.xml")

bcs = [DirichletBC(V, 5.0, boundaries, 1),# of course with your boundary
       DirichletBC(V, 0.0, boundaries, 0)]

Michael

Revision history for this message
Paul Scott (paul-scott-rocky) said : 		#2

Thank you, Michael for helping. I tried your suggestion but it still does not work. I will explain more so that you can understand my problem better.
So I have the mesh file from Gmsh " sphereinsphere.msh" and in FEniCS terminal, I convert it to .xml file by typing:

" dolfin-convert sphereinsphere.msh sphereinsphere.xml " and it has no problem at this process.

Now in my Python file, with your suggestions, I put this in my code:
"
mesh = Mesh("sphereinsphere.xml")
subdomains = MeshFunction("size_t", mesh, "sphereinsphere_physical_region.xml")
boundaries = MeshFunction("size_t", mesh, "sphereinsphere_facet_region.xml")

"
After running the Python file, it gives me this error:
"
Traceback (most recent call last):
  File "poisson3d_sphere.py", line 55, in <module>
    subdomains = MeshFunction("size_t", mesh, "sphere_physical_region.xml")
  File "/Applications/FEniCS.app/Contents/Resources/lib/python2.7/site-packages/dolfin/cpp/mesh.py", line 6418, in __new__
    return MeshFunctionSizet(*args)
  File "/Applications/FEniCS.app/Contents/Resources/lib/python2.7/site-packages/dolfin/cpp/mesh.py", line 5193, in __init__
    _mesh.MeshFunctionSizet_swiginit(self,_mesh.new_MeshFunctionSizet(*args))
RuntimeError:

*** -------------------------------------------------------------------------
*** DOLFIN encountered an error. If you are not able to resolve this issue
*** using the information listed below, you can ask for help at
***
*** https://answers.launchpad.net/dolfin
***
*** Remember to include the error message listed below and, if possible,
*** include a *minimal* running example to reproduce the error.
***
*** -------------------------------------------------------------------------
*** Error: Unable to read data from XML file.
*** Reason: Unable to open file "sphere_physical_region.xml".
*** Where: This error was encountered inside XMLFile.cpp.
*** Process: 0
*** -------------------------------------------------------------------------

"

Do you know what might cause the problem?

Revision history for this message
Best M Waste (michael-waste) said : 		#3

hi,

The error message refer to a file "sphere_physical_region.xml" but your python file call the file "sphereinsphere_physical_region.xml".

Either your xml file is not in your working directory or you made a typing error.

The geo file needs "physical surfaces", "physical volumes" and "physical lines" otherwise no "physical_region.xml" and "_facet_region.xml" were generated with dolfin convert.

Boundary labelling: physical lines, physical surfaces resp. Sart labelling with 1all other unused boundaries are then set to zero (natural boundaries).
plot your boundaries in order to see if all boundaries are correct.

HTH

Revision history for this message
Paul Scott (paul-scott-rocky) said : 		#4

Thanks M Waste, that solved my question.

Revision history for this message
Paul Scott (paul-scott-rocky) said : 		#5

In Gmsh .geo file, I relabeled/added physical lines,surfaces, volumes. It's all working now. Thanks again.