demo_automatic_chapman_enskog_analysis.ipynb 18.6 KB
Newer Older
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
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import sympy as sp\n",
    "from lbmpy.chapman_enskog import ChapmanEnskogAnalysis\n",
    "from lbmpy.creationfunctions import create_lb_method\n",
    "sp.init_printing()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Demo: Automatic Chapman Enskog Analysis \n",
    "\n",
    "\n",
    "First, we create a SRT lattice Boltzmann method. It is defined as the set of moments, together with one relaxation rate per moment."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "\n",
       "        <table style=\"border:none; width: 100%\">\n",
       "            <tr style=\"border:none\">\n",
       "                <th style=\"border:none\" >Moment</th>\n",
       "                <th style=\"border:none\" >Eq. Value </th>\n",
       "                <th style=\"border:none\" >Relaxation Rate</th>\n",
       "            </tr>\n",
       "            <tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$1$</td>\n",
       "                            <td style=\"border:none\">$\\rho$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{0}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$x$</td>\n",
       "                            <td style=\"border:none\">$u_{0}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{1}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$y$</td>\n",
       "                            <td style=\"border:none\">$u_{1}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{1}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$z$</td>\n",
       "                            <td style=\"border:none\">$u_{2}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{1}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$x^{2}$</td>\n",
       "                            <td style=\"border:none\">$\\frac{\\rho}{3} + u_{0}^{2}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{0}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$y^{2}$</td>\n",
       "                            <td style=\"border:none\">$\\frac{\\rho}{3} + u_{1}^{2}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{0}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$z^{2}$</td>\n",
       "                            <td style=\"border:none\">$\\frac{\\rho}{3} + u_{2}^{2}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{0}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$x y$</td>\n",
       "                            <td style=\"border:none\">$u_{0} u_{1}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{0}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$x z$</td>\n",
       "                            <td style=\"border:none\">$u_{0} u_{2}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{0}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$y z$</td>\n",
       "                            <td style=\"border:none\">$u_{1} u_{2}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{0}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$x^{2} y$</td>\n",
       "                            <td style=\"border:none\">$\\frac{u_{1}}{3}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{1}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$x^{2} z$</td>\n",
       "                            <td style=\"border:none\">$\\frac{u_{2}}{3}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{1}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$x y^{2}$</td>\n",
       "                            <td style=\"border:none\">$\\frac{u_{0}}{3}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{1}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$x z^{2}$</td>\n",
       "                            <td style=\"border:none\">$\\frac{u_{0}}{3}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{1}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$y^{2} z$</td>\n",
       "                            <td style=\"border:none\">$\\frac{u_{2}}{3}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{1}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$y z^{2}$</td>\n",
       "                            <td style=\"border:none\">$\\frac{u_{1}}{3}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{1}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$x^{2} y^{2}$</td>\n",
       "                            <td style=\"border:none\">$\\frac{\\rho}{9} + \\frac{u_{0}^{2}}{3} + \\frac{u_{1}^{2}}{3}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{0}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$x^{2} z^{2}$</td>\n",
       "                            <td style=\"border:none\">$\\frac{\\rho}{9} + \\frac{u_{0}^{2}}{3} + \\frac{u_{2}^{2}}{3}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{0}$</td>\n",
       "                         </tr>\n",
       "<tr style=\"border:none\">\n",
       "                            <td style=\"border:none\">$y^{2} z^{2}$</td>\n",
       "                            <td style=\"border:none\">$\\frac{\\rho}{9} + \\frac{u_{1}^{2}}{3} + \\frac{u_{2}^{2}}{3}$</td>\n",
       "                            <td style=\"border:none\">$\\omega_{0}$</td>\n",
       "                         </tr>\n",
       "\n",
       "        </table>\n",
       "        "
      ],
      "text/plain": [
       "<lbmpy.methods.momentbased.MomentBasedLbMethod at 0x7f9bbc455a90>"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "method = create_lb_method(method='trt', stencil='D3Q19', compressible=False)\n",
    "method"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, the Chapman Enskog analysis object is created. This may take a while..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "analysis = ChapmanEnskogAnalysis(method)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This object now information about the method, e.g. the relation of relaxation rate to viscosities, if the method approximates the compressible or incompressible continuity equation ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "False"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "analysis.compressible"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAABwAAAAmBAMAAAAlwuZsAAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEO+Zu3ZEIqsyiWbdVM2WrI0hAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA2UlEQVQoFWMQUlJmgAIjJUUGITCbzbQTRMO4ixmSkLjcBxjqExCyTBsY8i8guPkCDPEFCK4NA8MbIA9m1GEGzp8ILuMXBpkNCC7v91BLEA+qmPkBmAPjMgWgcOsnoHBzoDyYRRRw/4PAB6ABMO+jWAQzmCxZ9lVLLoAMgJpczMD7G4kbf4HhDxI3KYHxHxKXAVUxA4PMAmTZqsMCyFwGLjUULoMSSBpqbzkDgz/IHVDufwEUrgoDw/kCsKy8iysDw24G3q8MDCEuHxnY0tIZGNiMTRIYGMrSEgBPnDRDf1mzMgAAAABJRU5ErkJggg==\n",
      "text/latex": [
       "$$\\left [ \\frac{\\rho}{3}\\right ]$$"
      ],
      "text/plain": [
       "⎡ρ⎤\n",
       "⎢─⎥\n",
       "⎣3⎦"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "analysis.pressure_equation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAEwAAAAuBAMAAABwheJJAAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEM3dMrsi70RmdpmJVKuALsSwAAAACXBIWXMAAA7EAAAOxAGVKw4bAAABkUlEQVQ4Ed2UP0vDQBiHf5L+TYIpBXGzHYrgVMHFzXwDO+jkUnBwcag4CILQ0bEfwCGbgw5xE+rgVmg7ZBdE0VGQqpQKYrw0196b3NHWtTe973O/3OVC7gGmjVRz156WYfMrMAczxJ6Bkxlix8C7PT1XzkVjFfbIQV313FZO0GyQKDh6+16wsMr8EGIEteld4JHAYak9EXIOFJGtbsIKtshvBKM0nG+SFBwkXJgLX0h7FAO6S3sXuofDRB9ag2LgAamqIDvYdrGa6EFzBGQV20OvCtLqnP22Gmy1tCMgq/a7nRcClj7sy1Owd7PqhAJl3/+MgGHDTypPRMmV/N2igbAz725UeA7Y+FfwJ405OKh8hO61LUOJ5OtGUYIy2EOyItM4yfTjRNkna0och1apfQsE26p1wvOFNXaduU7iS5C+0INxxHWisM4oaa0jO+A6UVhnFNOYbr5DnWDCXVyssdVCnSisM1oMzCFGjetEts44hjcse1wnknVECmbzFaFOZOuQmChl64g5Wk04KY393zp/Qi2ByxcRfecAAAAASUVORK5CYII=\n",
      "text/latex": [
       "$$- \\frac{\\omega_{0} - 2}{6 \\omega_{0}}$$"
      ],
      "text/plain": [
       "-(ω₀ - 2) \n",
       "──────────\n",
       "   6⋅ω₀   "
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "analysis.get_kinematic_viscosity()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAFwAAAAuBAMAAABXK2OhAAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEM3dMlTvq5l2Zokiu0Rn3bgMAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAB3UlEQVRIDd2UzysEYRjHv+/sj9lfMjmIE5bauFAOwmWPSlknFwfbFoUNJynKnKwSbaJolT0pJe1Fthx2T3JwcHA3/wErCmG872rWO+/MaPfg4j3M+zzf72eeed/nnV6ginG1eWtDjSk2IpXIAgqaaJGlVQdcVuCNizhw6IB7eyA/VY9LpZpwWthVqr46Je9zNeF9Vtpxq4BPrQnfs6Gdq3tUXFtfcOo7hoCD6vHQ4EWqx4Kn5iNZi8gESdd1K26L/oXoSddU9X/iDb1sdAD0LHT9gXYkGA63J8LhBRrWl0XuYd+w/9mZn702XawYSWWrZOO8aIjmmczg0rCIETRmQ0kzZmRyHLJqJMa8CCwZsXmub0WQnZJpPAIFk1BJWm7gFy9O8kEv2mgF4QMvrf7JCyxO0OpaHQ1Cc4JFVyK9CRpbSXf2lKrkVbQ6Mf4salLaUxhVmXqGTMbkBvaHxbUDxzuFEQ0ScBLIedhPyg3Z5hbHajPYNZZzaRBsSXi9XCniA3YRyrVEMcmVDiQRK3J5OTwqut/rou4XSEVqDnC2v410cel3uK00ZIOzaxPLM4gp6Of9/LrGp+U4uL4Fctfhn0qzT09bfGeBrt1yis40rJ35BQb8qtj3X3Hk8wq+APCLgMe5ZQ0JAAAAAElFTkSuQmCC\n",
      "text/latex": [
       "$$- \\frac{1}{9} + \\frac{2}{9 \\omega_{0}}$$"
      ],
      "text/plain": [
       "  1    2  \n",
       "- ─ + ────\n",
       "  9   9⋅ω₀"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "analysis.get_bulk_viscosity()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "But also details of the analysis are available:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}{\\partial_{t} \\rho} + {\\partial_{0} u_{0}} + {\\partial_{1} u_{1}} + {\\partial_{2} u_{2}}\\\\- \\frac{\\epsilon \\omega_{0}}{2} {\\partial_{2} {\\Pi_{02}^{(1)}}} - \\frac{\\epsilon \\omega_{0}}{2} {\\partial_{1} {\\Pi_{01}^{(1)}}} - \\frac{\\epsilon \\omega_{0}}{2} {\\partial_{0} {\\Pi_{00}^{(1)}}} + \\epsilon {\\partial_{2} {\\Pi_{02}^{(1)}}} + \\epsilon {\\partial_{1} {\\Pi_{01}^{(1)}}} + \\epsilon {\\partial_{0} {\\Pi_{00}^{(1)}}} + \\frac{{\\partial_{0} \\rho}}{3} + {\\partial_{t} u_{0}} + {\\partial_{0} (u_{0}^{2}) } + {\\partial_{1} (u_{0} u_{1}) } + {\\partial_{2} (u_{0} u_{2}) }\\\\- \\frac{\\epsilon \\omega_{0}}{2} {\\partial_{2} {\\Pi_{12}^{(1)}}} - \\frac{\\epsilon \\omega_{0}}{2} {\\partial_{1} {\\Pi_{11}^{(1)}}} - \\frac{\\epsilon \\omega_{0}}{2} {\\partial_{0} {\\Pi_{01}^{(1)}}} + \\epsilon {\\partial_{2} {\\Pi_{12}^{(1)}}} + \\epsilon {\\partial_{1} {\\Pi_{11}^{(1)}}} + \\epsilon {\\partial_{0} {\\Pi_{01}^{(1)}}} + \\frac{{\\partial_{1} \\rho}}{3} + {\\partial_{t} u_{1}} + {\\partial_{1} (u_{1}^{2}) } + {\\partial_{0} (u_{0} u_{1}) } + {\\partial_{2} (u_{1} u_{2}) }\\\\- \\frac{\\epsilon \\omega_{0}}{2} {\\partial_{2} {\\Pi_{22}^{(1)}}} - \\frac{\\epsilon \\omega_{0}}{2} {\\partial_{1} {\\Pi_{12}^{(1)}}} - \\frac{\\epsilon \\omega_{0}}{2} {\\partial_{0} {\\Pi_{02}^{(1)}}} + \\epsilon {\\partial_{2} {\\Pi_{22}^{(1)}}} + \\epsilon {\\partial_{1} {\\Pi_{12}^{(1)}}} + \\epsilon {\\partial_{0} {\\Pi_{02}^{(1)}}} + \\frac{{\\partial_{2} \\rho}}{3} + {\\partial_{t} u_{2}} + {\\partial_{2} (u_{2}^{2}) } + {\\partial_{0} (u_{0} u_{2}) } + {\\partial_{1} (u_{1} u_{2}) }\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "⎡                                                                             \n",
       "⎢                                                                             \n",
       "⎢  ε⋅ω₀⋅D(\\Pi_(1)_(1, 0, 1))   ε⋅ω₀⋅D(\\Pi_(1)_(1, 1, 0))   ε⋅ω₀⋅D(\\Pi_(1)_(2, \n",
       "⎢- ───────────────────────── - ───────────────────────── - ───────────────────\n",
       "⎢              2                           2                           2      \n",
       "⎢                                                                             \n",
       "⎢  ε⋅ω₀⋅D(\\Pi_(1)_(0, 1, 1))   ε⋅ω₀⋅D(\\Pi_(1)_(0, 2, 0))   ε⋅ω₀⋅D(\\Pi_(1)_(1, \n",
       "⎢- ───────────────────────── - ───────────────────────── - ───────────────────\n",
       "⎢              2                           2                           2      \n",
       "⎢                                                                             \n",
       "⎢  ε⋅ω₀⋅D(\\Pi_(1)_(0, 0, 2))   ε⋅ω₀⋅D(\\Pi_(1)_(0, 1, 1))   ε⋅ω₀⋅D(\\Pi_(1)_(1, \n",
       "⎢- ───────────────────────── - ───────────────────────── - ───────────────────\n",
       "⎣              2                           2                           2      \n",
       "\n",
       "             D(rho) + D(u_0) + D(u_1) + D(u_2)                                \n",
       "                                                                              \n",
       "0, 0))                                                                        \n",
       "────── + ε⋅D(\\Pi_(1)_(1, 0, 1)) + ε⋅D(\\Pi_(1)_(1, 1, 0)) + ε⋅D(\\Pi_(1)_(2, 0, \n",
       "                                                                              \n",
       "                                                                              \n",
       "1, 0))                                                                        \n",
       "────── + ε⋅D(\\Pi_(1)_(0, 1, 1)) + ε⋅D(\\Pi_(1)_(0, 2, 0)) + ε⋅D(\\Pi_(1)_(1, 1, \n",
       "                                                                              \n",
       "                                                                              \n",
       "0, 1))                                                                        \n",
       "────── + ε⋅D(\\Pi_(1)_(0, 0, 2)) + ε⋅D(\\Pi_(1)_(0, 1, 1)) + ε⋅D(\\Pi_(1)_(1, 0, \n",
       "                                                                              \n",
       "\n",
       "                                                           ⎤\n",
       "                                                           ⎥\n",
       "      D(rho)                                               ⎥\n",
       "0)) + ────── + D(u_0) + D(u_0**2) + D(u_0*u_1) + D(u_0*u_2)⎥\n",
       "        3                                                  ⎥\n",
       "                                                           ⎥\n",
       "      D(rho)                                               ⎥\n",
       "0)) + ────── + D(u_1) + D(u_1**2) + D(u_0*u_1) + D(u_1*u_2)⎥\n",
       "        3                                                  ⎥\n",
       "                                                           ⎥\n",
       "      D(rho)                                               ⎥\n",
       "1)) + ────── + D(u_2) + D(u_2**2) + D(u_0*u_2) + D(u_1*u_2)⎥\n",
       "        3                                                  ⎦"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sp.Matrix(analysis.get_macroscopic_equations())"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
Frederik Hennig's avatar
Frederik Hennig committed
361
   "version": "3.8.2"
362
363
364
365
366
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}