35 |
\frac{\partial^2 p}{\partial x^2\hackscore{i}} |
\frac{\partial^2 p}{\partial x^2\hackscore{i}} |
36 |
\label{eqn:laplacian} |
\label{eqn:laplacian} |
37 |
\end{equation} |
\end{equation} |
38 |
For the two dimensional case in Cartesian coordinates \refEq{eqn:laplacian} |
For the two dimensional case in Cartesian coordinates \autoref{eqn:laplacian} |
39 |
becomes; |
becomes; |
40 |
\begin{equation} |
\begin{equation} |
41 |
\nabla^2 p = \frac{\partial^2 p}{\partial x^2} |
\nabla^2 p = \frac{\partial^2 p}{\partial x^2} |
66 |
\label{eqn:gradrank1} |
\label{eqn:gradrank1} |
67 |
\end{equation} |
\end{equation} |
68 |
|
|
69 |
\refEq{eqn:grad} corresponds to the Linear PDE general form value |
\autoref{eqn:grad} corresponds to the Linear PDE general form value |
70 |
$X$. Notice however that the general form contains the term $X |
$X$. Notice however that the general form contains the term $X |
71 |
\hackscore{i,j}$\footnote{This is the first derivative in the $j^{th}$ |
\hackscore{i,j}$\footnote{This is the first derivative in the $j^{th}$ |
72 |
direction for the $i^{th}$ component of the solution.}, |
direction for the $i^{th}$ component of the solution.}, |
73 |
hence for a rank 0 object there is no need to do more than calculate the |
hence for a rank 0 object there is no need to do more than calculate the |
74 |
gradient and submit it to the solver. In the case of the rank 1 or greater |
gradient and submit it to the solver. In the case of the rank 1 or greater |
75 |
object, it is nesscary to calculate the trace also. This is the sum of the |
object, it is nesscary to calculate the trace also. This is the sum of the |
76 |
diagonal in \refeq{eqn:gradrank1}. |
diagonal in \autoref{eqn:gradrank1}. |
77 |
|
|
78 |
Thus when solving for equations containing the Laplacian one of two things must |
Thus when solving for equations containing the Laplacian one of two things must |
79 |
be completed. If the object \verb p is less than rank 1 the gradient is |
be completed. If the object \verb!p! is less than rank 1 the gradient is |
80 |
calculated via; |
calculated via; |
81 |
\begin{python} |
\begin{python} |
82 |
gradient=grad(p) |
gradient=grad(p) |
83 |
\end{python} |
\end{python} |
84 |
and if the object is greater thank or equal to a rank 1 tensor, the trace of |
and if the object is greater thank or equal to a rank 1 tensor, the trace of |
85 |
the gradient is calculated. |
the gradient is calculated. |
148 |
f''(x) \approx \frac{f(x+h - 2f(x) + f(x-h)}{h^2} |
f''(x) \approx \frac{f(x+h - 2f(x) + f(x-h)}{h^2} |
149 |
\label{eqn:centdiff} |
\label{eqn:centdiff} |
150 |
\end{equation} |
\end{equation} |
151 |
substituting \refEq{eqn:centdiff} for $\frac{\partial ^2 p }{\partial t ^2}$ |
substituting \autoref{eqn:centdiff} for $\frac{\partial ^2 p }{\partial t ^2}$ |
152 |
in \refEq{eqn:acswave}; |
in \autoref{eqn:acswave}; |
153 |
\begin{equation} |
\begin{equation} |
154 |
\nabla ^2 p - \frac{1}{c^2h^2} \left[p\hackscore{(t+1)} - 2p\hackscore{(t)} + |
\nabla ^2 p - \frac{1}{c^2h^2} \left[p\hackscore{(t+1)} - 2p\hackscore{(t)} + |
155 |
p\hackscore{(t-1)} \right] |
p\hackscore{(t-1)} \right] |
206 |
suface in mayavi. The second is to take a cross section of the model. |
suface in mayavi. The second is to take a cross section of the model. |
207 |
|
|
208 |
For the later, we will require the \textit{Locator} function. |
For the later, we will require the \textit{Locator} function. |
209 |
First \verb Locator must be imported; |
First \verb!Locator! must be imported; |
210 |
\begin{python} |
\begin{python} |
211 |
from esys.escript.pdetools import Locator |
from esys.escript.pdetools import Locator |
212 |
\end{python} |
\end{python} |
222 |
cut_loc.append(xstep*i) |
cut_loc.append(xstep*i) |
223 |
src_cut.append([xstep*i,xc[1]]) |
src_cut.append([xstep*i,xc[1]]) |
224 |
\end{python} |
\end{python} |
225 |
We then submit the output to \verb Locator and finally return the appropriate |
We then submit the output to \verb!Locator! and finally return the appropriate |
226 |
values using the \verb getValue function. |
values using the \verb!getValue! function. |
227 |
\begin{python} |
\begin{python} |
228 |
src=Locator(mydomain,src_cut) |
src=Locator(mydomain,src_cut) |
229 |
src_cut=src.getValue(u) |
src_cut=src.getValue(u) |
236 |
\end{python} |
\end{python} |
237 |
\begin{figure}[h] |
\begin{figure}[h] |
238 |
\centering |
\centering |
239 |
FIXME PLEASE! |
\includegraphics[width=6in]{figures/sourceline.png} |
|
% \includegraphics[width=6in]{figures/sourceline.png} |
|
240 |
\caption{Cross section of the source function.} |
\caption{Cross section of the source function.} |
241 |
\label{fig:cxsource} |
\label{fig:cxsource} |
242 |
\end{figure} |
\end{figure} |
254 |
\begin{python} |
\begin{python} |
255 |
u_rec=rec.getValue(u) |
u_rec=rec.getValue(u) |
256 |
\end{python} |
\end{python} |
257 |
For consecutive time steps one can record the values from \verb u_rec in an |
For consecutive time steps one can record the values from \verb!u_rec! in an |
258 |
array initialised as \verb u_rec0=[] with; |
array initialised as \verb!u_rec0=[]! with; |
259 |
\begin{python} |
\begin{python} |
260 |
u_rec0.append(rec.getValue(u)) |
u_rec0.append(rec.getValue(u)) |
261 |
\end{python} |
\end{python} |
262 |
|
|
263 |
It can be useful to monitor the value at a single or multiple individual points |
It can be useful to monitor the value at a single or multiple individual points |
264 |
in the model during the modelling process. This is done using |
in the model during the modelling process. This is done using |
265 |
the \verb Locator function. |
the \verb!Locator! function. |
266 |
|
|
267 |
|
|
268 |
\section{Acceleration Solution} |
\section{Acceleration Solution} |
270 |
|
|
271 |
An alternative method is to solve for the acceleration $\frac{\partial ^2 |
An alternative method is to solve for the acceleration $\frac{\partial ^2 |
272 |
p}{\partial t^2}$ directly, and derive the displacement solution from the |
p}{\partial t^2}$ directly, and derive the displacement solution from the |
273 |
PDE solution. \refEq{eqn:waveu} is thus modified; |
PDE solution. \autoref{eqn:waveu} is thus modified; |
274 |
\begin{equation} |
\begin{equation} |
275 |
\nabla ^2 p - \frac{1}{c^2} a = 0 |
\nabla ^2 p - \frac{1}{c^2} a = 0 |
276 |
\label{eqn:wavea} |
\label{eqn:wavea} |
309 |
frequency is $f\hackscore{r} = \frac{380.0}{5} = 76.0 Hz$. This is a worst case |
frequency is $f\hackscore{r} = \frac{380.0}{5} = 76.0 Hz$. This is a worst case |
310 |
scenario with a small source and the models maximum velocity. |
scenario with a small source and the models maximum velocity. |
311 |
|
|
312 |
Furthermore, we know from \refSec{sec:nsstab}, that the spatial sampling |
Furthermore, we know from \autoref{sec:nsstab}, that the spatial sampling |
313 |
frequency must be at least twice this value to ensure stability. If we assume |
frequency must be at least twice this value to ensure stability. If we assume |
314 |
the model mesh is a square equispaced grid, |
the model mesh is a square equispaced grid, |
315 |
then the sampling interval is the side length divided by the number of samples, |
then the sampling interval is the side length divided by the number of samples, |
316 |
given by $\Delta x = \frac{1000.0m}{400} = 2.5m$ and the maximum sampling |
given by $\Delta x = \frac{1000.0m}{400} = 2.5m$ and the maximum sampling |
317 |
frequency capable at this interval is |
frequency capable at this interval is |
318 |
$f\hackscore{s}=\frac{380.0ms^{-1}}{2.5m}=152Hz$ this is just equal to the |
$f\hackscore{s}=\frac{380.0ms^{-1}}{2.5m}=152Hz$ this is just equal to the |
319 |
required rate satisfying \refeq{eqn:samptheorem}. |
required rate satisfying \autoref{eqn:samptheorem}. |
320 |
|
|
321 |
\reffig{fig:ex07sampth} depicts three examples where the grid has been |
\autoref{fig:ex07sampth} depicts three examples where the grid has been |
322 |
undersampled, sampled correctly, and over sampled. The grids used had |
undersampled, sampled correctly, and over sampled. The grids used had |
323 |
200, 400 and 800 nodes per side respectively. Obviously, the oversampled grid |
200, 400 and 800 nodes per side respectively. Obviously, the oversampled grid |
324 |
retains the best resolution of the modelled wave. |
retains the best resolution of the modelled wave. |