hmax on a gmsh mesh

Asked by Vladimir

Hello,

I have a problem with computing hmax on a mesh generated by gmsh. In particular, I created a unit sphere mesh in gmsh. Here is the geo-file:

lc = 0.2 ; //characteristic length

Point(1) = {0.0, 0.0, 0.0, lc} ;
Point(2) = {1.0, 0.0, 0.0, lc} ;
Point(3) = {0.0, 1.0, 0.0, lc} ;
Point(4) = {-1.0, 0.0, 0.0, lc} ;
Point(5) = {0.0, -1.0, 0.0, lc} ;
Point(6) = {0.0, 0.0, 1.0, lc} ;
Point(7) = {0.0, 0.0, -1.0, lc} ;

Circle(1) = {2, 1, 3} ;
Circle(2) = {3, 1, 4} ;
Circle(3) = {4, 1, 5} ;
Circle(4) = {5, 1, 2} ;
Circle(5) = {6, 1, 3} ;
Circle(6) = {3, 1, 7} ;
Circle(7) = {7, 1, 5} ;
Circle(8) = {5, 1, 6} ;

Line Loop(11) = {4, 1, -5, -8} ;
Line Loop(12) = {2, 3, 8, 5} ;
Line Loop(13) = {-2, -3, 7, 6} ;
Line Loop(14) = {-1, -4, -7, -6} ;

Ruled Surface(21) = {11} ;
Ruled Surface(22) = {12} ;
Ruled Surface(23) = {13} ;
Ruled Surface(24) = {14} ;

Surface Loop(25) = {21, 22, 23, 24} ;
Volume(31) = {25} ;
Physical Volume(301) = {31} ;

I used hmax to get maximum diameter over tetrahedra

from dolfin import *
import shlex, subprocess

# Create mesh
args = shlex.split('gmsh sphere.geo -3')
subprocess.call(args)
# Convert .msh to .xml
args = shlex.split('dolfin-convert sphere.msh sphere.xml')
subprocess.call(args)

# Load mesh
domain = Mesh('sphere.xml')

# Print h
print domain.hmax()

I got a number which is even greater then the diameter of unit sphere. I got the same result when I computed the parameter h on my own by reading the corresponding sphere.msh file.
There wasn't any problem with a cube mesh. Where is the problem?

I also used dolfin to make a spherical mesh

domain = UnitSphere(20)

The result for hmax is always 2. Is it correct?

Thank you in advance for your help
Vladimir Vrabel

Question information

Language:
English Edit question
Status:
Answered
For:
DOLFIN Edit question
Assignee:
No assignee Edit question
Last query:
Last reply:
Revision history for this message
Anders Logg (logg) said :
#1

On Fri, Mar 02, 2012 at 02:41:00PM -0000, Vladimir wrote:
> New question #189443 on DOLFIN:
> https://answers.launchpad.net/dolfin/+question/189443
>
> Hello,
>
> I have a problem with computing hmax on a mesh generated by gmsh. In particular, I created a unit sphere mesh in gmsh. Here is the geo-file:
>
> lc = 0.2 ; //characteristic length
>
> Point(1) = {0.0, 0.0, 0.0, lc} ;
> Point(2) = {1.0, 0.0, 0.0, lc} ;
> Point(3) = {0.0, 1.0, 0.0, lc} ;
> Point(4) = {-1.0, 0.0, 0.0, lc} ;
> Point(5) = {0.0, -1.0, 0.0, lc} ;
> Point(6) = {0.0, 0.0, 1.0, lc} ;
> Point(7) = {0.0, 0.0, -1.0, lc} ;
>
> Circle(1) = {2, 1, 3} ;
> Circle(2) = {3, 1, 4} ;
> Circle(3) = {4, 1, 5} ;
> Circle(4) = {5, 1, 2} ;
> Circle(5) = {6, 1, 3} ;
> Circle(6) = {3, 1, 7} ;
> Circle(7) = {7, 1, 5} ;
> Circle(8) = {5, 1, 6} ;
>
> Line Loop(11) = {4, 1, -5, -8} ;
> Line Loop(12) = {2, 3, 8, 5} ;
> Line Loop(13) = {-2, -3, 7, 6} ;
> Line Loop(14) = {-1, -4, -7, -6} ;
>
> Ruled Surface(21) = {11} ;
> Ruled Surface(22) = {12} ;
> Ruled Surface(23) = {13} ;
> Ruled Surface(24) = {14} ;
>
> Surface Loop(25) = {21, 22, 23, 24} ;
> Volume(31) = {25} ;
> Physical Volume(301) = {31} ;
>
> I used hmax to get maximum diameter over tetrahedra
>
> from dolfin import *
> import shlex, subprocess
>
> # Create mesh
> args = shlex.split('gmsh sphere.geo -3')
> subprocess.call(args)
> # Convert .msh to .xml
> args = shlex.split('dolfin-convert sphere.msh sphere.xml')
> subprocess.call(args)
>
> # Load mesh
> domain = Mesh('sphere.xml')
>
> # Print h
> print domain.hmax()
>
> I got a number which is even greater then the diameter of unit sphere. I got the same result when I computed the parameter h on my own by reading the corresponding sphere.msh file.
> There wasn't any problem with a cube mesh. Where is the problem?
>
> I also used dolfin to make a spherical mesh
>
> domain = UnitSphere(20)
>
> The result for hmax is always 2. Is it correct?
>
> Thank you in advance for your help
> Vladimir Vrabel

Thanks for reporting.

This is a bug in the UnitSphere class and possibly the function that
computes the diameter of tetrahedral cells. See these bug reports:

https://bugs.launchpad.net/dolfin/+bug/948361
https://bugs.launchpad.net/dolfin/+bug/948379

--
Anders

Can you help with this problem?

Provide an answer of your own, or ask Vladimir for more information if necessary.

To post a message you must log in.