File: t14.jl

package info (click to toggle)
gmsh 4.14.0%2Bds1-1
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid
  • size: 96,556 kB
  • sloc: cpp: 438,695; ansic: 114,912; f90: 15,477; python: 14,025; yacc: 7,333; java: 3,491; lisp: 3,194; lex: 631; perl: 571; makefile: 497; sh: 439; xml: 414; javascript: 113; pascal: 35; modula3: 32
file content (137 lines) | stat: -rw-r--r-- 4,563 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
# ------------------------------------------------------------------------------
#
#  Gmsh Julia tutorial 14
#
#  Homology and cohomology computation
#
# ------------------------------------------------------------------------------

# Homology computation in Gmsh finds representative chains of (relative)
# (co)homology space bases using a mesh of a model.  The representative basis
# chains are stored in the mesh as physical groups of Gmsh, one for each chain.

import gmsh

gmsh.initialize(append!(["gmsh"], ARGS))

# Create an example geometry
gmsh.model.add("t14")

m = 0.5  # mesh size
h = 2  # geometry height in the z-direction

gmsh.model.geo.addPoint(0, 0, 0, m, 1)
gmsh.model.geo.addPoint(10, 0, 0, m, 2)
gmsh.model.geo.addPoint(10, 10, 0, m, 3)
gmsh.model.geo.addPoint(0, 10, 0, m, 4)

gmsh.model.geo.addPoint(4, 4, 0, m, 5)
gmsh.model.geo.addPoint(6, 4, 0, m, 6)
gmsh.model.geo.addPoint(6, 6, 0, m, 7)
gmsh.model.geo.addPoint(4, 6, 0, m, 8)

gmsh.model.geo.addPoint(2, 0, 0, m, 9)
gmsh.model.geo.addPoint(8, 0, 0, m, 10)
gmsh.model.geo.addPoint(2, 10, 0, m, 11)
gmsh.model.geo.addPoint(8, 10, 0, m, 12)

gmsh.model.geo.addLine(1, 9, 1)
gmsh.model.geo.addLine(9, 10, 2)
gmsh.model.geo.addLine(10, 2, 3)

gmsh.model.geo.addLine(2, 3, 4)
gmsh.model.geo.addLine(3, 12, 5)
gmsh.model.geo.addLine(12, 11, 6)

gmsh.model.geo.addLine(11, 4, 7)
gmsh.model.geo.addLine(4, 1, 8)
gmsh.model.geo.addLine(5, 6, 9)

gmsh.model.geo.addLine(6, 7, 10)
gmsh.model.geo.addLine(7, 8, 11)
gmsh.model.geo.addLine(8, 5, 12)

gmsh.model.geo.addCurveLoop([6, 7, 8, 1, 2, 3, 4, 5], 13)
gmsh.model.geo.addCurveLoop([11, 12, 9, 10], 14)
gmsh.model.geo.addPlaneSurface([13, 14], 15)

e = gmsh.model.geo.extrude([(2, 15)], 0, 0, h)

gmsh.model.geo.synchronize()

# Create physical groups, which are used to define the domain of the
# (co)homology computation and the subdomain of the relative (co)homology
# computation.

# Whole domain
domain_tag = e[2][2]
domain_physical_tag = 1001
gmsh.model.addPhysicalGroup(3, [domain_tag], domain_physical_tag,
                            "Whole domain")

# Four "terminals" of the model
terminal_tags = [e[4][2], e[6][2], e[8][2], e[10][2]]
terminals_physical_tag = 2001
gmsh.model.addPhysicalGroup(2, terminal_tags, terminals_physical_tag,
                            "Terminals")

# Find domain boundary tags
boundary_dimtags = gmsh.model.getBoundary([(3, domain_tag)], false, false)
boundary_tags = []
complement_tags = []

for tag in boundary_dimtags
    push!(complement_tags, tag[2])
    push!(boundary_tags, tag[2])
end

for tag in terminal_tags
    deleteat!(complement_tags, findall(x -> x == tag, complement_tags))
end

# Whole domain surface
boundary_physical_tag = 2002
gmsh.model.addPhysicalGroup(2, boundary_tags, boundary_physical_tag,
                            "Boundary")

# Complement of the domain surface with respect to the four terminals
complement_physical_tag = 2003
gmsh.model.addPhysicalGroup(2, complement_tags, complement_physical_tag,
                            "Complement")

# Find bases for relative homology spaces of the domain modulo the four
# terminals.
gmsh.model.mesh.addHomologyRequest("Homology", [domain_physical_tag],
                                   [terminals_physical_tag], [0, 1, 2, 3])

# Find homology space bases isomorphic to the previous bases: homology spaces
# modulo the non-terminal domain surface, a.k.a the thin cuts.
gmsh.model.mesh.addHomologyRequest("Homology", [domain_physical_tag],
                                   [complement_physical_tag], [0, 1, 2, 3])

# Find cohomology space bases isomorphic to the previous bases: cohomology
# spaces of the domain modulo the four terminals, a.k.a the thick cuts.
gmsh.model.mesh.addHomologyRequest("Cohomology", [domain_physical_tag],
                                   [terminals_physical_tag], [0, 1, 2, 3])

# more examples
# gmsh.model.mesh.addHomologyRequest()
# gmsh.model.mesh.addHomologyRequest("Homology", [domain_physical_tag])
# gmsh.model.mesh.addHomologyRequest("Homology", [domain_physical_tag],
#                                    [boundary_physical_tag], [0,1,2,3])

# Generate the mesh and perform the requested homology computations
gmsh.model.mesh.generate(3)

# For more information, see M. Pellikka, S. Suuriniemi, L. Kettunen and
# C. Geuzaine. Homology and cohomology computation in finite element
# modeling. SIAM Journal on Scientific Computing 35(5), pp. 1195-1214, 2013.

gmsh.write("t14.msh")

# Launch the GUI to see the results:
if !("-nopopup" in ARGS)
    gmsh.fltk.run()
end

gmsh.finalize()