Coverage for pyxfem/core/gmsh_mesh_generator.py: 24%
419 statements
« prev ^ index » next coverage.py v7.17.0, created at 2026-10-24 07:17 +0000
1#!/usr/bin/env python3
2# -*- coding: utf-8 -*-
3"""
4Created on Tue Dec 5 14:08:33 2023
6@author: nchevaug
7"""
8from typing import Tuple
10import gmsh
11def square(length=2., relh=1., transfinite = False)-> Tuple[str, dict]:
12 """ generate a gmsh mesh model of a square of lenght L and mesh density lc """
13 modelname = 'ref_square_lc_'+str(relh)
14 gmsh.model.add(modelname)
15 gmsh.model.geo.addPoint(-length/2., -length/2., 0, relh, 1)
16 gmsh.model.geo.addPoint(length/2., -length/2., 0, relh, 2)
17 gmsh.model.geo.addPoint(length/2., length/2., 0, relh, 3)
18 gmsh.model.geo.addPoint(-length/2., length/2., 0, relh, 4)
19 gmsh.model.geo.addLine(1, 2, 1)
20 gmsh.model.geo.addLine(2, 3, 2)
21 gmsh.model.geo.addLine(3, 4, 3)
22 gmsh.model.geo.addLine(4, 1, 4)
23 gmsh.model.geo.addCurveLoop([1, 2, 3, 4], 1)
24 gmsh.model.geo.addPlaneSurface([1], 1)
25 gmsh.model.geo.synchronize()
26 if transfinite :
27 gmsh.model.mesh.setTransfiniteSurface(1, arrangement="Alternate Left")
28 phys_dict = dict()
29 phys_dict['Omega'] = (2, gmsh.model.addPhysicalGroup(2, [1], -1, name='Omega'))
30 phys_dict['bottom'] = (1, gmsh.model.addPhysicalGroup(1, [1], -1, name='bottom'))
31 phys_dict['right'] = (1, gmsh.model.addPhysicalGroup(1, [2], -1, name='right'))
32 phys_dict['top'] = (1, gmsh.model.addPhysicalGroup(1, [3], -1, name='top'))
33 phys_dict['left'] = (1, gmsh.model.addPhysicalGroup(1, [4], -1, name='left'))
34 phys_dict['corner'] = (0, gmsh.model.addPhysicalGroup(0, [1], -1, name='corner'))
35 gmsh.model.mesh.generate()
36 return modelname, phys_dict
38def reg_rectangle(L, H, Nx, Ny)-> Tuple[str, dict]:
39 """ generate a gmsh mesh model of a rectangle with regular mesh of lenght L and mesh density lc """
40 modelname = 'reg_rectangle'
41 gmsh.model.add(modelname)
42 gmsh.model.geo.addPoint(-L/2., -H/2., 0, 1., 1)
43 gmsh.model.geo.addPoint(L/2., -H/2., 0, 1., 2)
44 gmsh.model.geo.addPoint(L/2., H/2., 0, 1., 3)
45 gmsh.model.geo.addPoint(-L/2., H/2., 0, 1., 4)
46 gmsh.model.geo.addLine(1, 2, 1)
47 gmsh.model.geo.addLine(2, 3, 2)
48 gmsh.model.geo.addLine(3, 4, 3)
49 gmsh.model.geo.addLine(4, 1, 4)
50 gmsh.model.geo.addCurveLoop([1, 2, 3, 4], 1)
51 gmsh.model.geo.addPlaneSurface([1], 1)
52 gmsh.model.geo.synchronize()
53 gmsh.model.mesh.setTransfiniteCurve(1, Nx)
54 gmsh.model.mesh.setTransfiniteCurve(2, Ny)
55 gmsh.model.mesh.setTransfiniteCurve(3, Nx)
56 gmsh.model.mesh.setTransfiniteCurve(4, Ny)
58 gmsh.model.mesh.setTransfiniteSurface(1, arrangement="AlternateLeft")
59 #gmsh.model.mesh.setTransfiniteSurface(1)
61 phys_dict = dict()
62 phys_dict['Omega'] = (2, gmsh.model.addPhysicalGroup(2, [1], -1, name='Omega'))
63 phys_dict['bottom'] = (1, gmsh.model.addPhysicalGroup(1, [1], -1, name='bottom'))
64 phys_dict['right'] = (1, gmsh.model.addPhysicalGroup(1, [2], -1, name='right'))
65 phys_dict['top'] = (1, gmsh.model.addPhysicalGroup(1, [3], -1, name='top'))
66 phys_dict['left'] = (1, gmsh.model.addPhysicalGroup(1, [4], -1, name='left'))
67 phys_dict['corner'] = (0, gmsh.model.addPhysicalGroup(0, [1], -1, name='corner'))
68 gmsh.model.mesh.generate()
69 return modelname, phys_dict
72def square_with_inner_seg(length:float=2., relh1:float=1., relh2:float=0.1, xc0=(-0.5, 0.), xc1=(0.5,0.))-> Tuple[str, dict]:
73 """ generate a gmsh mesh model of a square of lenght L and mesh density relh
74 that embed a segment between xc0 ans xc1
75 """
76 modelname = 'square__with_inner_seg'
77 gmsh.model.add(modelname)
78 gmsh.model.geo.addPoint(-length/2., -length/2., 0, relh1, 1)
79 gmsh.model.geo.addPoint(length/2., -length/2., 0, relh1, 2)
80 gmsh.model.geo.addPoint(length/2., length/2., 0, relh1, 3)
81 gmsh.model.geo.addPoint(-length/2., length/2., 0, relh1, 4)
82 gmsh.model.geo.addLine(1, 2, 1)
83 gmsh.model.geo.addLine(2, 3, 2)
84 gmsh.model.geo.addLine(3, 4, 3)
85 gmsh.model.geo.addLine(4, 1, 4)
86 gmsh.model.geo.addPoint( xc0[0], xc0[1], 0., relh2, 5)
87 gmsh.model.geo.addPoint( xc1[0], xc1[1], 0., relh2, 6)
89 # gmsh.model.geo.addPoint( 0.1, -0.2, 0., relh/10., 7)
90 # gmsh.model.geo.addPoint( 0.2, -0.5, 0., relh/10., 8)
92 gmsh.model.geo.addLine(5, 6, 5)
93 # gmsh.model.geo.addLine(5, 7, 5)
94 # gmsh.model.geo.addLine(6, 7, 6)
95 # gmsh.model.geo.addLine(8, 7, 7)
97 gmsh.model.geo.addCurveLoop([1, 2, 3, 4], 1)
98 gmsh.model.geo.addPlaneSurface([1], 1)
99 gmsh.model.geo.synchronize()
100 gmsh.model.mesh.embed(1, [5], 2, 1)
102 #gmsh.model.mesh.embed(1, [5, 6, 7], 2, 1)
104 phys_dict = dict()
105 phys_dict['Omega'] = (2, gmsh.model.addPhysicalGroup(2, [1], -1, name='Omega'))
106 phys_dict['Crack'] = (1, gmsh.model.addPhysicalGroup(1, [5], -1, name='Crack'))
107 phys_dict['Top'] = (1, gmsh.model.addPhysicalGroup(1, [3], -1, name='Top'))
108 phys_dict['Bottom'] = (1, gmsh.model.addPhysicalGroup(1, [1], -1, name='Bottom'))
109 gmsh.model.mesh.generate()
110 return modelname, phys_dict
112def two_square(length=2., relh=1.)-> Tuple[str, dict]:
113 """ generate a gmsh mesh model of a square of lenght L and mesh density lc """
114 modelname = 'ref_square_lc_'+str(relh)
115 gmsh.model.add(modelname)
116 gmsh.model.geo.addPoint(-length, -length/2., 0, relh, 1)
117 gmsh.model.geo.addPoint(0, -length/2., 0, relh, 2)
118 gmsh.model.geo.addPoint(0, length/2., 0, relh, 3)
119 gmsh.model.geo.addPoint(-length, length/2., 0, relh, 4)
120 gmsh.model.geo.addPoint(length, -length/2., 0, relh, 5)
121 gmsh.model.geo.addPoint(length, length/2., 0, relh, 6)
122 gmsh.model.geo.addLine(1, 2, 1)
123 gmsh.model.geo.addLine(2, 3, 2)
124 gmsh.model.geo.addLine(3, 4, 3)
125 gmsh.model.geo.addLine(4, 1, 4)
126 gmsh.model.geo.addLine(2, 5, 5)
127 gmsh.model.geo.addLine(5, 6, 6)
128 gmsh.model.geo.addLine(6, 3, 7)
129 gmsh.model.geo.addCurveLoop([1, 2, 3, 4], 1)
130 gmsh.model.geo.addCurveLoop([5, 6, 7, -2], 2)
131 gmsh.model.geo.addPlaneSurface([1], 1)
132 gmsh.model.geo.addPlaneSurface([2], 2)
133 gmsh.model.geo.synchronize()
134 phys_dict = dict()
135 phys_dict['Omega1'] = (2, gmsh.model.addPhysicalGroup(2, [1], -1, name='Omega1'))
136 phys_dict['Omega2'] = (2, gmsh.model.addPhysicalGroup(2, [2], -1, name='Omega2'))
137 phys_dict['midline'] = (2, gmsh.model.addPhysicalGroup(1, [2], -1, name='midline'))
138 phys_dict['leftline'] = (2, gmsh.model.addPhysicalGroup(1, [4], -1, name='leftline'))
139 phys_dict['rightline'] = (2, gmsh.model.addPhysicalGroup(1, [6], -1, name='rightline'))
140 phys_dict['topline'] = (2, gmsh.model.addPhysicalGroup(1, [3,7], -1, name='topline'))
141 phys_dict['bottomline'] = (2, gmsh.model.addPhysicalGroup(1, [1,5], -1, name='bottomline'))
142 gmsh.model.mesh.generate()
143 return modelname, phys_dict
145def disk_in_ref_square(length=2., radius=.5, relh=1.)-> Tuple[str, dict]:
146 """ generate a gmsh mesh model of a square of lenght L containing a centerd disk of radius r
147 and mesh density lc """
148 modelname = 'disk_in_ref_square'
149 gmsh.model.add(modelname)
150 gmsh.model.geo.addPoint(-length/2., -length/2., 0, relh, 1)
151 gmsh.model.geo.addPoint(length/2., -length/2., 0, relh, 2)
152 gmsh.model.geo.addPoint(length/2., length/2., 0, relh, 3)
153 gmsh.model.geo.addPoint(-length/2., length/2., 0, relh, 4)
154 lc2 = 1.*relh
155 gmsh.model.geo.addPoint(0., 0., 0, lc2, 5)
156 gmsh.model.geo.addPoint(0., -radius, 0, lc2, 6)
157 gmsh.model.geo.addPoint(radius, 0., 0, lc2, 7)
158 gmsh.model.geo.addPoint(0., radius, 0, lc2, 8)
159 gmsh.model.geo.addPoint(-radius, 0., 0, lc2, 9)
160 gmsh.model.geo.addLine(1, 2, 1)
161 gmsh.model.geo.addLine(2, 3, 2)
162 gmsh.model.geo.addLine(3, 4, 3)
163 gmsh.model.geo.addLine(4, 1, 4)
164 gmsh.model.geo.addCircleArc(6, 5, 7, 5)
165 gmsh.model.geo.addCircleArc(7, 5, 8, 6)
166 gmsh.model.geo.addCircleArc(8, 5, 9, 7)
167 gmsh.model.geo.addCircleArc(9, 5, 6, 8)
168 gmsh.model.geo.addCurveLoop([1, 2, 3, 4], 1)
169 gmsh.model.geo.addCurveLoop([5, 6, 7, 8], 2)
170 gmsh.model.geo.addPlaneSurface([1, -2], 1)
171 gmsh.model.geo.addPlaneSurface([2], 2)
172 gmsh.model.geo.synchronize()
173 phys_dict = dict()
174 phys_dict['Omega'] = (2, gmsh.model.addPhysicalGroup(2, [1, 2], -1, name='Omega'))
175 phys_dict['bottom'] = (1, gmsh.model.addPhysicalGroup(1, [1], -1, name='bottom'))
176 phys_dict['right'] = (1, gmsh.model.addPhysicalGroup(1, [2], -1, name='right'))
177 phys_dict['top'] = (1, gmsh.model.addPhysicalGroup(1, [3], -1, name='top'))
178 phys_dict['left'] = (1, gmsh.model.addPhysicalGroup(1, [4], -1, name='left'))
179 gmsh.model.mesh.generate()
180 return modelname, phys_dict
182def hole_in_square(length=2., radius=.5, relh=.1)-> Tuple[str, dict]:
183 """ generate a gmsh mesh model of a square of lenght L containing a
184 centered circular hole of radius r. mesh density is lc """
185 modelname = 'CircularHoleInRefSquare'
186 gmsh.model.add(modelname)
187 gmsh.model.geo.addPoint(-length/2., -length/2., 0, relh, 1)
188 gmsh.model.geo.addPoint(length/2., -length/2., 0, relh, 2)
189 gmsh.model.geo.addPoint(length/2., length/2., 0, relh, 3)
190 gmsh.model.geo.addPoint(-length/2., length/2., 0, relh, 4)
191 ratio = (3.14*radius)/length
192 lc2 = ratio*relh
193 gmsh.model.geo.addPoint(0., 0., 0, lc2, 5)
194 gmsh.model.geo.addPoint(0., -radius, 0, lc2, 6)
195 gmsh.model.geo.addPoint(radius, 0., 0, lc2, 7)
196 gmsh.model.geo.addPoint(0., radius, 0, lc2, 8)
197 gmsh.model.geo.addPoint(-radius, 0., 0, lc2, 9)
198 gmsh.model.geo.addLine(1, 2, 1)
199 gmsh.model.geo.addLine(2, 3, 2)
200 gmsh.model.geo.addLine(3, 4, 3)
201 gmsh.model.geo.addLine(4, 1, 4)
202 gmsh.model.geo.addCircleArc(7, 5, 6, 5)
203 gmsh.model.geo.addCircleArc(6, 5, 9, 6)
204 gmsh.model.geo.addCircleArc(9, 5, 8, 7)
205 gmsh.model.geo.addCircleArc(8, 5, 7, 8)
206 gmsh.model.geo.addCurveLoop([1, 2, 3, 4], 1)
207 gmsh.model.geo.addCurveLoop([5, 6, 7, 8], 2)
208 gmsh.model.geo.addPlaneSurface([1, 2], 1)
209 gmsh.model.geo.synchronize()
210 phys_dict = dict()
211 phys_dict['Omega'] = (2, gmsh.model.addPhysicalGroup(2, [1], -1, name='Omega'))
212 phys_dict['bottom'] = (1, gmsh.model.addPhysicalGroup(1, [1], -1, name='bottom'))
213 phys_dict['right'] = (1, gmsh.model.addPhysicalGroup(1, [2], -1, name='right'))
214 phys_dict['top'] = (1, gmsh.model.addPhysicalGroup(1, [3], -1, name='top'))
215 phys_dict['left'] = (1, gmsh.model.addPhysicalGroup(1, [4], -1, name='left'))
216 phys_dict['circle'] = (1, gmsh.model.addPhysicalGroup(1, [5, 6, 7, 8], -1, name='circle'))
217 gmsh.model.mesh.generate()
218 return modelname, phys_dict
220def hole_in_square_2(length=2., radius=.5, relh=.1)-> Tuple[str, dict]:
221 """ generate a gmsh mesh model of a square of lenght L containing a
222 centered circular hole of radius r. mesh density is lc """
223 modelname = 'CircularHoleInRefSquare'
224 gmsh.model.add(modelname)
225 gmsh.model.geo.addPoint(-length/2., -length/2., 0, relh, 1)
226 gmsh.model.geo.addPoint(0, -length/2., 0, relh, 10)
228 gmsh.model.geo.addPoint(length/2., -length/2., 0, relh, 2)
230 gmsh.model.geo.addPoint(length/2., length/2., 0, relh, 3)
231 gmsh.model.geo.addPoint(-length/2., length/2., 0, relh, 4)
233 gmsh.model.geo.addPoint(0, length/2., 0, relh, 20)
235 ratio = (3.14*radius)/length
236 lc2 = ratio*relh
237 gmsh.model.geo.addPoint(0., 0., 0, lc2, 5)
238 gmsh.model.geo.addPoint(0., -radius, 0, lc2, 6)
239 gmsh.model.geo.addPoint(radius, 0., 0, lc2, 7)
240 gmsh.model.geo.addPoint(0., radius, 0, lc2, 8)
241 gmsh.model.geo.addPoint(-radius, 0., 0, lc2, 9)
242 gmsh.model.geo.addLine(1, 10, 1)
243 gmsh.model.geo.addLine(10, 2, 20)
244 gmsh.model.geo.addLine(10, 6, 40)
245 gmsh.model.geo.addLine(8, 20, 50)
246 gmsh.model.geo.addLine(2, 3, 2)
247 gmsh.model.geo.addLine(3, 20, 3)
249 gmsh.model.geo.addLine(20, 4, 30)
252 gmsh.model.geo.addLine(4, 1, 4)
253 gmsh.model.geo.addCircleArc(7, 5, 6, 5)
254 gmsh.model.geo.addCircleArc(6, 5, 9, 6)
255 gmsh.model.geo.addCircleArc(9, 5, 8, 7)
256 gmsh.model.geo.addCircleArc(8, 5, 7, 8)
257 gmsh.model.geo.addCurveLoop([1, 40, 6, 7, 50, 30, 4], 1)
258 gmsh.model.geo.addCurveLoop([20, 2, 3, -50, 8, 5, -40], 2)
259 gmsh.model.geo.addPlaneSurface([1], 1)
260 gmsh.model.geo.addPlaneSurface([2], 2)
262 gmsh.model.geo.synchronize()
263 phys_dict = dict()
264 phys_dict['Omega'] = (2, gmsh.model.addPhysicalGroup(2, [1,2], -1, name='Omega'))
265 phys_dict['bottom'] = (1, gmsh.model.addPhysicalGroup(1, [1,20], -1, name='bottom'))
266 phys_dict['right'] = (1, gmsh.model.addPhysicalGroup(1, [2], -1, name='right'))
267 phys_dict['top'] = (1, gmsh.model.addPhysicalGroup(1, [3,30], -1, name='top'))
268 phys_dict['left'] = (1, gmsh.model.addPhysicalGroup(1, [4], -1, name='left'))
269 phys_dict['circle'] = (1, gmsh.model.addPhysicalGroup(1, [5, 6, 7, 8], -1, name='circle'))
270 gmsh.model.mesh.generate()
271 return modelname, phys_dict
273def crack_in_square(relh=1., L=1.,)-> Tuple[str, dict]:
274 """ generate a gmsh mesh model of a square of lenght 2L
275 containing an horizontal crack from left to square center """
276 modelname = 'square_with_crack'
277 #relh2 = relh*np.pi
278 gmsh.model.add(modelname)
279 gmsh.model.geo.addPoint(-L, -L, 0, relh, 1)
280 gmsh.model.geo.addPoint(L, -L, 0, relh, 2)
281 gmsh.model.geo.addPoint(L, L, 0, relh, 3)
282 gmsh.model.geo.addPoint(-L, L, 0, relh, 4)
283 gmsh.model.geo.addPoint(-L, 1.e-6, 0, relh, 5)
284 gmsh.model.geo.addPoint(0., 0., 0., relh, 6)
285 gmsh.model.geo.addPoint(-L, -1.e-6, 0, relh, 7)
286 gmsh.model.geo.addLine(1, 2, 1)
287 gmsh.model.geo.addLine(2, 3, 2)
288 gmsh.model.geo.addLine(3, 4, 3)
289 gmsh.model.geo.addLine(4, 5, 4)
290 gmsh.model.geo.addLine(5, 6, 5)
291 gmsh.model.geo.addLine(6, 7, 6)
292 gmsh.model.geo.addLine(7, 1, 7)
293 gmsh.model.geo.addCurveLoop([1, 2, 3, 4, 5, 6, 7], 1)
294 gmsh.model.geo.addPlaneSurface([1], 1)
295 gmsh.model.geo.synchronize()
296 gmsh.model.addPhysicalGroup(2, [1], -1, name='Omega')
297 phys_dict = dict()
298 phys_dict['Omega'] = (2, gmsh.model.addPhysicalGroup(2, [1], -1, name='Omega'))
299 phys_dict['dOmega_ext'] = (1, gmsh.model.addPhysicalGroup(1,[1,2,3,4,7], -1, name='dOmega_ext'))
300 phys_dict['dOmega_crack'] = (1, gmsh.model.addPhysicalGroup(1,[5,6], -1, name='dOmega_crack'))
301 phys_dict['top'] = gmsh.model.addPhysicalGroup(1, [3], -1, name='top')
302 phys_dict['bottom'] = gmsh.model.addPhysicalGroup(1, [1], -1, name='bottom')
303 #phys_dict['dOmega_left'] = (1, gmsh.model.addPhysicalGroup(1, [7], -1, name='dOmega_left'))
304 #phys_dict['dOmega_right'] = (1, gmsh.model.addPhysicalGroup(1, [3, 4], -1, name='dOmega_right'))
305 #phys_dict['dOmega_left'] = (1, gmsh.model.addPhysicalGroup(1, [7], -1, name='dOmega_left'))
306 gmsh.model.mesh.generate()
307 return modelname, phys_dict
309def crack_in_square_contour(relh=1., l=1., r=0.5)-> Tuple[str, dict]:
310 """ generate a gmsh mesh model of a square of lenght 2L
311 containing an horizontal crack from left to square center
312 and a circular contour for inter-*actions integral of radius R
313 """
314 modelname = 'square_with_crack'
315 eps = 1.e-10
317 gmsh.model.add(modelname)
318 gmsh.model.geo.addPoint(-l, -l, 0, relh, 1)
319 gmsh.model.geo.addPoint(l, -l, 0, relh, 2)
320 gmsh.model.geo.addPoint(l, l, 0, relh, 3)
321 gmsh.model.geo.addPoint(-l, l, 0, relh, 4)
322 gmsh.model.geo.addPoint(-l, +eps, 0, relh, 5)
323 gmsh.model.geo.addPoint(0., 0., 0., relh, 6)
324 gmsh.model.geo.addPoint(-l, -eps, 0, relh, 7)
326 gmsh.model.geo.addPoint(-r, +eps, 0, relh, 8)
327 gmsh.model.geo.addPoint(0, r, 0, relh, 9)
328 gmsh.model.geo.addPoint(r, 0, 0, relh, 10)
329 gmsh.model.geo.addPoint(0, -r, 0, relh, 11)
330 gmsh.model.geo.addPoint(-r, -eps, 0, relh, 12)
332 gmsh.model.geo.addLine(1, 2, 1)
333 gmsh.model.geo.addLine(2, 3, 2)
334 gmsh.model.geo.addLine(3, 4, 3)
335 gmsh.model.geo.addLine(4, 5, 4)
336 gmsh.model.geo.addLine(5, 8, 5)
337 gmsh.model.geo.addLine(8, 6, 6)
338 gmsh.model.geo.addLine(6, 12, 7)
339 gmsh.model.geo.addLine(12, 7, 8)
340 gmsh.model.geo.addLine(7, 1, 9)
342 gmsh.model.geo.addCircleArc(8, 6, 9, 11)
343 gmsh.model.geo.addCircleArc(9, 6, 10, 12)
344 gmsh.model.geo.addCircleArc(10, 6, 11, 13)
345 gmsh.model.geo.addCircleArc(11, 6, 12, 14)
347 gmsh.model.geo.addCurveLoop([1, 2, 3, 4, 5, 11, 12, 13, 14, 8, 9], 1)
348 gmsh.model.geo.addPlaneSurface([1], 1)
349 #gmsh.model.geo.addCurveLoop([-6, 11, 12, 13, 14, -7], 2)
350 gmsh.model.geo.addCurveLoop([ -11, -12, -13, -14, 7, 6], 2)
351 gmsh.model.geo.addPlaneSurface([2], 2)
352 gmsh.model.geo.synchronize()
353 phys_dict={}
354 phys_dict['Omega'] = (2, gmsh.model.addPhysicalGroup(2, [1,2], -1, name='Omega'))
355 phys_dict['Omega_ext'] = (2, gmsh.model.addPhysicalGroup(2, [1], -1, name='Omega_ext'))
356 phys_dict['Omega_in'] = (2, gmsh.model.addPhysicalGroup(2, [2], -1, name='Omega_in'))
357 phys_dict['bottom'] = (1, gmsh.model.addPhysicalGroup(1, [1], -1, name='bottom'))
358 phys_dict['right'] = (1, gmsh.model.addPhysicalGroup(1, [2], -1, name='right'))
359 phys_dict['top'] = (1, gmsh.model.addPhysicalGroup(1, [3], -1, name='top'))
360 phys_dict['left_top'] = (1, gmsh.model.addPhysicalGroup(1, [4], -1, name='left_top'))
361 phys_dict['left_bottom'] = (1, gmsh.model.addPhysicalGroup(1, [9], -1, name='left_bottom'))
362 phys_dict['crack'] = (1, gmsh.model.addPhysicalGroup(1, [5,6,7,8], -1, name='crack'))
363 phys_dict['contour'] = (1, gmsh.model.addPhysicalGroup(1, [11,12,13,14], -1, name='contour'))
364 gmsh.model.mesh.generate()
365 return modelname, phys_dict
367def l_shape(relh=1.)-> Tuple[str, dict]:
368 """ generate a gmsh mesh model of a L shaped domain, divided in 3 squares """
369 modelname = 'L'
370 gmsh.model.add(modelname)
371 gmsh.model.geo.addPoint(0., 0., 0., relh, 1)
372 gmsh.model.geo.addPoint(0., -1., 0., relh, 2)
373 gmsh.model.geo.addPoint(1., -1., 0., relh, 3)
374 gmsh.model.geo.addPoint(1., 0., 0., relh, 4)
375 gmsh.model.geo.addPoint(1., 1., 0., relh, 5)
376 gmsh.model.geo.addPoint(0., 1., 0., relh, 6)
377 gmsh.model.geo.addPoint(-1., 1., 0., relh, 7)
378 gmsh.model.geo.addPoint(-1., 0., 0., relh, 8)
379 gmsh.model.geo.addLine(1, 2, 1)
380 gmsh.model.geo.addLine(2, 3, 2)
381 gmsh.model.geo.addLine(3, 4, 3)
382 gmsh.model.geo.addLine(4, 5, 4)
383 gmsh.model.geo.addLine(5, 6, 5)
384 gmsh.model.geo.addLine(6, 7, 6)
385 gmsh.model.geo.addLine(7, 8, 7)
386 gmsh.model.geo.addLine(8, 1, 8)
387 gmsh.model.geo.addLine(1, 6, 9)
388 gmsh.model.geo.addLine(1, 4, 10)
389 gmsh.model.geo.addCurveLoop([1, 2, 3, -10], 1)
390 gmsh.model.geo.addCurveLoop([10, 4, 5, -9], 2)
391 gmsh.model.geo.addCurveLoop([9, 6, 7, 8], 3)
392 gmsh.model.geo.addPlaneSurface([1], 1)
393 gmsh.model.geo.addPlaneSurface([2], 2)
394 gmsh.model.geo.addPlaneSurface([3], 3)
395 gmsh.model.geo.synchronize()
396 phys_dict = dict()
397 phys_dict['dOmega_left'] = (1, gmsh.model.addPhysicalGroup(1, [7], -1, name='dOmega_left'))
398 phys_dict['dOmega_right'] = (1, gmsh.model.addPhysicalGroup(1, [3, 4], -1, name='dOmega_right'))
399 phys_dict['dOmega_other'] = (1, gmsh.model.addPhysicalGroup(1, [1, 2, 5, 6, 8], -1,
400 name='dOmega_other'))
401 phys_dict['dOmega8'] = (1, gmsh.model.addPhysicalGroup(1, [8], -1,
402 name='dOmega8'))
403 phys_dict['Omega'] = (2, gmsh.model.addPhysicalGroup(2, [1, 2, 3], -1, name='Omega'))
404 phys_dict['Omega1'] = (2, gmsh.model.addPhysicalGroup(2, [1], -1, name='Omega1'))
405 phys_dict['Omega2'] = (2, gmsh.model.addPhysicalGroup(2, [2], -1, name='Omega2'))
406 phys_dict['Omega3'] = (2, gmsh.model.addPhysicalGroup(2, [3], -1, name='Omega3'))
407 gmsh.model.mesh.generate()
408 #gmsh.model.mesh.setNode(7, np.array([-1,0.,0.]), np.array([], dtype=np.float64))
409 return modelname, phys_dict
411def ref_tri(relh=1)-> Tuple[str, dict]:
412 """ generate a gmsh mesh model of the reference triangle """
413 modelname = 'ref_tri'
414 gmsh.model.add(modelname)
415 gmsh.model.geo.addPoint(0., 0., 0., relh, 1)
416 gmsh.model.geo.addPoint(1., 0., 0., relh, 2)
417 gmsh.model.geo.addPoint(0., 1., 0., relh, 3)
418 gmsh.model.geo.addLine(1, 2, 1)
419 gmsh.model.geo.addLine(2, 3, 2)
420 gmsh.model.geo.addLine(3, 1, 3)
421 gmsh.model.geo.addCurveLoop([1, 2, 3], 1)
422 gmsh.model.geo.addPlaneSurface([1], 1)
423 gmsh.model.geo.synchronize()
424 phys_dict = dict()
425 phys_dict['Omega'] = gmsh.model.addPhysicalGroup(2, [1], -1, name='Omega')
426 phys_dict['Omega'] = gmsh.model.addPhysicalGroup(1, [1], -1, name='bottom')
427 phys_dict['Omega'] = gmsh.model.addPhysicalGroup(1, [2], -1, name='right')
428 phys_dict['Omega'] = gmsh.model.addPhysicalGroup(1, [3], -1, name='left')
429 gmsh.model.mesh.generate()
430 return modelname, phys_dict
432def cylinder(radius:float, height:float, relh:float=1.)-> Tuple[str, dict]:
433 ''' Construct a cylinder in the Z direction (just the surface)'''
434 modelname = 'cylinder'
435 gmsh.model.add(modelname)
436 gmsh.model.geo.addPoint(0., 0., 0., relh, 1)
437 gmsh.model.geo.addPoint(radius, 0., 0., relh, 2)
438 gmsh.model.geo.addPoint(0., radius, 0., relh, 3)
439 gmsh.model.geo.addPoint(-radius, 0., 0., relh, 4)
440 gmsh.model.geo.addPoint(0., -radius, 0., relh, 5)
441 gmsh.model.geo.addCircleArc(2, 1, 3, 1)
442 gmsh.model.geo.addCircleArc(3, 1, 4, 2)
443 gmsh.model.geo.addCircleArc(4, 1, 5, 3)
444 gmsh.model.geo.addCircleArc(5, 1, 2, 4)
445 gmsh.model.geo.extrude([(1, 1), (1, 2), (1, 3), (1,4)],0.,0., height)
447 #print(out_dimtags)
448 bottom_geo = [1,2,3,4]
449 top_geo = [5,9,13,17]
450 surf_geo = [8, 12, 16, 20]
451 gmsh.model.geo.synchronize()
452 phys_dict = dict()
453 phys_dict['omega'] = gmsh.model.addPhysicalGroup(2, surf_geo, -1, name='omega')
454 phys_dict['bottom'] = gmsh.model.addPhysicalGroup(1, bottom_geo, -1, name='bottom')
455 phys_dict['top'] = gmsh.model.addPhysicalGroup(1, top_geo, -1, name='top')
456 # phys_dict['Omega'] = gmsh.model.addPhysicalGroup(1, [3], -1, name='left')
457 gmsh.model.mesh.generate()
458 return modelname, phys_dict
460def structured_cylinder(
461 radius: float,
462 height: float,
463 n_theta_quarter: int = 10,
464 n_z: int = 20,
465 relh: float = 1.0,
466) -> Tuple[str, dict]:
467 modelname = "cylinder"
468 gmsh.model.add(modelname)
470 # 1. Base points
471 gmsh.model.geo.addPoint(0.0, 0.0, 0.0, relh, 1)
472 gmsh.model.geo.addPoint(radius, 0.0, 0.0, relh, 2)
473 gmsh.model.geo.addPoint(0.0, radius, 0.0, relh, 3)
474 gmsh.model.geo.addPoint(-radius, 0.0, 0.0, relh, 4)
475 gmsh.model.geo.addPoint(0.0, -radius, 0.0, relh, 5)
477 # 2. Base arcs
478 arc_tags = [1, 2, 3, 4]
479 gmsh.model.geo.addCircleArc(2, 1, 3, arc_tags[0])
480 gmsh.model.geo.addCircleArc(3, 1, 4, arc_tags[1])
481 gmsh.model.geo.addCircleArc(4, 1, 5, arc_tags[2])
482 gmsh.model.geo.addCircleArc(5, 1, 2, arc_tags[3])
484 # 3. Transfinite base curves
485 for tag in arc_tags:
486 gmsh.model.geo.mesh.setTransfiniteCurve(tag, n_theta_quarter + 1)
488 # 4. Extrude one curve at a time to accurately isolate top vs. lateral lines
489 top_geo = []
490 surf_geo = []
491 lat_lines = set()
493 for arc in arc_tags:
494 out = gmsh.model.geo.extrude(
495 [(1, arc)],
496 0.0,
497 0.0,
498 height,
499 numElements=[n_z],
500 recombine=False,
501 )
502 # out layout for a curve extrusion:
503 # out[0] -> (2, surface_tag)
504 # out[1] -> (1, top_curve_tag)
505 # out[2:] -> (1, lateral_edge_tags)
506 print(out)
507 surf_geo.append(out[1][1])
508 top_geo.append(out[0][1])
509 for dim, tag in out[2:]:
510 if dim == 1:
511 lat_lines.add(tag)
513 bottom_geo = arc_tags
515 # 5. Transfinite surfaces for structured triangle mesh
516 for s_tag in surf_geo:
517 gmsh.model.geo.mesh.setTransfiniteSurface(s_tag, arrangement="Alternate")
519 gmsh.model.geo.synchronize()
521 # 6. Physical groups
522 phys_dict = {}
523 # phys_dict["omega"] = gmsh.model.addPhysicalGroup(
524 # 2, surf_geo, -1, name="omega"
525 # )
526 phys_dict["bottom"] = gmsh.model.addPhysicalGroup(
527 1, bottom_geo, -1, name="bottom"
528 )
529 phys_dict["top"] = gmsh.model.addPhysicalGroup(1, top_geo, -1, name="top")
531 gmsh.model.mesh.generate()
532 return modelname, phys_dict
534if __name__ == '__main__':
535 gmsh.initialize()
536 #crack_in_square_contour(relh=0.05, L=1., R=0.5)
537 crack_in_square(relh=0.05, L=1.)
540 gmsh.fltk.run()
541 gmsh.finalize()