import salome, smesh, geompy

quad  = geompy.MakeSketcher("Sketcher:F 5 0:TT 30 0:TT 0 30:TT 0 5:WF")
OZ    = geompy.MakeVectorDXDYDZ(0,0,1)
quad4 = geompy.MultiRotate1D( quad, OZ, 4 )

quaIn = geompy.MakeSketcher("Sketcher:F 5 0:TT 0 5:TT -5 0:TT 0 -5:WF")

circ  = geompy.MakeCircleR( 15 )

edge  = geompy.GetEdgeNearPoint( quad4, geompy.MakeVertex( 100, 100, 0 ))
rect  = geompy.MakePrismDXDYDZ( edge, 100, 100, 0 )

all2D = geompy.Partition( [quad4, quaIn, circ, rect ] )

state = geompy.GEOM.ST_IN
circSubs = geompy.GetShapesOnCylinder( all2D, geompy.ShapeType["FACE"], OZ, 15, state )
circSubs = geompy.MakeCompound( circSubs )

cyl = geompy.MakePrismDXDYDZ( circSubs, 0,0,-100)
box = geompy.MakePrismDXDYDZ( all2D,    0,0,100)

shape = geompy.MakeCompound( [cyl, box ] )
geompy.addToStudy( shape, "shape" )


mesh = smesh.Mesh( shape )
mesh.Segment().LocalLength( 2 )
mesh.Quadrangle()
mesh.Hexahedron()
subAlgo = mesh.Segment(circSubs)
subAlgo.NumberOfSegments(8)
subAlgo.Propagation()

mesh.Compute()