[Gmsh] Volume{} creates duplicate Surfaces: gmsh w/o OpenCascade

Paul Fulmek fulmek at tuwien.ac.at
Fri Jan 19 16:23:58 CET 2018


Dear gmsh-list,

I noticed a different or strange behaviour of gmsh-3.0.1/3.0.6 w/o 
SetFactory("OpenCASCADE") when creating a Volume{}.
Usually I rely upon my list of surfaces, which I use in the Surface 
Loop{} command to create a Volume, accordingly I declare the Physical-ID 
Volumes/Surfaces from my private book-keeping of surfaces and volumes. 
If I do use SetFactory("OpenCASCADE") in my script, however, this 
assumption fails, at least in my simple test-case. The Volume-command 
creates new Surfaces, ignoring my original surface-IDs, which leads to 
physical volumes which are not bounded by my physical surfaces, leading 
to unexpected behaviour in the subsequent GetDP simulations.

I have attached a simple test-case geo-script which shows this 
behaviour, this script contains some more comments on what has been tried:

If using gmsh-built-in everything is OK, as expected

If using openCASCADE a strange thing happens:
       Volume(P_Vol) = { surface_loop };
creates some duplicate surfaces from the complete and correct 
surface-list!  I.e. my original surface-loop does no longer describe the 
boundaries of the new volume!
gmsh-3.0.6 shows this behavior always
gmsh-3.0.1 shows this only if my original surface-ID-list has been 
modified in the process!

Whats wrong with my script/assumptions?

Thank you very much,

kind regards

Paul


-------------- next part --------------
/***************************
 * model_1.geo
 * dummy-test-file to show a strange gmsh behaviour if I modify a volume-surface of a simple block:
 * 
  * if using gmsh-built-in everything is OK, as expected
 * 
 * if using openCASCADE a strange thing happens:
 *	 Volume(P_Vol) = { surface_list };
 * creates duplicate surfaces from surface-list!
 * i.e. my original surface-list does no longer describe the boundaries of the volume!
 * I must not use my stored surface-list-IDs as physical-boundary of the volume!
 * gmsh-3.0.6 shows this behavior always
 * gmsh-3.0.1 shows this only if the P_Surf[]-list has been modified!
 * 
 * The output-comments of this script are designed for processing with gmsg-3.0.1!
***************************/

use_openCascade = 1;
choose_new_ID = 0;

/***************************
 * test-case w/o SetFactory("OpenCASCADE");
 
1)  create lines, create surface, create volume by extrusion, store volume-boundaries in my P_surf[]-list
2)  delete the volume, delete one surface with the ID P_Surf[0]
3A) create a new surface with the ID P_surf[0]		-> behaves OK, as expected in gmsh-3.0.1
3B) create a new surface with a NEW ID P_surf[0]
4)  create the surface-loop from P_surf[]
5)  create the volume again							-> behaves STRANGE, unexpected! 
6)  compare surface-loop P_surf[] to the queried volume-boundaries

***********************************************************
*** OK:  choose_new_ID = 0  *******************************

OK: all as expected, if I simply re-use my ID P_surf[0] for creating the new surface and the new volume.
	The boundary of the new volume exactly reflects my surface-list.
	I do receive 1 volume and 6 surfaces, as expected my 6 items of P_surf[]

*** OUTPUT:
*** model_1.geo: test case geometry design ***
*** OK: keep my private surface-list       ***
***** create simple Plane Surface *****
***** extrude surface to a volume and store the boundary-surface-list *****
***** resulting P_Vol, P_Surf[] *****
  P_Vol #1: 1, 
  P_Surf #6: 6, 11, 7, 8, 9, 10, 
  bnd P_Surf [0] #4: 1, 2, 3, 4, 
  bnd P_Surf [1] #4: 7, 9, 11, 12, 
  bnd P_Surf [2] #4: 5, 7, -6, -1, 
  bnd P_Surf [3] #4: 6, 9, -8, -2, 
  bnd P_Surf [4] #4: 8, 11, -10, -3, 
  bnd P_Surf [5] #4: 10, 12, -5, -4, 
***** query volume boundary *****
  P_Vbnd #6: 7, 8, 9, 10, 6, 11, 
  bnd P_Vbnd [0] #4: 5, 7, -6, -1, 
  bnd P_Vbnd [1] #4: 6, 9, -8, -2, 
  bnd P_Vbnd [2] #4: 8, 11, -10, -3, 
  bnd P_Vbnd [3] #4: 10, 12, -5, -4, 
  bnd P_Vbnd [4] #4: 1, 2, 3, 4, 
  bnd P_Vbnd [5] #4: 7, 9, 11, 12, 
***** delete P_Vol 1
***** delete P_Surf[0] 6
***** keep the original  P_Surf-list *****
***** create new surface from the original line-loop *****
  P_S0bd #4: 1, 2, 3, 4, 
***** create new surface-loop and volume *****
***** resulting P_Vol, P_Surf[] *****
  P_Vol #1: 15, 
  P_Surf #6: 6, 11, 7, 8, 9, 10, 
  bnd P_Surf [0] #4: 1, 2, 3, 4, 
  bnd P_Surf [1] #4: 7, 9, 11, 12, 
  bnd P_Surf [2] #4: 5, 7, -6, -1, 
  bnd P_Surf [3] #4: 6, 9, -8, -2, 
  bnd P_Surf [4] #4: 8, 11, -10, -3, 
  bnd P_Surf [5] #4: 10, 12, -5, -4, 
***** query volume boundary *****
  P_Vbnd #6: 6, 7, 11, 8, 9, 10, 
  bnd P_Vbnd [0] #4: 1, 2, 3, 4, 
  bnd P_Vbnd [1] #4: 5, 7, -6, -1, 
  bnd P_Vbnd [2] #4: 7, 9, 11, 12, 
  bnd P_Vbnd [3] #4: 6, 9, -8, -2, 
  bnd P_Vbnd [4] #4: 8, 11, -10, -3, 
  bnd P_Vbnd [5] #4: 10, 12, -5, -4, 

***********************************************************
*** STRANGE:  choose_new_ID = 1  **************************

STRANGE: after deleting the Volume P_Vol and the surface P_Surf[0] 
	a new ID for my new surface P_surf[0] is chosen, e.g. by P_Surf[0] = news.
	Then create the new Surface P_surf[0], create a surface-loop from my modified P_Surf[]-list,
	and then create the new Volume.
	Obviously the Volume{}-command creates new surfaces, which are exact duplicates of my original
	surfaces declared in the surface-loop?
	What is the reason for the discrepancy between the Boundary of the volume and my Surface-list,
	which is the root for the Volume creation!
	I do receive 1 volume and 10 surfaces!
	The Volume-command has created 4 new surfaces, duplicating 4 of my surface-loop-items!
	Why?

*** OUTPUT:
*** model_1.geo: test case geometry design ***
*** STRANGE: modify my private surface-list       ***
***** create simple Plane Surface *****
***** extrude surface to a volume and store the boundary-surface-list *****
***** resulting P_Vol, P_Surf[] *****
  P_Vol #1: 1, 
  P_Surf #6: 6, 11, 7, 8, 9, 10, 
  bnd P_Surf [0] #4: 1, 2, 3, 4, 
  bnd P_Surf [1] #4: 7, 9, 11, 12, 
  bnd P_Surf [2] #4: 5, 7, -6, -1, 
  bnd P_Surf [3] #4: 6, 9, -8, -2, 
  bnd P_Surf [4] #4: 8, 11, -10, -3, 
  bnd P_Surf [5] #4: 10, 12, -5, -4, 
***** query volume boundary *****
  P_Vbnd #6: 7, 8, 9, 10, 6, 11, 
  bnd P_Vbnd [0] #4: 5, 7, -6, -1, 
  bnd P_Vbnd [1] #4: 6, 9, -8, -2, 
  bnd P_Vbnd [2] #4: 8, 11, -10, -3, 
  bnd P_Vbnd [3] #4: 10, 12, -5, -4, 
  bnd P_Vbnd [4] #4: 1, 2, 3, 4, 
  bnd P_Vbnd [5] #4: 7, 9, 11, 12, 
***** delete P_Vol 1
***** delete P_Surf[0] 6
***** choose a new ID for P_Surf[0] *****
  P_Surf #6: 13, 11, 7, 8, 9, 10, 
***** create new surface from the original line-loop *****
  P_S0bd #4: 1, 2, 3, 4, 
***** create new surface-loop and volume *****
***** resulting P_Vol, P_Surf[] *****
  P_Vol #1: 18, 
  P_Surf #6: 13, 11, 7, 8, 9, 10, 
  bnd P_Surf [0] #4: 1, 2, 3, 4, 
  bnd P_Surf [1] #4: 7, 9, 11, 12, 
  bnd P_Surf [2] #4: 5, 7, -6, -1, 
  bnd P_Surf [3] #4: 6, 9, -8, -2, 
  bnd P_Surf [4] #4: 8, 11, -10, -3, 
  bnd P_Surf [5] #4: 10, 12, -5, -4, 
***** query volume boundary *****
  P_Vbnd #6: 13, 14, 11, 15, 16, 17, 
  bnd P_Vbnd [0] #4: 1, 2, 3, 4, 
  bnd P_Vbnd [1] #4: 5, 7, -6, -1, 
  bnd P_Vbnd [2] #4: 7, 9, 11, 12, 
  bnd P_Vbnd [3] #4: 6, 9, -8, -2, 
  bnd P_Vbnd [4] #4: 8, 11, -10, -3, 
  bnd P_Vbnd [5] #4: 10, 12, -5, -4, 

*********************************/

Printf("*** model_1.geo: test case geometry design ***");

If( use_openCascade )
	SetFactory("OpenCASCADE");
	Printf("*** gmsh OpenCASCADE:");
Else
	Printf("*** gmsh built-in:");
EndIf

If( choose_new_ID )
	If( use_openCascade )
		Printf("*** STRANGE: modify my private surface-list       ***");
	Else
		Printf("*** OK: keep my private surface-list       ***");
	EndIf
Else	
	Printf("*** OK: keep my private surface-list       ***");
EndIf

Pxsize = 1.0;	Pysize = 1.0;	Pzsize = 0.3;	Pmesh = 0.1;
Cxsize = 0.5;	Cysize = 0.5;	Czsize = 0.2;	Cmesh = 0.1;

/**** print my lists... ****/
Macro print_list
/* print contents of variable "VName" */
    tmpList[] = StringToName[ VName ][];
    N = #tmpList[] - 1;
    tmpString = StrCat["  ",StrCat[ VName, Sprintf(" #%g: ",N+1) ]];
    For t In {0:N}
        tmpString = StrCat[ tmpString, Sprintf("%g, ",tmpList[t])];
    EndFor
    Printf( tmpString);
Return
/**** das wars mit: Macro print my lists... ****/

/**** print the surface boundaries ****/
Macro print_surf_bnd
/* print surface-boundaries of surface-list "VName" */
    tmpSurf[] = StringToName[ VName ][];
    Ns = #tmpSurf[] - 1;
	For ts In {0:Ns}
		tmpList = Boundary{ Surface{ tmpSurf[ts] }; };
		N = #tmpList[] - 1;
		tmpString = StrCat["  bnd ",StrCat[ VName, Sprintf(" [%g] #%g: ", ts, N+1 ) ]];
		For t In {0:N}
			tmpString = StrCat[ tmpString, Sprintf("%g, ",tmpList[t])];
		EndFor
		Printf( tmpString);
	EndFor
Return
/**** das wars mit: Macro print my lists... ****/

// PLATE corners: pp[]
xx = Pxsize/2.0;	yy = Pysize/2.0;	zz = 0;	ms = Pmesh;
pp[0] = newp;	Point(pp[0]) = { xx, yy, 0, ms };
pp[1] = newp;	Point(pp[1]) = {-xx, yy, 0, ms };
pp[2] = newp;	Point(pp[2]) = {-xx,-yy, 0, ms };
pp[3] = newp;	Point(pp[3]) = { xx,-yy, 0, ms };

// we have all points/lists!
// create lines and loops: PLATE
pl[0] = newl;	Line(pl[0]) = { pp[0], pp[1] };
pl[1] = newl;	Line(pl[1]) = { pp[1], pp[2] };
pl[2] = newl;	Line(pl[2]) = { pp[2], pp[3] };
pl[3] = newl;	Line(pl[3]) = { pp[3], pp[0] };
pll = newll;	Line Loop(pll) = { pl[] };

Printf("***** create simple Plane Surface *****");
P_Surf[0] = news;	Plane Surface(P_Surf[0]) = { pll };	// bottom
Printf("***** extrude surface to a volume and store the boundary-surface-list *****");
tmp[] = Extrude{0,0,-Pzsize}{ Surface{P_Surf[0]}; };
P_Vol = tmp[1];
P_Surf[1] = tmp[0];	// top
P_Surf[2] = tmp[2];
P_Surf[3] = tmp[3];
P_Surf[4] = tmp[4];
P_Surf[5] = tmp[5];
Printf("***** resulting P_Vol, P_Surf[] *****");
VName =  Str["P_Vol"];	Call print_list;
VName =  Str["P_Surf"];	Call print_surf_bnd;	Call print_list;
Printf("***** query volume boundary *****");
P_Vbnd[] = Boundary{ Volume{ P_Vol }; };
VName = Str["P_Vbnd"];	Call print_surf_bnd;	Call print_list;

// modify the top-surface O_Surf[0]
Delete{ Volume{ P_Vol }; }
Printf("***** delete P_Vol %g",P_Vol);
Delete{ Surface{ P_Surf[0] }; }
Printf("***** delete P_Surf[0] %g",P_Surf[0]);

// NOW CHOOSE A NEW SURFACE-ID:
If( choose_new_ID )
	Printf("***** choose a new ID for P_Surf[0] *****");
	P_Surf[0] = news;		VName =  Str["P_Surf"];	Call print_list;
Else
	Printf("***** keep the original  P_Surf-list *****");
EndIf

Printf("***** create new surface from the original line-loop *****");
Plane Surface(P_Surf[0]) = { pll };
P_S0bd[] = Boundary{ Surface{ P_Surf[0] }; };
VName = Str["P_S0bd"];		Call print_list;

Printf("***** create new surface-loop and volume *****");
tmp = newsl;	Surface Loop(tmp) = { P_Surf[] };
P_Vol = newv;	Volume(P_Vol) = { tmp };
Printf("***** resulting P_Vol, P_Surf[] *****");
VName =  Str["P_Vol"];	Call print_list;
VName =  Str["P_Surf"];	Call print_surf_bnd;	Call print_list;
Printf("***** query volume boundary *****");
P_Vbnd[] = Boundary{ Volume{ P_Vol }; };
VName = Str["P_Vbnd"];	Call print_surf_bnd;	Call print_list;



More information about the gmsh mailing list