1 # -*- coding: latin-1 -*-
2 # Copyright (C) 2009-2015 CEA/DEN, EDF R&D
4 # This library is free software; you can redistribute it and/or
5 # modify it under the terms of the GNU Lesser General Public
6 # License as published by the Free Software Foundation; either
7 # version 2.1 of the License, or (at your option) any later version.
9 # This library is distributed in the hope that it will be useful,
10 # but WITHOUT ANY WARRANTY; without even the implied warranty of
11 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 # Lesser General Public License for more details.
14 # You should have received a copy of the GNU Lesser General Public
15 # License along with this library; if not, write to the Free Software
16 # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
18 # See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
21 # Francis KLOSS : 2012 : CEA-Saclay, DEN, DM2S, SFME, LGLS, F-91191 Gif-sur-Yvette, France
22 # ========================================================================================
27 geompy = hexablock.geompy
29 # Charger la géométrie
30 # ====================
34 geometrie = geompy.ImportBREP(nom+".brep")
36 # Sélectionner des sous-parties de la géométrie
37 # ---------------------------------------------
39 sommets = geompy.SubShapeAllSortedCentres(geometrie, geompy.ShapeType["VERTEX"])
40 aretes = geompy.SubShapeAllSortedCentres(geometrie, geompy.ShapeType["EDGE" ])
41 faces = geompy.SubShapeAllSortedCentres(geometrie, geompy.ShapeType["FACE" ])
43 petit_arri_v = sommets[0]
44 petit_arri_x, petit_arri_y, petit_arri_z = geompy.PointCoordinates(petit_arri_v)
45 petit_avan_x, petit_avan_y, petit_avan_z = geompy.PointCoordinates(sommets[1])
47 moyen_arri_v = sommets[3]
48 moyen_arri_x, moyen_arri_y, moyen_arri_z = geompy.PointCoordinates(moyen_arri_v)
49 moyen_avan_x, moyen_avan_y, moyen_avan_z = geompy.PointCoordinates(sommets[4])
51 petit0, petit_x, petit_y, petit_z, petit1, petit2, petit3, petit_rayon = geompy.KindOfShape(aretes[0])
53 moyen0, moyen_x, moyen_y, moyen_z, moyen1, moyen2, moyen3, moyen_rayon = geompy.KindOfShape(aretes[4])
55 grand0, grand_xb, grand_yb, grand_zb, grand1, grand2, grand3, grand_rayon = geompy.KindOfShape(aretes[7])
56 grand4, grand_xh, grand_yh, grand_zh, grand5, grand6, grand7, grand8 = geompy.KindOfShape(aretes[8])
58 grand_generat = aretes[11]
60 arc0, arc_x, arc_y, arc_z, arc_dx, arc_dy, arc_dz, arc_rayon, arc_x1, arc_y1, arc_z1, arc_x2, arc_y2, arc_z2 = geompy.KindOfShape(aretes[12])
62 face_cylindre = faces[5]
64 # Construire des géométries intermédiaires
65 # ----------------------------------------
67 arc_sc = geompy.MakeVertex(arc_x , arc_y , arc_z )
68 arc_s1 = geompy.MakeVertex(arc_x1, arc_y1, arc_z1)
69 arc_s2 = geompy.MakeVertex(arc_x2, arc_y2, arc_z2)
71 arc_a1 = geompy.MakeEdge(arc_sc, arc_s1)
72 arc_a2 = geompy.MakeEdge(arc_sc, arc_s2)
74 # Construire le modèle de bloc
75 # ============================
77 doc = hexablock.addDocument(nom)
79 # Construire le grand cylindre en 2 parties
80 # -----------------------------------------
82 grand_base = doc.addVertex(grand_xb, grand_yb, grand_zb)
83 grand_oppo = doc.addVertex(grand_xh, grand_yh, grand_zh)
85 grand_dir_a = doc.addVectorVertices(grand_base, grand_oppo)
86 grand_dir_b = doc.addVectorVertices(grand_oppo, grand_base)
88 grand_hauteur_a = grand_dir_a.getNorm()*0.20
89 grand_hauteur_b = grand_dir_b.getNorm()*0.60
91 #### grand_cylindre_a = doc.addCylinder(grand_base, grand_dir_a, grand_rayon, grand_hauteur_a)
92 #### grand_cylindre_b = doc.addCylinder(grand_oppo, grand_dir_b, grand_rayon, grand_hauteur_b)
94 # Construire le moyen cylindre
95 # ----------------------------
97 moyen_base = doc.addVertex(moyen_x, moyen_y, moyen_z)
99 moyen_arri = doc.addVertex(moyen_arri_x, moyen_arri_y, moyen_arri_z)
100 moyen_avan = doc.addVertex(moyen_avan_x, moyen_avan_y, moyen_avan_z)
102 moyen_dir = doc.addVectorVertices(moyen_arri, moyen_avan)
104 moyen_hauteur = geompy.MinDistance(moyen_arri_v, grand_generat) - grand_rayon
107 # Construire le petit cylindre
108 # ----------------------------
110 petit_base = doc.addVertex(petit_x, petit_y, petit_z)
112 petit_arri = doc.addVertex(petit_arri_x, petit_arri_y, petit_arri_z)
113 petit_avan = doc.addVertex(petit_avan_x, petit_avan_y, petit_avan_z)
115 petit_dir = doc.addVectorVertices(petit_arri, petit_avan)
117 petit_hauteur = geompy.MinDistance(petit_arri_v, grand_generat)
120 # Construire les 2 cylindres qui s'intersectent
121 # ---------------------------------------------
123 moyen_rayon = moyen_rayon*0.7
124 cylindres_gm = doc.makeCylinders (grand_base, grand_dir_a, grand_rayon, grand_hauteur_a,
125 moyen_base, moyen_dir, moyen_rayon, moyen_hauteur)
126 cylindres_gp = doc.makeCylinders (grand_oppo, grand_dir_b, grand_rayon, grand_hauteur_b,
127 petit_base, petit_dir, petit_rayon, petit_hauteur)
129 doc.saveVtk ("tuyauterie0.vtk")
131 print "grand_rayon = ", grand_rayon
132 print "moyen_rayon = ", moyen_rayon
133 print "petit_rayon = ", petit_rayon
136 # Joindre les 2 croix
137 # -------------------
143 qstart = cylindres_gm.getQuadIJ (CYL_BIG, 0,0,CYL_KMAX)
145 qdest = cylindres_gp.getQuadIJ (CYL_BIG, 0,0,CYL_KMAX)
147 va1 = qstart.getVertex (0)
148 va2 = qstart.getVertex (1)
149 vb1 = qdest .nearestVertex (va1)
150 vb2 = qdest .nearestVertex (va2)
154 for nj in range (CYL_JMAX) :
155 quad = cylindres_gm.getQuadIJ (CYL_BIG, 1,nj, CYL_KMAX)
156 gm_quads.append (quad)
158 prisme = doc.joinQuadsUni (gm_quads, qdest, va1,vb1,va2,vb2, 1)
161 # Ajouter le coude au grand cylindre
162 # ----------------------------------
164 coude_quads = [ cylindres_gp.getQuadIJ(CYL_BIG, 0, 0, 0) ]
166 for nj in range (CYL_JMAX) :
167 quad = cylindres_gp.getQuadIJ(CYL_BIG, 1, nj, 0)
168 coude_quads.append(quad)
170 coude_centre = doc.addVertex(arc_x , arc_y , arc_z )
171 coude_dir = doc.addVector(arc_dx, arc_dy, arc_dz)
172 coude_angle = geompy.GetAngle(arc_a1, arc_a2)
174 coude = doc.revolutionQuads (coude_quads, coude_centre, coude_dir, [coude_angle])
175 ##coude = doc.revolutionQuadsUni (coude_quads, coude_centre, coude_dir, coude_angle, 8)
177 doc.saveVtk ("tuyauterie1.vtk")
178 doc.addLaws (0.02, True)
179 blocs = hexablock.mesh(doc)
180 muv, mue, muq, muh = hexablock.dump(doc, blocs)
184 # Associer la géométrie au modèle de bloc
185 # =======================================
187 doc.addShape(geometrie, nom)
189 # Associer les arêtes du modèle issues du joinQuads
190 # -------------------------------------------------
192 for i in xrange(hexablock.CV_MAXI_EXT):
193 asso_h = prisme.getHexa(i)
194 asso_a = [ hexablock.E_CE, hexablock.E_CF, hexablock.E_CE, hexablock.E_CF, hexablock.E_CE, hexablock.E_CF, hexablock.E_CE, hexablock.E_CF ][i]
195 asso_e = asso_h.getEdge(asso_a)
196 asso_s1 = asso_e.getVertex(0)
197 asso_s2 = asso_e.getVertex(1)
198 asso_v1 = asso_s1.getAssociation()
199 asso_v2 = asso_s2.getAssociation()
200 x1, y1, z1 = geompy.PointCoordinates(asso_v1)
201 x2, y2, z2 = geompy.PointCoordinates(asso_v2)
202 asso_pts = [ asso_v1 ]
204 pas = math.sqrt( (x2-x1)**2 + (y2-y1)**2 + (z2-z1)**2 ) / n
205 for i in xrange(1, n):
206 x = x1 + (x2 - x1) * i * pas
207 y = y1 + (y2 - y1) * i * pas
208 z = z1 + (z2 - z1) * i * pas
209 asso_v = geompy.MakeVertexOnSurfaceByCoord(face_cylindre, x, y, z)
210 asso_pts.append(asso_v)
211 asso_pts.append(asso_v2)
212 asso_curve = geompy.MakeInterpol(asso_pts)
213 asso_e.addAssociation(asso_curve, 0, 1)
215 # Associer les quadrangles du modèle issues du joinQuads
216 # ------------------------------------------------------
218 for i in xrange(hexablock.CV_MAXI_EXT):
219 asso_h = prisme.getHexa(i)
220 asso_f = [ hexablock.Q_C, hexablock.Q_F, hexablock.Q_E, hexablock.Q_F, hexablock.Q_E, hexablock.Q_F, hexablock.Q_E, hexablock.Q_F ][i]
221 asso_q = asso_h.getQuad(asso_f)
222 asso_q.addAssociation(face_cylindre)
224 # Mailler le modèle de bloc
225 # =========================
227 # Définir 4 groupes de faces
228 # --------------------------
230 groupe_d_grand = doc.addQuadGroup("grand:disque")
231 groupe_d_moyen = doc.addQuadGroup("moyen:disque")
232 groupe_d_petit = doc.addQuadGroup("petit:disque")
233 groupe_d_coude = doc.addQuadGroup("coude:disque")
236 for j in xrange( [hexablock.CV_MAXI_INT, hexablock.CV_MAXI_EXT][i] ):
237 quad = cylindres_gm.getQuadIJ(hexablock.CYL_BIG , i, j, 0)
238 groupe_d_grand.addElement(quad)
240 quad = cylindres_gm.getQuadIJ(hexablock.CYL_SMALL, i, j, 0)
241 groupe_d_moyen.addElement(quad)
243 quad = cylindres_gp.getQuadIJ(hexablock.CYL_SMALL, i, j, 0)
244 groupe_d_petit.addElement(quad)
246 for i in xrange( hexablock.CV_MAXI_INT + hexablock.CV_MAXI_EXT ):
248 quad = h.getQuad(hexablock.Q_B)
249 groupe_d_coude.addElement(quad)
251 # Définir 3 groupes de volumes
252 # ----------------------------
254 groupe_grand = doc.addHexaGroup("grand")
255 groupe_moyen = doc.addHexaGroup("moyen")
256 groupe_petit = doc.addHexaGroup("petit")
258 for i in xrange( doc.countUsedHexa() ):
259 h = doc.getUsedHexa(i)
260 groupe_grand.addElement(h)
263 for j in xrange( [hexablock.CV_MAXI_INT, hexablock.CV_MAXI_EXT][i] ):
264 h = cylindres_gm.getHexaIJK(hexablock.CYL_SMALL, i, j, 0)
265 groupe_moyen.addElement(h)
266 groupe_grand.removeElement(h)
268 h = cylindres_gp.getHexaIJK(hexablock.CYL_SMALL, i, j, 0)
269 groupe_petit.addElement(h)
270 groupe_grand.removeElement(h)
272 # Mailler le modèle de bloc avec ses associations
273 # -----------------------------------------------
275 hexablock.addLaws(doc, 0.02, True)
277 blocs = hexablock.mesh(doc)
279 muv, mue, muq, muh = hexablock.dump(doc, blocs)