39 |
$-(A\hackscore{jl} u\hackscore{,l})\hackscore{,j}+D u =Y$} |
$-(A\hackscore{jl} u\hackscore{,l})\hackscore{,j}+D u =Y$} |
40 |
is described by; |
is described by; |
41 |
\begin{equation}\label{eqn:commonform nabla} |
\begin{equation}\label{eqn:commonform nabla} |
42 |
-\nabla.(A.\nabla u) + Du = f |
-\nabla\cdot(A\cdot\nabla u) + Du = f |
43 |
\end{equation} |
\end{equation} |
44 |
where $A$, $D$ and $f$ are known values. The symbol $\nabla$ which is called the \textit{Nabla operator} or \textit{del operator} represents |
where $A$, $D$ and $f$ are known values. The symbol $\nabla$ which is called the \textit{Nabla operator} or \textit{del operator} represents |
45 |
the spatial derivative of its subject - in this case $u$. Lets assume for a moment that we deal with a one-dimensional problem then ; |
the spatial derivative of its subject - in this case $u$. Lets assume for a moment that we deal with a one-dimensional problem then ; |
80 |
\frac{\rho c\hackscore p}{h} T^{(n)} - \kappa \frac{\partial^2 T^{(n)}}{\partial x^2} = q\hackscore H + \frac{\rho c\hackscore p}{h} T^{(n-1)} |
\frac{\rho c\hackscore p}{h} T^{(n)} - \kappa \frac{\partial^2 T^{(n)}}{\partial x^2} = q\hackscore H + \frac{\rho c\hackscore p}{h} T^{(n-1)} |
81 |
\label{eqn:hdgenf} |
\label{eqn:hdgenf} |
82 |
\end{equation} |
\end{equation} |
83 |
The PDE is now in a form that statisfies \refEq{eqn:commonform nabla} which is required for \ESCRIPT to solve our PDE. This can be done by generating a solution for sucessive increments in the time nodes $t^{(n)}$ where |
The PDE is now in a form that satisfies \refEq{eqn:commonform nabla} which is required for \ESCRIPT to solve our PDE. This can be done by generating a solution for successive increments in the time nodes $t^{(n)}$ where |
84 |
$t^{(0)}=0$ and $t^{(n)}=t^{(n-1)}+h$ where $h>0$ is the step size and assumed to be constant. |
$t^{(0)}=0$ and $t^{(n)}=t^{(n-1)}+h$ where $h>0$ is the step size and assumed to be constant. |
85 |
In the following the upper index ${(n)}$ refers to a value at time $t^{(n)}$. Finally, by comparing \ref{eqn:hdgenf} with \ref{eqn:commonform} it can be seen that; |
In the following the upper index ${(n)}$ refers to a value at time $t^{(n)}$. Finally, by comparing \ref{eqn:hdgenf} with \ref{eqn:commonform} it can be seen that; |
86 |
\begin{equation} |
\begin{equation} |
96 |
\subsection{Boundary Conditions} |
\subsection{Boundary Conditions} |
97 |
With the PDE sufficiently modified, consideration must now be given to the boundary conditions of our model. Typically there are two main types of boundary conditions known as Neumann and Dirichlet\footnote{More information on Boundary Conditions is available at Wikipedia \url{http://en.wikipedia.org/wiki/Boundary_conditions}}. In this example, we have utilised both conditions. Dirichlet is conceptually simpler and is used to prescribe a known value to the model on its boundary. This is like holding a snake by the tail; we know where the tail will be as we hold it however, we have no control over the rest of the snake. Dirichlet boundary conditions exist where we have applied our heat source. As the heat source is a constant, we can simulate its presence on that boundary. This is done by continuously resetting the temperature of the boundary, so that is is the same as the heat source. |
With the PDE sufficiently modified, consideration must now be given to the boundary conditions of our model. Typically there are two main types of boundary conditions known as Neumann and Dirichlet\footnote{More information on Boundary Conditions is available at Wikipedia \url{http://en.wikipedia.org/wiki/Boundary_conditions}}. In this example, we have utilised both conditions. Dirichlet is conceptually simpler and is used to prescribe a known value to the model on its boundary. This is like holding a snake by the tail; we know where the tail will be as we hold it however, we have no control over the rest of the snake. Dirichlet boundary conditions exist where we have applied our heat source. As the heat source is a constant, we can simulate its presence on that boundary. This is done by continuously resetting the temperature of the boundary, so that is is the same as the heat source. |
98 |
|
|
99 |
Neumann boundary conditions describe the radiation or flux normal to the boundayor surface. This aptly describes our insulation conditions as we do not want to exert a constant temperature as with the heat source. However, we do want to prevent any loss of energy from the system. These natural boundary conditions can be described by specifying a radiation condition which prescribes the normal component of the flux $\kappa T\hackscore{,i}$ to be proportional |
Neumann boundary conditions describe the radiation or flux normal to the boundary surface. This aptly describes our insulation conditions as we do not want to exert a constant temperature as with the heat source. However, we do want to prevent any loss of energy from the system. These natural boundary conditions can be described by specifying a radiation condition which prescribes the normal component of the flux $\kappa T\hackscore{,i}$ to be proportional |
100 |
to the difference of the current temperature to the surrounding temperature $T\hackscore{ref}$; in general terms this is; |
to the difference of the current temperature to the surrounding temperature $T\hackscore{ref}$; in general terms this is; |
101 |
\begin{equation} |
\begin{equation} |
102 |
\kappa T\hackscore{,i} n\hackscore i = \eta (T\hackscore{ref}-T) |
\kappa T\hackscore{,i} n\hackscore i = \eta (T\hackscore{ref}-T) |
231 |
plx = np.array(plx) #convert to tuple to numpy array |
plx = np.array(plx) #convert to tuple to numpy array |
232 |
plx = plx[:,0] #extract x locations |
plx = plx[:,0] #extract x locations |
233 |
\end{verbatim} |
\end{verbatim} |
234 |
As there are two solution outputs, we will generate two plots and save each to a file for every time step in the solution. The following is appended to the end of the \verb while loop and creates two figures. The first figure is for the temperature distribution, and the second the total temperature in the bar. Both cases are similar with a few minor changes for scale and labeling. We start by converting the solution to a tuple and then plotting this against our \textit{x coordinates} \verb plx from before. The axis is then standardised and a title applied. The figure is then saved to a *.png file and cleared for the following iteration. |
As there are two solution outputs, we will generate two plots and save each to a file for every time step in the solution. The following is appended to the end of the \verb while loop and creates two figures. The first figure is for the temperature distribution, and the second the total temperature in the bar. Both cases are similar with a few minor changes for scale and labelling. We start by converting the solution to a tuple and then plotting this against our \textit{x coordinates} \verb plx from before. The axis is then standardised and a title applied. The figure is then saved to a *.png file and cleared for the following iteration. |
235 |
\begin{verbatim} |
\begin{verbatim} |
236 |
#establish figure 1 for temperature vs x plots |
#establish figure 1 for temperature vs x plots |
237 |
tempT = T.toListOfTuples(scalarastuple=False) |
tempT = T.toListOfTuples(scalarastuple=False) |
261 |
\end{center} |
\end{center} |
262 |
\end{figure} |
\end{figure} |
263 |
|
|
264 |
|
\subsubsection{Parallel scripts (MPI)} |
265 |
|
In some of the example files for this cookbook the plot part of the script looks a little different. |
266 |
|
For example, |
267 |
|
\begin{verbatim} |
268 |
|
pl.title("Total temperature accross Rod") |
269 |
|
if getMPIRankWorld() == 0: |
270 |
|
pl.savefig(os.path.join(save_path+"/totT","ttrodpyplot%03d.png")%i) |
271 |
|
pl.clf() |
272 |
|
\end{verbatim} |
273 |
|
|
274 |
|
The additional \verb if statement is not necessary for normal desktop use. |
275 |
|
It becomes important for scripts run on parallel computers. |
276 |
|
Its purpose is to ensure that only one copy of the file is written. |
277 |
|
For more details on writing scripts for parallel computers please consult the \emph{user's guide}. |
278 |
|
|
279 |
\subsection{Make a video} |
\subsection{Make a video} |
280 |
Our saved plots from the previous section can be cast into a video using the following command appended to the end of the script. \verb mencoder is linux only however, and other platform users will need to use an alternative video encoder. |
Our saved plots from the previous section can be cast into a video using the following command appended to the end of the script. \verb mencoder is linux only however, and other platform users will need to use an alternative video encoder. |
281 |
\begin{verbatim} |
\begin{verbatim} |
282 |
# compile the *.png files to create two *.avi videos that show T change |
# compile the *.png files to create two *.avi videos that show T change |
283 |
# with time. This opperation uses linux mencoder. For other operating |
# with time. This operation uses linux mencoder. For other operating |
284 |
# systems it is possible to use your favourite video compiler to |
# systems it is possible to use your favourite video compiler to |
285 |
# convert image files to videos. |
# convert image files to videos. |
286 |
|
|