-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathManual on cell_geom.py.html
executable file
·270 lines (250 loc) · 26.2 KB
/
Manual on cell_geom.py.html
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
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN"
"http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd">
<html xmlns="http://www.w3.org/1999/xhtml">
<head>
<meta http-equiv="Content-Type" content="text/html; charset=utf-8" />
<title>cell_geom.py — Unit Cell Module Documentation 1.0 documentation</title>
<link rel="stylesheet" href="_static/default.css" type="text/css" />
<link rel="stylesheet" href="_static/pygments.css" type="text/css" />
<script type="text/javascript">
var DOCUMENTATION_OPTIONS = {
URL_ROOT: './',
VERSION: '1.0',
COLLAPSE_INDEX: false,
FILE_SUFFIX: '.html',
HAS_SOURCE: true
};
</script>
<script type="text/javascript" src="_static/jquery.js"></script>
<script type="text/javascript" src="_static/underscore.js"></script>
<script type="text/javascript" src="_static/doctools.js"></script>
<link rel="top" title="Unit Cell Module Documentation 1.0 documentation" href="index.html" />
<link rel="next" title="cell_material.py" href="Manual on cell_material.py.html" />
<link rel="prev" title="Welcome to Unit Cell Module’s documentation!" href="index.html" />
</head>
<body>
<div class="related">
<h3>Navigation</h3>
<ul>
<li class="right" style="margin-right: 10px">
<a href="genindex.html" title="General Index"
accesskey="I">index</a></li>
<li class="right" >
<a href="py-modindex.html" title="Python Module Index"
>modules</a> |</li>
<li class="right" >
<a href="Manual on cell_material.py.html" title="cell_material.py"
accesskey="N">next</a> |</li>
<li class="right" >
<a href="index.html" title="Welcome to Unit Cell Module’s documentation!"
accesskey="P">previous</a> |</li>
<li><a href="index.html">Unit Cell Module Documentation 1.0 documentation</a> »</li>
</ul>
</div>
<div class="document">
<div class="documentwrapper">
<div class="bodywrapper">
<div class="body">
<div class="section" id="cell-geom-py">
<h1>cell_geom.py<a class="headerlink" href="#cell-geom-py" title="Permalink to this headline">¶</a></h1>
<ul class="simple">
<li><a class="reference external" href="#overview">Overview</a></li>
<li><a class="reference external" href="#inclusions">Inclusions</a><ul>
<li><a class="reference external" href="#d-case">2D Case</a></li>
<li><a class="reference external" href="#id1">3D Case</a></li>
</ul>
</li>
<li><a class="reference external" href="#peirodic-boundary-condition">Peirodic Boundary Condition</a></li>
</ul>
<div class="section" id="overview">
<h2>Overview<a class="headerlink" href="#overview" title="Permalink to this headline">¶</a></h2>
<p>In this file <tt class="docutils literal"><span class="pre">class</span> <span class="pre">UnitCell</span></tt> is defined, where possible inclusions
can be added to the unit cell. The member methods of this class are <tt class="docutils literal"><span class="pre">set_append_inclusion</span></tt>, <tt class="docutils literal"><span class="pre">add_mark_boundary</span></tt>,
<tt class="docutils literal"><span class="pre">view_mesh</span></tt>, and <tt class="docutils literal"><span class="pre">view_domain</span></tt>. The instance of this method is
instantiated with a <tt class="docutils literal"><span class="pre">Mesh</span></tt> object in <em>FEniCS</em>. A <tt class="docutils literal"><span class="pre">UnitCell</span></tt> instance
can be either two dimensional or three dimensional.</p>
<p>Classes for creation of inclusions are included in the current file,
namely <tt class="docutils literal"><span class="pre">InclusionCircle</span></tt> and <tt class="docutils literal"><span class="pre">InclusionRectangle</span></tt>. Besides,
<tt class="docutils literal"><span class="pre">PeriodicBoundary_no_corner</span></tt> is a class specifying the periodic map
for <em>periodic boundary condition</em> in homogenization problem.</p>
</div>
<div class="section" id="inclusions">
<h2>Inclusions<a class="headerlink" href="#inclusions" title="Permalink to this headline">¶</a></h2>
<p>Setting a unit cell and its inclusions is introduced in this part. We
first import modules</p>
<div class="code python highlight-python"><div class="highlight"><pre><span class="kn">from</span> <span class="nn">dolfin</span> <span class="kn">import</span> <span class="o">*</span>
<span class="kn">import</span> <span class="nn">sys</span>
<span class="c"># Include the module directory</span>
<span class="n">sys</span><span class="o">.</span><span class="n">path</span><span class="o">.</span><span class="n">append</span><span class="p">(</span><span class="s">'./'</span><span class="p">)</span>
<span class="kn">import</span> <span class="nn">cell_geom</span> <span class="kn">as</span> <span class="nn">geom</span>
</pre></div>
</div>
<div class="section" id="d-case">
<h3>2D Case<a class="headerlink" href="#d-case" title="Permalink to this headline">¶</a></h3>
<p><strong>Import mesh and instantiation</strong></p>
<div class="code python highlight-python"><div class="highlight"><pre><span class="n">mesh</span> <span class="o">=</span> <span class="n">Mesh</span><span class="p">(</span><span class="s">r"../m.xml"</span><span class="p">)</span>
<span class="c"># Generate Inclusion</span>
<span class="n">inc1</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">InclusionCircle</span><span class="p">(</span><span class="mi">2</span><span class="p">,</span> <span class="p">(</span><span class="mf">0.5</span><span class="p">,</span> <span class="mf">0.5</span><span class="p">),</span> <span class="mf">0.25</span><span class="p">)</span>
<span class="n">inc_group</span> <span class="o">=</span> <span class="p">{</span><span class="s">'circle_inc1'</span><span class="p">:</span> <span class="n">inc1</span><span class="p">}</span>
<span class="c"># Initiate UnitCell Instance with Inclusion</span>
<span class="n">cell</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">UnitCell</span><span class="p">(</span><span class="n">mesh</span><span class="p">,</span> <span class="n">inc_group</span><span class="p">)</span>
<span class="n">cell</span><span class="o">.</span><span class="n">view_domain</span><span class="p">()</span>
</pre></div>
</div>
<p><strong>Multiple inclusions and append inclusionis</strong></p>
<div class="code python highlight-python"><div class="highlight"><pre><span class="n">mesh</span> <span class="o">=</span> <span class="n">UnitSquareMesh</span><span class="p">(</span><span class="mi">40</span><span class="p">,</span> <span class="mi">40</span><span class="p">,</span> <span class="s">'crossed'</span><span class="p">)</span>
<span class="c"># Instantiation with inclusions</span>
<span class="n">inc1</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">InclusionCircle</span><span class="p">(</span><span class="mi">2</span><span class="p">,</span> <span class="p">(</span><span class="mf">0.1</span><span class="p">,</span> <span class="mf">0.1</span><span class="p">),</span> <span class="mf">0.5</span><span class="p">)</span>
<span class="n">inc2</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">InclusionCircle</span><span class="p">(</span><span class="mi">2</span><span class="p">,</span> <span class="p">(</span><span class="mf">0.9</span><span class="p">,</span> <span class="mf">0.9</span><span class="p">),</span> <span class="mf">0.5</span><span class="p">)</span>
<span class="n">inc_group_1</span> <span class="o">=</span> <span class="p">{</span><span class="s">'circle_inc1'</span><span class="p">:</span> <span class="n">inc1</span><span class="p">,</span> <span class="s">'circle_inc2'</span><span class="p">:</span> <span class="n">inc2</span><span class="p">,}</span>
<span class="n">cell</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">UnitCell</span><span class="p">(</span><span class="n">mesh</span><span class="p">,</span> <span class="n">inc_group_1</span><span class="p">)</span>
<span class="n">cell</span><span class="o">.</span><span class="n">view_domain</span><span class="p">()</span>
<span class="c"># Another group of inlusions</span>
<span class="n">inc3</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">InclusionRectangle</span><span class="p">(</span><span class="mi">2</span><span class="p">,</span> <span class="mf">0.1</span><span class="p">,</span> <span class="mf">0.3</span><span class="p">,</span> <span class="mf">0.7</span><span class="p">,</span> <span class="mf">0.9</span><span class="p">)</span>
<span class="n">inc4</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">InclusionRectangle</span><span class="p">(</span><span class="mi">2</span><span class="p">,</span> <span class="mf">0.7</span><span class="p">,</span> <span class="mf">0.9</span><span class="p">,</span> <span class="mf">0.1</span><span class="p">,</span> <span class="mf">0.3</span><span class="p">)</span>
<span class="n">inc_group_2</span> <span class="o">=</span> <span class="p">{</span><span class="s">'rect_inc3'</span><span class="p">:</span> <span class="n">inc3</span><span class="p">,</span> <span class="s">'rect_inc4'</span><span class="p">:</span> <span class="n">inc4</span><span class="p">}</span>
<span class="c"># Append inclusions and view</span>
<span class="n">cell</span><span class="o">.</span><span class="n">set_append_inclusion</span><span class="p">(</span><span class="n">inc_group_2</span><span class="p">)</span>
<span class="n">cell</span><span class="o">.</span><span class="n">view_domain</span><span class="p">()</span>
</pre></div>
</div>
</div>
<div class="section" id="id1">
<h3>3D Case<a class="headerlink" href="#id1" title="Permalink to this headline">¶</a></h3>
<p><strong>Multiple inclusions and append inclusions</strong></p>
<div class="code python highlight-python"><div class="highlight"><pre><span class="n">mesh</span> <span class="o">=</span> <span class="n">UnitCubeMesh</span><span class="p">(</span><span class="mi">20</span><span class="p">,</span> <span class="mi">20</span><span class="p">,</span> <span class="mi">20</span><span class="p">)</span>
<span class="c"># 9 Inclusions with 8 corner inclusions and one sphere inclusion in the center</span>
<span class="n">inc</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">InclusionCircle</span><span class="p">(</span><span class="mi">3</span><span class="p">,</span> <span class="mf">0.5</span><span class="p">)</span>
<span class="n">inc1</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">InclusionRectangle</span><span class="p">(</span><span class="mi">3</span><span class="p">,</span> <span class="mf">0.</span><span class="p">,</span> <span class="mf">0.3</span><span class="p">,</span> <span class="mf">0.</span><span class="p">,</span> <span class="mf">0.3</span><span class="p">,</span> <span class="mf">0.</span><span class="p">,</span> <span class="mf">0.3</span><span class="p">)</span>
<span class="n">inc2</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">InclusionRectangle</span><span class="p">(</span><span class="mi">3</span><span class="p">,</span> <span class="mf">0.</span><span class="p">,</span> <span class="mf">0.3</span><span class="p">,</span> <span class="mf">0.</span><span class="p">,</span> <span class="mf">0.3</span><span class="p">,</span> <span class="mf">0.7</span><span class="p">,</span> <span class="mf">1.</span><span class="p">)</span>
<span class="n">inc3</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">InclusionRectangle</span><span class="p">(</span><span class="mi">3</span><span class="p">,</span> <span class="mf">0.</span><span class="p">,</span> <span class="mf">0.3</span><span class="p">,</span> <span class="mf">0.7</span><span class="p">,</span> <span class="mf">1.</span><span class="p">,</span> <span class="mf">0.</span><span class="p">,</span> <span class="mf">0.3</span><span class="p">)</span>
<span class="n">inc4</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">InclusionRectangle</span><span class="p">(</span><span class="mi">3</span><span class="p">,</span> <span class="mf">0.</span><span class="p">,</span> <span class="mf">0.3</span><span class="p">,</span> <span class="mf">0.7</span><span class="p">,</span> <span class="mf">1.</span><span class="p">,</span> <span class="mf">0.7</span><span class="p">,</span> <span class="mf">1.</span><span class="p">)</span>
<span class="n">inc5</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">InclusionRectangle</span><span class="p">(</span><span class="mi">3</span><span class="p">,</span> <span class="mf">0.7</span><span class="p">,</span> <span class="mf">1.</span><span class="p">,</span> <span class="mf">0.</span><span class="p">,</span> <span class="mf">0.3</span><span class="p">,</span> <span class="mf">0.</span><span class="p">,</span> <span class="mf">0.3</span><span class="p">)</span>
<span class="n">inc6</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">InclusionRectangle</span><span class="p">(</span><span class="mi">3</span><span class="p">,</span> <span class="mf">0.7</span><span class="p">,</span> <span class="mf">1.</span><span class="p">,</span> <span class="mf">0.</span><span class="p">,</span> <span class="mf">0.3</span><span class="p">,</span> <span class="mf">0.7</span><span class="p">,</span> <span class="mf">1.</span><span class="p">)</span>
<span class="n">inc7</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">InclusionRectangle</span><span class="p">(</span><span class="mi">3</span><span class="p">,</span> <span class="mf">0.7</span><span class="p">,</span> <span class="mf">1.</span><span class="p">,</span> <span class="mf">0.7</span><span class="p">,</span> <span class="mf">1.</span><span class="p">,</span> <span class="mf">0.</span><span class="p">,</span> <span class="mf">0.3</span><span class="p">)</span>
<span class="n">inc8</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">InclusionRectangle</span><span class="p">(</span><span class="mi">3</span><span class="p">,</span> <span class="mf">0.7</span><span class="p">,</span> <span class="mf">1.</span><span class="p">,</span> <span class="mf">0.7</span><span class="p">,</span> <span class="mf">1.</span><span class="p">,</span> <span class="mf">0.7</span><span class="p">,</span> <span class="mf">1.</span><span class="p">)</span>
<span class="n">inc_group</span> <span class="o">=</span> <span class="p">{</span><span class="s">'circle'</span><span class="p">:</span> <span class="n">inc</span><span class="p">,</span> <span class="s">'corner1'</span><span class="p">:</span> <span class="n">inc1</span><span class="p">,</span> <span class="s">'corner2'</span><span class="p">:</span> <span class="n">inc2</span><span class="p">,</span>
<span class="s">'corner3'</span><span class="p">:</span> <span class="n">inc3</span><span class="p">,</span> <span class="s">'corner4'</span><span class="p">:</span> <span class="n">inc4</span><span class="p">,</span> <span class="s">'corner5'</span><span class="p">:</span> <span class="n">inc5</span><span class="p">,</span>
<span class="s">'corner6'</span><span class="p">:</span> <span class="n">inc6</span><span class="p">,</span> <span class="s">'corner7'</span><span class="p">:</span> <span class="n">inc7</span><span class="p">,</span> <span class="s">'corner8'</span><span class="p">:</span> <span class="n">inc8</span><span class="p">}</span>
<span class="n">cell</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">UnitCell</span><span class="p">(</span><span class="n">mesh</span><span class="p">,</span> <span class="n">inc_group</span><span class="p">)</span>
<span class="n">cell</span><span class="o">.</span><span class="n">view_domain</span><span class="p">()</span>
</pre></div>
</div>
</div>
</div>
<div class="section" id="peirodic-boundary-condition">
<h2>Peirodic Boundary Condition<a class="headerlink" href="#peirodic-boundary-condition" title="Permalink to this headline">¶</a></h2>
<p>Periodic mapping for FunctionSpace initiallization. Both 2D case and 3D
case are covered. This periodic mapping excludes corners of unit cell.
In unit cell computation these corners are set fixed to prevent rigid
body movements.</p>
<div class="code python highlight-python"><div class="highlight"><pre><span class="c"># 2D</span>
<span class="n">a</span><span class="p">,</span> <span class="n">b</span> <span class="o">=</span> <span class="mi">3</span><span class="p">,</span> <span class="mi">6</span>
<span class="n">mesh_2d</span> <span class="o">=</span> <span class="n">UnitSquareMesh</span><span class="p">(</span><span class="n">a</span><span class="p">,</span> <span class="n">b</span><span class="p">)</span>
<span class="n">FS_2d</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">FunctionSpace</span><span class="p">(</span><span class="n">mesh_2d</span><span class="p">,</span> <span class="s">'CG'</span><span class="p">,</span> <span class="mi">1</span><span class="p">,</span>
<span class="n">constrained_domain</span><span class="o">=</span><span class="n">geom</span><span class="o">.</span><span class="n">PeriodicBoundary_no_corner</span><span class="p">(</span><span class="mi">2</span><span class="p">))</span>
<span class="n">f</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">Function</span><span class="p">(</span><span class="n">FS_2d</span><span class="p">)</span>
<span class="c"># DoF that are cancelled out</span>
<span class="k">print</span> <span class="s">'2D periodic map'</span>
<span class="k">print</span> <span class="s">'original DoF ='</span><span class="p">,</span> <span class="p">(</span><span class="n">a</span> <span class="o">+</span> <span class="mi">1</span><span class="p">)</span> <span class="o">*</span> <span class="p">(</span><span class="n">b</span> <span class="o">+</span> <span class="mi">1</span><span class="p">),</span> <span class="s">';'</span><span class="p">,</span>
<span class="k">print</span> <span class="s">'actual DoF ='</span><span class="p">,</span> <span class="n">f</span><span class="o">.</span><span class="n">vector</span><span class="p">()</span><span class="o">.</span><span class="n">size</span><span class="p">(),</span> <span class="s">';'</span><span class="p">,</span>
<span class="k">print</span> <span class="s">'the excluded DoF ='</span><span class="p">,</span> <span class="p">(</span><span class="n">a</span> <span class="o">-</span> <span class="mi">1</span> <span class="o">+</span> <span class="n">b</span> <span class="o">-</span> <span class="mi">1</span><span class="p">)</span>
</pre></div>
</div>
<p>Outputs for 2D case are as below,</p>
<div class="highlight-python"><div class="highlight"><pre>2D periodic map
original DoF = 28 ; actual DoF = 21 ; the excluded DoF = 7
</pre></div>
</div>
<p>If 3D problem is considered, the code above pass with little modification,</p>
<div class="code python highlight-python"><div class="highlight"><pre><span class="c"># 3D</span>
<span class="n">a</span><span class="p">,</span> <span class="n">b</span><span class="p">,</span> <span class="n">c</span> <span class="o">=</span> <span class="mi">3</span><span class="p">,</span> <span class="mi">6</span><span class="p">,</span> <span class="mi">9</span>
<span class="n">mesh_3d</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">UnitCubeMesh</span><span class="p">(</span><span class="n">a</span><span class="p">,</span> <span class="n">b</span><span class="p">,</span> <span class="n">c</span><span class="p">)</span>
<span class="n">FS_3d</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">FunctionSpace</span><span class="p">(</span><span class="n">mesh_3d</span><span class="p">,</span> <span class="s">'CG'</span><span class="p">,</span> <span class="mi">1</span><span class="p">,</span>
<span class="n">constrained_domain</span><span class="o">=</span><span class="n">geom</span><span class="o">.</span><span class="n">PeriodicBoundary_no_corner</span><span class="p">(</span><span class="mi">3</span><span class="p">))</span>
<span class="n">f</span> <span class="o">=</span> <span class="n">geom</span><span class="o">.</span><span class="n">Function</span><span class="p">(</span><span class="n">FS_3d</span><span class="p">)</span>
<span class="c"># DoF that are cancelled out</span>
<span class="k">print</span> <span class="s">'3D periodic map'</span>
<span class="k">print</span> <span class="s">'original DoF ='</span><span class="p">,</span> <span class="p">(</span><span class="n">a</span> <span class="o">+</span> <span class="mi">1</span><span class="p">)</span> <span class="o">*</span> <span class="p">(</span><span class="n">b</span> <span class="o">+</span> <span class="mi">1</span><span class="p">)</span> <span class="o">*</span> <span class="p">(</span><span class="n">c</span> <span class="o">+</span> <span class="mi">1</span><span class="p">),</span> <span class="s">';'</span><span class="p">,</span>
<span class="k">print</span> <span class="s">'actual DoF ='</span><span class="p">,</span> <span class="n">f</span><span class="o">.</span><span class="n">vector</span><span class="p">()</span><span class="o">.</span><span class="n">size</span><span class="p">(),</span> <span class="s">';'</span><span class="p">,</span>
<span class="k">print</span> <span class="s">'the excluded DoF ='</span><span class="p">,</span> <span class="p">(</span><span class="n">a</span> <span class="o">-</span> <span class="mi">1</span> <span class="o">+</span> <span class="n">b</span> <span class="o">-</span> <span class="mi">1</span> <span class="o">+</span> <span class="n">c</span> <span class="o">-</span> <span class="mi">1</span><span class="p">)</span> <span class="o">*</span> <span class="mi">3</span> <span class="o">+</span> \
<span class="p">(</span><span class="n">a</span> <span class="o">-</span> <span class="mi">1</span><span class="p">)</span> <span class="o">*</span> <span class="p">(</span><span class="n">b</span> <span class="o">-</span> <span class="mi">1</span><span class="p">)</span> <span class="o">+</span> <span class="p">(</span><span class="n">a</span> <span class="o">-</span> <span class="mi">1</span><span class="p">)</span> <span class="o">*</span> <span class="p">(</span><span class="n">c</span> <span class="o">-</span> <span class="mi">1</span><span class="p">)</span> <span class="o">+</span> <span class="p">(</span><span class="n">b</span> <span class="o">-</span> <span class="mi">1</span><span class="p">)</span> <span class="o">*</span> <span class="p">(</span><span class="n">c</span> <span class="o">-</span> <span class="mi">1</span><span class="p">)</span>
</pre></div>
</div>
<p>Outputs for 3D case are as below,</p>
<div class="highlight-python"><div class="highlight"><pre>3D periodic map
original DoF = 280 ; actual DoF = 169 ; the excluded DoF = 111
</pre></div>
</div>
</div>
</div>
</div>
</div>
</div>
<div class="sphinxsidebar">
<div class="sphinxsidebarwrapper">
<h3><a href="index.html">Table Of Contents</a></h3>
<ul>
<li><a class="reference internal" href="#">cell_geom.py</a><ul>
<li><a class="reference internal" href="#overview">Overview</a></li>
<li><a class="reference internal" href="#inclusions">Inclusions</a><ul>
<li><a class="reference internal" href="#d-case">2D Case</a></li>
<li><a class="reference internal" href="#id1">3D Case</a></li>
</ul>
</li>
<li><a class="reference internal" href="#peirodic-boundary-condition">Peirodic Boundary Condition</a></li>
</ul>
</li>
</ul>
<h4>Previous topic</h4>
<p class="topless"><a href="index.html"
title="previous chapter">Welcome to Unit Cell Module’s documentation!</a></p>
<h4>Next topic</h4>
<p class="topless"><a href="Manual on cell_material.py.html"
title="next chapter">cell_material.py</a></p>
<h3>This Page</h3>
<ul class="this-page-menu">
<li><a href="_sources/Manual on cell_geom.py.txt"
rel="nofollow">Show Source</a></li>
</ul>
<div id="searchbox" style="display: none">
<h3>Quick search</h3>
<form class="search" action="search.html" method="get">
<input type="text" name="q" />
<input type="submit" value="Go" />
<input type="hidden" name="check_keywords" value="yes" />
<input type="hidden" name="area" value="default" />
</form>
<p class="searchtip" style="font-size: 90%">
Enter search terms or a module, class or function name.
</p>
</div>
<script type="text/javascript">$('#searchbox').show(0);</script>
</div>
</div>
<div class="clearer"></div>
</div>
<div class="related">
<h3>Navigation</h3>
<ul>
<li class="right" style="margin-right: 10px">
<a href="genindex.html" title="General Index"
>index</a></li>
<li class="right" >
<a href="py-modindex.html" title="Python Module Index"
>modules</a> |</li>
<li class="right" >
<a href="Manual on cell_material.py.html" title="cell_material.py"
>next</a> |</li>
<li class="right" >
<a href="index.html" title="Welcome to Unit Cell Module’s documentation!"
>previous</a> |</li>
<li><a href="index.html">Unit Cell Module Documentation 1.0 documentation</a> »</li>
</ul>
</div>
<div class="footer">
© Copyright 2016, Yi Hu.
Created using <a href="http://sphinx-doc.org/">Sphinx</a> 1.2.2.
</div>
</body>
</html>