1 #@ MODIF ajout_quad_gmsh Macro DATE 14/09/2004 AUTEUR MCOURTOI M.COURTOIS
2 # -*- coding: iso-8859-1 -*-
3 # CONFIGURATION MANAGEMENT OF EDF VERSION
4 # ======================================================================
5 # COPYRIGHT (C) 1991 - 2002 EDF R&D WWW.CODE-ASTER.ORG
6 # THIS PROGRAM IS FREE SOFTWARE; YOU CAN REDISTRIBUTE IT AND/OR MODIFY
7 # IT UNDER THE TERMS OF THE GNU GENERAL PUBLIC LICENSE AS PUBLISHED BY
8 # THE FREE SOFTWARE FOUNDATION; EITHER VERSION 2 OF THE LICENSE, OR
9 # (AT YOUR OPTION) ANY LATER VERSION.
11 # THIS PROGRAM IS DISTRIBUTED IN THE HOPE THAT IT WILL BE USEFUL, BUT
12 # WITHOUT ANY WARRANTY; WITHOUT EVEN THE IMPLIED WARRANTY OF
13 # MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE. SEE THE GNU
14 # GENERAL PUBLIC LICENSE FOR MORE DETAILS.
16 # YOU SHOULD HAVE RECEIVED A COPY OF THE GNU GENERAL PUBLIC LICENSE
17 # ALONG WITH THIS PROGRAM; IF NOT, WRITE TO EDF R&D CODE_ASTER,
18 # 1 AVENUE DU GENERAL DE GAULLE, 92141 CLAMART CEDEX, FRANCE.
19 # ======================================================================
23 # script PYTHON de modification d'un fichier de maillage GMSH (*.msh)
24 # pour transformer les mailles lineiques en mailles quadratiques
25 # SEG2->SEG3 / TRIA3->TRIA6 / QUAD4->QUAD8,
26 # TETRA4 -> TETRA10 / PENTA6 -> PENTA15 / PYRAM5 -> PYRAM13
28 def ajout_quad_gmsh(texte):
33 # construction du dictionnaire nom du noeud / coordonnees : dict_no
34 # construction de la liste des noeuds : l_no
38 typ_mail['2D']=['2','3']
39 typ_mail['3D']=['4','5','6','7']
42 texte_eclate=texte.split('\n')
43 nbno=int(texte_eclate[texte_eclate.index('$NOD')+1])
44 l_no=texte_eclate[texte_eclate.index('$NOD')+2:texte_eclate.index('$ENDNOD')]
46 for i in texte_eclate[texte_eclate.index('$NOD')+2:texte_eclate.index('$ENDNOD')]:
49 for j in cc[1:]:dict_no[cc[0]].append(float(j))
51 l_el1=texte_eclate[texte_eclate.index('$ELM')+2:texte_eclate.index('$ENDELM')]
52 nbel =texte_eclate[texte_eclate.index('$ELM')+1]
55 # construction du dictionnaire : element / liste des aretes de l'element : elems_aretes
56 # et de son inverse : aretes / elements contenant cette arete : aretes_elems
62 connec0=elem.split(' ')[5:]
63 while '' in connec0 : connec0.remove('')
64 aa=elem.split(' ')[:4]
66 # attention : indicateur de type de maille : ajouter 7 pour passer
67 # de TRIA3 a TRIA6, de SEG2 a SEG3, de QUAD4 a QUAD8 (voir inigms.f)
69 # si maille POI1, on ne fait rien
71 if aa[1] not in typ_mail['PO'] : typel=str(int(aa[1])+7)
73 nom_elem=aa[0]+' '+typel+' '+aa[2]+' '+aa[3]
76 # traitement des mailles lineaires : un seul segment
77 if aa[1] in typ_mail['LI']:
78 segments.append([int(connec0[0]),int(connec0[1])])
80 # traitement des mailles planes (TRI3 QUA4) : on cree les segments en prenant les noeuds consecutivement, puis on ferme
81 elif aa[1] in typ_mail['2D']:
82 for i in range(len(connec0)-1) : segments.append([int(connec0[i]),int(connec0[i+1])])
83 segments.append([int(connec0[-1]),int(connec0[0])])
85 # traitement des mailles TETRA4 : pour comprendre, voir "fonctions de forme des elements isoparametriques" R3.01.01
87 for i in range(0,2) : segments.append([int(connec0[i]),int(connec0[i+1])])
88 segments.append([int(connec0[2]),int(connec0[0])])
89 for i in range(0,3) : segments.append([int(connec0[3]),int(connec0[i])])
91 # traitement des mailles HEXA8 : pour comprendre, voir "fonctions de forme des elements isoparametriques" R3.01.01
93 for i in range(0,3) : segments.append([int(connec0[i]),int(connec0[i+1])])
94 segments.append([int(connec0[3]),int(connec0[0])])
95 for i in range(4,7) : segments.append([int(connec0[i]),int(connec0[i+1])])
96 segments.append([int(connec0[7]),int(connec0[4])])
97 for i in range(0,4) : segments.append([int(connec0[i]),int(connec0[i+4])])
99 # traitement des mailles PENTA6 : pour comprendre, voir "fonctions de forme des elements isoparametriques" R3.01.01
101 for i in range(0,2) : segments.append([int(connec0[i]),int(connec0[i+1])])
102 segments.append([int(connec0[2]),int(connec0[0])])
103 for i in range(3,5) : segments.append([int(connec0[i]),int(connec0[i+1])])
104 segments.append([int(connec0[5]),int(connec0[3])])
105 for i in range(0,3) : segments.append([int(connec0[i]),int(connec0[i+3])])
107 # traitement des mailles PYRA5 : pour comprendre, voir "fonctions de forme des elements isoparametriques" R3.01.01
109 for i in range(0,3) : segments.append([int(connec0[i]),int(connec0[i+1])])
110 segments.append([int(connec0[3]),int(connec0[0])])
111 for i in range(0,4) : segments.append([int(connec0[4]),int(connec0[i])])
113 elems_aretes[nom_elem]=[]
114 for seg in segments :
115 elems_aretes[nom_elem].append(copy.copy(seg))
117 for seg in segments :
118 if aretes_elems.has_key(tuple(seg)):aretes_elems[tuple(seg)].append(nom_elem)
119 else :aretes_elems[tuple(seg)]=[nom_elem]
120 l_el.append(nom_elem)
122 # construction de la liste complete des aretes
124 l_ar=aretes_elems.keys()
127 # numero du plus grand noeud :
128 nmax=int(l_no[-1].split(' ')[0])
130 # nouveau nombre total de noeuds :
133 # insertion des noeuds milieux dans le dictionnaire "aretes" et la liste "l_no" :
137 x=(dict_no[str(i[0])][0]+dict_no[str(i[1])][0])/2
138 y=(dict_no[str(i[0])][1]+dict_no[str(i[1])][1])/2
139 z=(dict_no[str(i[0])][2]+dict_no[str(i[1])][2])/2
140 l_no.append(str(ind)+' '+str(x)+' '+str(y)+' '+str(z))
141 for elem in aretes_elems[i]:
142 for ar in elems_aretes[elem]:
149 # re-ecriture du fichier avec les noeuds milieux :
151 resu='$NOD\n'+str(ndtot)+'\n'
154 resu=resu+'$ENDNOD\n'+'$ELM\n'+nbel+'\n'
157 if aa[1] not in typ_mail['PO']:
158 typ=str(int(aa[1])-7)
160 if elems_aretes[i]!=[]:
163 # traitement des mailles lineaires : d'abord les noeuds sommets, puis le noeud milieu
164 if typ in typ_mail['LI']:
166 for j in elems_aretes[i]:
167 resu=resu+str(j[0])+' '+str(j[2])+' '+str(j[1])+' '
169 # traitement des mailles planes (TRI3 QUA4) : d'abord les noeuds sommets, puis, dans le meme ordre, les noeuds milieux
170 elif typ in typ_mail['2D']:
171 resu=resu+' '+str(len(elems_aretes[i])*2)+' '
172 for j in elems_aretes[i]:
173 resu=resu+str(j[0])+' '
174 for j in elems_aretes[i]:
175 resu=resu+str(j[1])+' '
177 # traitement des mailles TETRA4 : pour comprendre, voir "fonctions de forme des elements isoparametriques" R3.01.01
180 for j in elems_aretes[i][:4]:
181 resu=resu+str(j[0])+' '
182 resu=resu+str(elems_aretes[i][0][1])+' '
183 resu=resu+str(elems_aretes[i][1][1])+' '
184 resu=resu+str(elems_aretes[i][2][1])+' '
185 resu=resu+str(elems_aretes[i][3][1])+' '
186 resu=resu+str(elems_aretes[i][4][1])+' '
187 resu=resu+str(elems_aretes[i][5][1])+' '
189 # traitement des mailles HEXA8 : pour comprendre, voir "fonctions de forme des elements isoparametriques" R3.01.01
192 for j in elems_aretes[i][:8]:
193 resu=resu+str(j[0])+' '
194 resu=resu+str(elems_aretes[i][0][1])+' '
195 resu=resu+str(elems_aretes[i][1][1])+' '
196 resu=resu+str(elems_aretes[i][2][1])+' '
197 resu=resu+str(elems_aretes[i][3][1])+' '
199 resu=resu+str(elems_aretes[i][8][1])+' '
200 resu=resu+str(elems_aretes[i][9][1])+' '
201 resu=resu+str(elems_aretes[i][10][1])+' '
202 resu=resu+str(elems_aretes[i][11][1])+' '
204 resu=resu+str(elems_aretes[i][4][1])+' '
205 resu=resu+str(elems_aretes[i][5][1])+' '
206 resu=resu+str(elems_aretes[i][6][1])+' '
207 resu=resu+str(elems_aretes[i][7][1])+' '
209 # traitement des mailles PENTA6 : pour comprendre, voir "fonctions de forme des elements isoparametriques" R3.01.01
212 for j in elems_aretes[i][:6]:
213 resu=resu+str(j[0])+' '
214 resu=resu+str(elems_aretes[i][0][1])+' '
215 resu=resu+str(elems_aretes[i][1][1])+' '
216 resu=resu+str(elems_aretes[i][2][1])+' '
218 resu=resu+str(elems_aretes[i][6][1])+' '
219 resu=resu+str(elems_aretes[i][7][1])+' '
220 resu=resu+str(elems_aretes[i][8][1])+' '
222 resu=resu+str(elems_aretes[i][3][1])+' '
223 resu=resu+str(elems_aretes[i][4][1])+' '
224 resu=resu+str(elems_aretes[i][5][1])+' '
226 # traitement des mailles PYRA5 : pour comprendre, voir "fonctions de forme des elements isoparametriques" R3.01.01
229 for j in elems_aretes[i][:4]:
230 resu=resu+str(j[0])+' '
231 resu=resu+str(elems_aretes[i][0][1])+' '
232 resu=resu+str(elems_aretes[i][1][1])+' '
233 resu=resu+str(elems_aretes[i][2][1])+' '
235 resu=resu+str(elems_aretes[i][6][1])+' '
236 resu=resu+str(elems_aretes[i][7][1])+' '
237 resu=resu+str(elems_aretes[i][8][1])+' '
239 resu=resu+str(elems_aretes[i][3][1])+' '
240 resu=resu+str(elems_aretes[i][4][1])+' '
241 resu=resu+str(elems_aretes[i][5][1])+' '
244 # traitement des mailles POI1 : on recopie la connectivite telle quelle, sans ajouter de nouveau noeud
245 resu=resu+l_el1[l_el.index(i)]
249 resu=resu+'$ENDELM\n'