I am reading a really excellent book on numerical methods for PDE called Computational Seismology, but Igel. He has a really nice treatment of finite element methods, but I was
a bit confused about how he derives the off-diagonal elements for the mass matrix
on a time-dependent problem.
So the demo problem is a 2 dimensional wave equation, where there is just 1 spatial dimension $x$ and then the time dimension $t$.
$$
\rho \partial_{t}^{2} u=\partial_{x} \mu \partial_{x} u+f
$$
The author defines a coordinate system that is local to the basis functions. So $\xi = x – x_i$, for a basis discretized into $i$ intervals. And $h_i = x_i – x_{i-1}$.
Now using linear basis functions, we get the following definition for the basis.
$$
\varphi_{i}(\xi)= \begin{cases}\frac{\xi}{h_{i-1}}+1 & \text { for }-h_{i-1} \leq \xi \leq 0 \\ 1-\frac{\xi}{h_{i}} & \text { for } 0 \leq \xi \leq h_{i} \\ 0 & \text { elsewhere }\end{cases}
$$
The part that is confusing is the definition for the off-diagonal element of the mass matrix. The author defines that as.
So the mass matrix is defined as:
$$
M_{i j}=\int_{D} \rho \varphi_{i} \varphi_{j} d x
$$
The $M_{ij}$ will correspond to $M_{i,i-1}$ for the off-diagonal element. The author defined the off-diagonal element as below.
$$
M_{i, i-1}=\rho_{i-1} \int_{-h_{i-1}}^{0}\left(\frac{\xi}{h_{i-1}}+1\right) \frac{-\xi}{h_{i-1}} d \xi=\frac{1}{6} \rho_{i-1} h_{i-1}
$$
I was not sure why the author defined the integrand as:
$$
\left(\frac{\xi}{h_{i-1}}+1\right) \frac{-\xi}{h_{i-1}}
$$
instead of defining it as below:
$$
\left(\frac{\xi}{h_{i-1}}+1\right)\left(1 – \frac{\xi}{h_{i-1}} \right)
$$
It seems like the way the author defined it is not consistent with the definition of the basis function–as having either a $\pm 1$. The numbers seem to work out correctly as Igel specified the function.
Hence I was hoping someone could just explain why that integrand is setup with $ \frac{-\xi}{h_{i-1}} $ instead of $ \left(1 – \frac{\xi}{h_{i-1}} \right) $.
Best Answer
Okay, I think I figured it out. It is a simple thing, but it often helps to be explicit about these things.
So the short answer to the question is that it has to do with the reference for the value of $\xi$. If just look at the basis function where $\xi_i$ defined with reference to its own linear basis function $\phi_i$, then the definition:
$$ \varphi_{i}(\xi)= \begin{cases}\frac{\xi}{h_{i-1}}+1 & \text { for }-h_{i-1} \leq \xi \leq 0 \\ 1-\frac{\xi}{h_{i}} & \text { for } 0 \leq \xi \leq h_{i} \\ 0 & \text { elsewhere }\end{cases} $$
is correct.
But, since the integral is using basis function $\phi_i$, as the reference, then looking at the value of $\xi_i$ with respect to basis function $\phi_{i-1}$ is where the confusion happens. So if I were looking at $\xi_{i-1}$ with respect to $\phi_{i-1}$ then the correct function to use is $1 - \frac{\xi_{i-1}}{h_{i-1}}$.
BUT, if I am trying to compute the value of the basis function for $\xi_i$, using basis function $\phi_{i-1}$, then that would be $\frac{\xi_i}{h_i}$. This is basically just computing the value of the function starting from the left side of the interval versus the right side.