This repository has been archived by the owner on Oct 6, 2018. It is now read-only.
-
Notifications
You must be signed in to change notification settings - Fork 1
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
Added a tutorial for the 2D Heat Equation
- Loading branch information
Vincenzo Pandolfo
committed
Jul 26, 2016
1 parent
2a1ef85
commit f60b67a
Showing
3 changed files
with
60 additions
and
0 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,58 @@ | ||
Solving the 2D Heat Equation using Devito | ||
========================================= | ||
|
||
Consider the 2D Heat Equation defined by | ||
|
||
.. math:: | ||
u_t = a\left(u_{xx}+u_{yy}\right) | ||
This equation can be solved numerically using Finite Difference approximations. | ||
|
||
First of all, we need to allocate the grid and set its initial condition:: | ||
|
||
from devito import TimeData | ||
u = TimeData(name = 'u', shape = 'nx, ny', time_order = 1, space_order = 2) | ||
u.data[0, :] = ui[:] | ||
|
||
TimeData is a Devito data object used to store and manage time-varying data. | ||
|
||
We initialise our grid to be of size :obj:`nx * ny` for some :obj:`nx` | ||
and :obj:`ny`. :obj:`time_order` and :obj:`space_order` represent the discretization | ||
order for time and space respectively. The initial configuration is given as a | ||
:class:`numpy.array` in :obj:`u.data`. | ||
|
||
The next step is to generate the stencil to be solved by a | ||
:class:`devito.operator.Operator`:: | ||
from devito import Operator | ||
from sympy import Eq, solve, symbols | ||
a, h, s = symbols("a h s") | ||
eqn = Eq(u.dt, a * (u.dx2 + u.dy2)) | ||
stencil = solve(eqn, u.forward)[0] | ||
|
||
The stencil is generated according to Devito conventions. It uses a sympy | ||
equation to represent the 2D Heat equation and store it in :obj:`eqn`. | ||
Devito makes easy to represent the equation by providing properties :obj:`dt`, | ||
:obj:`dx2`, and :obj:`dx2` that represent the derivatives. | ||
|
||
We then generate the stencil by solving :obj:`eqn` for :obj:`u.forward`, a | ||
symbol for the time-forward state of the function. | ||
|
||
We plug the stencil in an Operator, as shown, and define the values of the | ||
thermal conductivity :obj:`a`, the spacing between cells :obj:`h` and the | ||
temporal spacing :obj:`s`.:: | ||
|
||
op = Operator( | ||
stencils = Eq(u.forward, stencil), | ||
substitutions = {a: tc, h: dx, s: dt}, | ||
nt = timesteps, | ||
shape = (nx, ny), | ||
spc_border = 1, | ||
time_border = 1) | ||
|
||
|
||
To execute the generated Operator, we simply call :samp:`op.apply()`. The | ||
results will then be found in :obj:`u.data[1, :]`. | ||
|
||
For a complete example of this code, check file `test_diffusion.py` in the | ||
`tests` folder. |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters