/[escript]/trunk/doc/cookbook/example03.tex
ViewVC logotype

Diff of /trunk/doc/cookbook/example03.tex

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 2982 by caltinay, Wed Mar 10 04:27:47 2010 UTC revision 3308 by jfenwick, Tue Oct 26 03:24:54 2010 UTC
# Line 40  To expand upon our 1D problem, the domai Line 40  To expand upon our 1D problem, the domai
40  have solved a two dimensional problem already but essentially ignored the second  have solved a two dimensional problem already but essentially ignored the second
41  dimension. In our definition phase we create a square domain in $x$ and  dimension. In our definition phase we create a square domain in $x$ and
42  $y$\footnote{In \esc the notation  $y$\footnote{In \esc the notation
43  $x\hackscore{0}$ and $x\hackscore{1}$ is used for $x$ and $y$, respectively.}  $x_{0}$ and $x_{1}$ is used for $x$ and $y$, respectively.}
44  that is $600$ meters along each side \reffig{fig:twodhdmodel}. Now we set the  that is $600$ meters along each side \reffig{fig:twodhdmodel}. Now we set the
45  number of discrete spatial cells to $150$ in both directions and the radius of  number of discrete spatial cells to $150$ in both directions and the radius of
46  the intrusion to $200$ meters with the centre located at the $300$ meter mark  the intrusion to $200$ meters with the centre located at the $300$ meter mark
# Line 72  In fact, the temperature model \refEq{eq Line 72  In fact, the temperature model \refEq{eq
72  1D case with a constant material parameter only. For the more general case  1D case with a constant material parameter only. For the more general case
73  examined in this chapter, the correct model equation is  examined in this chapter, the correct model equation is
74  \begin{equation}  \begin{equation}
75  \rho c\hackscore p \frac{\partial T}{\partial t} -  \frac{\partial }{\partial x}  \rho c_p \frac{\partial T}{\partial t} -  \frac{\partial }{\partial x}
76  \kappa \frac{\partial T}{\partial x} -  \frac{\partial }{\partial y} \kappa  \kappa \frac{\partial T}{\partial x} -  \frac{\partial }{\partial y} \kappa
77  \frac{\partial T}{\partial y} = q\hackscore H  \frac{\partial T}{\partial y} = q_H
78  \label{eqn:hd2}  \label{eqn:hd2}
79  \end{equation}  \end{equation}
80  Notice that for the 1D case we have $\frac{\partial T}{\partial y}=0$ and  Notice that for the 1D case we have $\frac{\partial T}{\partial y}=0$ and
# Line 84  with a simplified model equation for thi Line 84  with a simplified model equation for thi
84  to write this equation using the $\nabla$ notation as we have already seen in  to write this equation using the $\nabla$ notation as we have already seen in
85  \refEq{eqn:commonform nabla};  \refEq{eqn:commonform nabla};
86  \begin{equation}\label{eqn:Tform nabla}  \begin{equation}\label{eqn:Tform nabla}
87  \rho c\hackscore p \frac{\partial T}{\partial t}  \rho c_p \frac{\partial T}{\partial t}
88  -\nabla \cdot \kappa \nabla T = q\hackscore H  -\nabla \cdot \kappa \nabla T = q_H
89  \end{equation}  \end{equation}
90  We can easily apply the backward Euler scheme as before to end up with  We can easily apply the backward Euler scheme as before to end up with
91  \begin{equation}  \begin{equation}
92  \frac{\rho c\hackscore p}{h} T^{(n)} -\nabla \cdot \kappa \nabla T^{(n)}  =  \frac{\rho c_p}{h} T^{(n)} -\nabla \cdot \kappa \nabla T^{(n)}  =
93  q\hackscore H +  \frac{\rho c\hackscore p}{h} T^{(n-1)}  q_H +  \frac{\rho c_p}{h} T^{(n-1)}
94  \label{eqn:hdgenf2}  \label{eqn:hdgenf2}
95  \end{equation}  \end{equation}
96  which is very similar to \refEq{eqn:hdgenf} used to define the temperature  which is very similar to \refEq{eqn:hdgenf} used to define the temperature
# Line 100  $\nabla \cdot \kappa \nabla T^{(n)}$. Un Line 100  $\nabla \cdot \kappa \nabla T^{(n)}$. Un
100  we need to revisit the \esc PDE form as shown in \ref{eqn:commonform2D}.  we need to revisit the \esc PDE form as shown in \ref{eqn:commonform2D}.
101  For the 2D case with variable PDE coefficients the form needs to be read as  For the 2D case with variable PDE coefficients the form needs to be read as
102  \begin{equation}\label{eqn:commonform2D 2}  \begin{equation}\label{eqn:commonform2D 2}
103  -\frac{\partial }{\partial x} A\hackscore{00}\frac{\partial u}{\partial x}  -\frac{\partial }{\partial x} A_{00}\frac{\partial u}{\partial x}
104  -\frac{\partial }{\partial x} A\hackscore{01}\frac{\partial u}{\partial y}  -\frac{\partial }{\partial x} A_{01}\frac{\partial u}{\partial y}
105  -\frac{\partial }{\partial y} A\hackscore{10}\frac{\partial u}{\partial x}  -\frac{\partial }{\partial y} A_{10}\frac{\partial u}{\partial x}
106  -\frac{\partial }{\partial x} A\hackscore{11}\frac{\partial u}{\partial y}  -\frac{\partial }{\partial x} A_{11}\frac{\partial u}{\partial y}
107  + Du = f  + Du = f
108  \end{equation}  \end{equation}
109  So besides the settings $u=T^{(n)}$, $D = \frac{\rho c \hackscore{p}}{h}$ and  So besides the settings $u=T^{(n)}$, $D = \frac{\rho c _{p}}{h}$ and
110  $f = q \hackscore{H} + \frac{\rho c\hackscore p}{h} T^{(n-1)}$ as we have used  $f = q _{H} + \frac{\rho c_p}{h} T^{(n-1)}$ as we have used
111  before (see \refEq{ESCRIPT SET}) we need to set  before (see \refEq{ESCRIPT SET}) we need to set
112  \begin{equation}\label{eqn: kappa general}  \begin{equation}\label{eqn: kappa general}
113  A\hackscore{00}=A\hackscore{11}=\kappa; A\hackscore{01}=A\hackscore{10}=0  A_{00}=A_{11}=\kappa; A_{01}=A_{10}=0
114  \end{equation}  \end{equation}
115  The fundamental difference to the 1D case is that $A\hackscore{11}$ is not set  The fundamental difference to the 1D case is that $A_{11}$ is not set
116  to zero but $\kappa$,  to zero but $\kappa$,
117  which brings in the second dimension. It is important to note that the  which brings in the second dimension. It is important to note that the
118  coefficients of the PDE may depend on their location in the domain which does  coefficients of the PDE may depend on their location in the domain which does
# Line 146  we set $\kappa$. As pointed out \refEq{N Line 146  we set $\kappa$. As pointed out \refEq{N
146  the boundary conditions. In our case these assumptions translate to;  the boundary conditions. In our case these assumptions translate to;
147  \begin{equation}  \begin{equation}
148  -n \cdot \kappa \nabla T^{(n)} =  -n \cdot \kappa \nabla T^{(n)} =
149  -n\hackscore{0} \cdot \kappa \frac{\partial T^{(n)}}{\partial x} -  -n_{0} \cdot \kappa \frac{\partial T^{(n)}}{\partial x} -
150  n\hackscore{1} \cdot  \kappa \frac{\partial T^{(n)}}{\partial y} = 0  n_{1} \cdot  \kappa \frac{\partial T^{(n)}}{\partial y} = 0
151  \label{eq:hom flux}  \label{eq:hom flux}
152  \end{equation}  \end{equation}
153  which sets the normal component of the heat flux $- \kappa \cdot (\frac{\partial  which sets the normal component of the heat flux $- \kappa \cdot (\frac{\partial
154  T^{(n)}}{\partial x}, \frac{\partial T^{(n)}}{\partial y})$ to zero. This means  T^{(n)}}{\partial x}, \frac{\partial T^{(n)}}{\partial y})$ to zero. This means
155  that the region is insulated which is what we want.  that the region is insulated which is what we want.
156  On the left and right face of the domain where we have  On the left and right face of the domain where we have
157  $(n\hackscore{0},n\hackscore{1} ) = (\mp 1,0)$  $(n_{0},n_{1} ) = (\mp 1,0)$
158  this means $\frac{\partial T^{(n)}}{\partial x}=0$ and on the top and bottom on  this means $\frac{\partial T^{(n)}}{\partial x}=0$ and on the top and bottom on
159  the domain  the domain
160  where we have  $(n\hackscore{0},n\hackscore{1} ) = (\pm 1,0)$ this is  where we have  $(n_{0},n_{1} ) = (\pm 1,0)$ this is
161  $\frac{\partial T^{(n)}}{\partial y}=0$.  $\frac{\partial T^{(n)}}{\partial y}=0$.
162    
163  \section{Setting variable PDE Coefficients}  \section{Setting variable PDE Coefficients}
164  Now we need to look into the problem of how we define the material coefficients  Now we need to look into the problem of how we define the material coefficients
165  $\kappa$ and $\rho c\hackscore p$ depending on their location in the domain.  $\kappa$ and $\rho c_p$ depending on their location in the domain.
166  We can make use of the technique used in the granite block example in  We can make use of the technique used in the granite block example in
167  \refSec{Sec:1DHDv00}  \refSec{Sec:1DHDv00}
168  to set up the initial temperature. However,  to set up the initial temperature. However,
# Line 172  curved interface between the two sub-dom Line 172  curved interface between the two sub-dom
172  Prior to setting up the PDE, the interface between the two materials must be  Prior to setting up the PDE, the interface between the two materials must be
173  established.  established.
174  The distance $s\ge 0$ between two points $[x,y]$ and  The distance $s\ge 0$ between two points $[x,y]$ and
175  $[x\hackscore{0},y\hackscore{0}]$ in Cartesian coordinates is defined as:  $[x_{0},y_{0}]$ in Cartesian coordinates is defined as:
176  \begin{equation}  \begin{equation}
177   (x-x\hackscore{0})^{2}+(y-y\hackscore{0})^{2} = s^{2}   (x-x_{0})^{2}+(y-y_{0})^{2} = s^{2}
178  \end{equation}  \end{equation}
179  If we define the point $[x\hackscore{0},y\hackscore{0}]$ as $ic$ which denotes  If we define the point $[x_{0},y_{0}]$ as $ic$ which denotes
180  the centre of the semi-circle of our intrusion, then for all the points $[x,y]$  the centre of the semi-circle of our intrusion, then for all the points $[x,y]$
181  in our model we can calculate a distance to $ic$.  in our model we can calculate a distance to $ic$.
182  All the points that fall within the given radius $r$ of our intrusion will have  All the points that fall within the given radius $r$ of our intrusion will have
# Line 204  mypde.setValue(A=kappa*kronecker(model)) Line 204  mypde.setValue(A=kappa*kronecker(model))
204  Notice that we are using the sample points of the \verb|Function| function space  Notice that we are using the sample points of the \verb|Function| function space
205  as expected for the PDE coefficient \verb|A|\footnote{For the experienced user: use  as expected for the PDE coefficient \verb|A|\footnote{For the experienced user: use
206  \texttt{x=mypde.getFunctionSpace("A").getX()}.}.  \texttt{x=mypde.getFunctionSpace("A").getX()}.}.
207  The corresponding statements are used to set $\rho c\hackscore p$.  The corresponding statements are used to set $\rho c_p$.
208    
209  Our PDE has now been properly established. The initial conditions for  Our PDE has now been properly established. The initial conditions for
210  temperature are set out in a similar manner:  temperature are set out in a similar manner:

Legend:
Removed from v.2982  
changed lines
  Added in v.3308

  ViewVC Help
Powered by ViewVC 1.1.26