///
/// This file is part of Rheolef.
///
/// Copyright (C) 2000-2009 Pierre Saramito <Pierre.Saramito@imag.fr>
///
/// Rheolef is free software; you can redistribute it and/or modify
/// it under the terms of the GNU General Public License as published by
/// the Free Software Foundation; either version 2 of the License, or
/// (at your option) any later version.
///
/// Rheolef is distributed in the hope that it will be useful,
/// but WITHOUT ANY WARRANTY; without even the implied warranty of
/// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
/// GNU General Public License for more details.
///
/// You should have received a copy of the GNU General Public License
/// along with Rheolef; if not, write to the Free Software
/// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
/// 
/// =========================================================================
// tested with gmsh version 2.8.3
R = 1;   // large radius
r = 0.6; // small radius
n = 10;  // density of discretisation
h = 1.0/n;
n1 = Ceil(2*Pi*r/h);
n2 = Ceil(2*Pi*R/(3*h));
Printf("n = %g, n1 = %g, n2 = %g", n, n1, n2);
Point(1) = {  R, 0, 0, h};
Point(2) = {R+r, 0, 0, h};
Point(3) = {  R, 0, r, h};
Point(4) = {  R, 0,-r, h};
Point(5) = {R-r, 0, 0, h};
Circle(1) = {2,1,3};
Circle(2) = {3,1,5};
Circle(3) = {5,1,4};
Circle(4) = {4,1,2};
Line Loop(5) = {1,2,3,4};
Plane Surface(6) = {5};
Geometry.AutoCoherence = 0;
nn1 = n;
third1[]=Extrude {{0,0,1}, {0,0,0}, 2*Pi/3}{Surface{6}; Layers{nn1}; };
nn2 = n;
third2[]=Extrude {{0,0,1}, {0,0,0}, 2*Pi/3}{Surface{third1[0]}; Layers{nn2}; };
nn3 = n;
third3[]=Extrude {{0,0,1}, {0,0,0}, 2*Pi/3}{Surface{third2[0]}; Layers{nn3}; };
Delete {
  Surface{6};
  Surface{third1[0]};
  Surface{third2[0]};
}
Geometry.AutoCoherence = 1;
Coherence;
Physical Surface("boundary") = {
 -third1[2],-third1[3],-third1[4],-third1[5],
 -third2[2],-third2[3],-third2[4],-third2[5],
 -third3[2],-third3[3],-third3[4],-third3[5]
};
