Skip to content
Snippets Groups Projects
Commit e5d7ddd3 authored by Erik Strand's avatar Erik Strand
Browse files

Start pset 6

parent 72893da0
No related branches found
No related tags found
No related merge requests found
from sympy import *
from sympy.assumptions.assume import global_assumptions
x = Symbol("x", real=True)
h = Symbol("h", real=True, positive=True)
# not sure why I did this
def linear_hats():
n = 5
A = zeros(5, 5)
for i in range(0, n):
if i - 1 >= 0:
A[i, i - 1] = h / 6
A[i, i] = 2 * h / 3
if i + 1 < n:
A[i, i + 1] = h / 6
print(latex(A))
B = zeros(5, 5)
for i in range(0, n):
if i - 1 >= 0:
B[i, i - 1] = 1 / h
B[i, i] = -2 / h
if i + 1 < n:
B[i, i + 1] = 1 / h
print(latex(B))
if __name__ == "__main__":
linear_hats()
---
title: Problem Set 5 (FEM)
---
## 1
{:.question}
Consider the damped wave equation
$$
\frac{\partial^2 u}{\partial t^2}
= v^2 \frac{\partial^2 u}{\partial x^2}
+ \gamma \frac{\partial}{\partial t} \frac{\partial^2 u}{\partial x^2}
$$
{:.question}
Take the solution domain to be the interval [0, 1].
### (a)
{:.question}
Use the Galerkin method to find an approximating system of differential equations.
We will search for an approximate solution of the form
$$
u(x, t) = \sum_{i=1}^n a_i(t) \varphi_i(x)
$$
(where $$n$$ is yet to be determined). From here on out I will just write $$a_i$$ and $$\varphi_i$$; it
is understood that these are functions of time and space, respectively. Thus the damped wave
equation becomes
$$
\sum_{i=1}^n \left(
\frac{d^2 a_i}{d t^2} \varphi_i
- v^2 a_i \frac{d^2 \varphi_i}{d x^2}
- \gamma \frac{d a_i}{d t} \frac{d^2 \varphi_i}{d x^2}
\right)
= 0
$$
For each basis function $$\varphi_i$$, we construct a Galerkin constraint
$$
\begin{aligned}
0 &=
\int_0^1
\sum_{j=1}^n \left(
\frac{d^2 a_j}{d t^2} \varphi_j
- v^2 a_j \frac{d^2 \varphi_j}{d x^2}
- \gamma \frac{d a_j}{d t} \frac{d^2 \varphi_j}{d x^2}
\right) \varphi_i dx \\
&=
\sum_{j=1}^n
\left(
\frac{d^2 a_j}{d t^2} \int_0^1 \varphi_j \varphi_i dx
- \left(v^2 a_j + \gamma \frac{d a_j}{d t} \right)
\int_0^1 \frac{d^2 \varphi_j}{d x^2} \varphi_i dx
\right) \\
\end{aligned}
$$
This may be written as
$$
\bold{A} \frac{d^2 \overrightharpoon{a}}{d t^2}
- \gamma \bold{B} \frac{d \overrightharpoon{a}}{d t}
- v^2 \bold{B} \overrightharpoon{a}
= 0
$$
where
$$
\begin{aligned}
\bold{A}_{i, j} &= \int_0^1 \varphi_i \varphi_j dx \\
\bold{B}_{i, j} &= \int_0^1 \varphi_i \frac{d^2 \varphi_j}{d x^2} dx \\
&= \left[ \varphi_i \frac{d \varphi_j}{d x} \right]_ {0}^1
- \int_0^1 \frac{d \varphi_i}{dx} \frac{d \varphi_j}{dx} dx \\
\end{aligned}
$$
This gives a system of second order linear homogeneous ODEs. Note that I integrated $$\bold{B}$$ by
parts in order to reduce the order of the highest order derivative of $$\varphi_j$$ that appears.
This lets us use lower order basis functions.
### (b)
{:.question}
Evaluate the matrix coefficients for linear hat basis functions, using elements with a fixed size of
h.
I will assume both boundaries are fixed at zero. Thus it's fine if all hat functions are zero at
each boundary. So I will place $$n$$ hat functions evenly spaced at $$1/h, \ldots, n/h$$, where $$h
= 1/(n + 1)$$. Each hat ranges linearly from 0 to 1, then back down to zero. So for $$1 \leq i \leq
n$$,
$$
\begin{aligned}
\varphi_i(x) &=
\begin{cases}
\Large{\frac{x - (i - 1)h}{h}} &\text{ for } (i - 1)h < x \leq ih \\
\Large{\frac{(i + 1)h - x}{h}} &\text{ for } ih < x < (i + 1)h \\
0 &\text{ otherwise}
\end{cases} \\
\frac{d \varphi_i(x)}{dx} &=
\begin{cases}
\Large{\frac{1}{h}} &\text{ for } (i - 1)h < x \leq ih \\
\Large{\frac{-1}{h}} &\text{ for } ih < x < (i + 1)h \\
0 &\text{ otherwise}
\end{cases}
\end{aligned}
$$
First, let's compute $$\bold{A}$$. Since the hat functions only overlap their immediate neighbors,
the only nonzero elements will be those on the diagonal, and one off it.
First let's compute the diagonal elements. The dependence on $$i$$ drops out after a change of
variables.
$$
\begin{aligned}
\bold{A}_{i, i} &= \int_{(i - 1)h}^{(i + 1)h} \varphi_i(x)^2 dx \\
&= \int_{(i - 1)h}^{ih} \left( \frac{x - (i - 1)h}{h} \right)^2 dx
+ \int_{ih}^{(i + 1)h} \left( \frac{(i + 1)h - x}{h} \right)^2 dx \\
&= \int_{0}^{h} \left( \frac{x}{h} \right)^2 dx
+ \int_{-h}^{0} \left( \frac{-x}{h} \right)^2 dx \\
&= \frac{2}{h^2} \int_{0}^{h} x^2 dx \\
&= \frac{2h}{3}
\end{aligned}
$$
Now we compute the off diagonal elements.
$$
\begin{aligned}
\bold{A}_{i, i+1} &= \int_{ih}^{(i + 1)h} \varphi_i(x) \varphi_{i+1}(x) dx \\
&= \int_{ih}^{(i + 1)h} \left( \frac{(i + 1)h - x}{h} \right) \left( \frac{x - ih}{h} \right) dx \\
&= \int_{0}^{h} \left( \frac{h - x}{h} \right) \left( \frac{x}{h} \right) dx \\
&= \frac{1}{h^2} \int_{0}^{h} \left( hx - x^2 \right) dx \\
&= \frac{1}{h^2} \left( \frac{h^3}{2} - \frac{h^3}{3} \right) dx \\
&= \frac{h}{6}
\end{aligned}
$$
By symmetry, $$\bold{A}_{i, i+1} = \bold{A}_{i, i-1}$$. All other elements are zero.
Now we compute the diagonal elements of $$\bold{B}$$. The boundary term drops out since all
$$\varphi_i$$ are zero at 0 and 1.
$$
\begin{aligned}
\bold{B}_{i, i}
&= \left[ \varphi_i(x) \frac{d \varphi_i(x)}{d x} \right]_ {0}^1
- \int_0^1 \left( \frac{d \varphi_i(x)}{dx} \right)^2 dx \\
&= - \int_{(i-1)h}^{ih} \left( \frac{1}{h} \right)^2 dx
- \int_{ih}^{(i+1)h} \left( \frac{-1}{h} \right)^2 dx \\
&= - \frac{2}{h^2} \int_{0}^{h} dx \\
&= \frac{-2}{h}
\end{aligned}
$$
Finally, the off diagonal elements.
$$
\begin{aligned}
\bold{B}_{i, i + 1}
&= \left[ \varphi_i(x) \frac{d \varphi_i(x)}{d x} \right]_ {0}^1
- \int_0^1 \frac{d \varphi_i(x)}{dx} \frac{d \varphi_{i + 1}(x)}{dx} dx \\
&= - \int_{ih}^{(i+1)h} \left( \frac{-1}{h} \right) \left( \frac{1}{h} \right) dx \\
&= \frac{1}{h^2} \int_{0}^{h} dx \\
&= \frac{1}{h}
\end{aligned}
$$
All other elements are zero.
### (c)
{:.question}
Now find the matrix coefficients for Hermite polynomial interpolation basis functions, once again
using elements with a fixed size of h. A symbolic math environment is useful for this problem.
## 2
{:.question}
Model the bending of a beam (equation 9.29) under an applied load. Use Hermite polynomial
interpolation, and boundary conditions fixing the displacement and slope at one end, and applying a
force at the other end.
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment