You can not select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
2886 lines
98 KiB
2886 lines
98 KiB
|
|
\documentclass[12pt,ragged2e,oneside,a4paper]{book} |
|
|
|
% For embeddeing the source file |
|
\usepackage{attachfile} |
|
\usepackage{caption} |
|
|
|
\renewcommand{\captionfont}{\small} |
|
|
|
|
|
\usepackage[english]{babel} |
|
\usepackage{tocloft} |
|
\usepackage{titlesec} |
|
\usepackage{enumitem} |
|
\usepackage{marginnote} |
|
\usepackage{graphicx} |
|
% for different enumerate styles |
|
\usepackage{multicol} |
|
\usepackage[utf8]{inputenc} |
|
\usepackage[T1]{fontenc} |
|
\usepackage{amsmath} |
|
\usepackage{amsthm} |
|
\usepackage{bookmark} |
|
\usepackage{pgfplots} |
|
\usepackage{import} |
|
\usepackage{xcolor} |
|
\usepackage{soul} |
|
\newcommand{\mathcolorbox}[2]{\colorbox{#1}{$\displaystyle #2$}} |
|
|
|
\usepgfplotslibrary{external} |
|
\tikzexternalize |
|
\pgfplotsset{width=8cm,compat=1.9} |
|
|
|
\usepackage[pdf]{graphviz} |
|
% for unnumbered theorems |
|
|
|
\usepackage{mathtools} |
|
% for := |
|
|
|
\usepackage{amssymb} |
|
% for mathbb |
|
|
|
|
|
|
|
\usepackage[top=2cm, bottom=3cm, outer=7cm, inner=1cm, heightrounded, |
|
marginparwidth=5cm, marginparsep=0.3cm]{geometry} |
|
|
|
\title{MATH3102 Lecture Notes} |
|
\author{Alistair Michael} |
|
|
|
\titleformat{\chapter}{\normalfont\huge\bf}{\thechapter.}{20pt}{\huge\bf} |
|
|
|
\graphicspath{{images/}} |
|
|
|
\newcommand{\sforall}{\enspace\forall} |
|
|
|
\renewcommand{\labelenumi}{\roman{enumi})} |
|
\newcommand{\dsum}{\displaystyle\sum} |
|
\newcommand{\ds}{\displaystyle} |
|
\newcommand{\lsum}{\underline{\int_a^b}} |
|
\newcommand{\usum}{\overline{\int_a^b}} |
|
\newcommand{\rieman}{\int_a^b} |
|
\newcommand{\Arg}{\text{Arg}\,} |
|
\newcommand{\Log}{\text{Log}\,} |
|
\newcommand{\Int}{\text{Int}\,} |
|
\newcommand{\Ext}{\text{Ext}\,} |
|
\newcommand{\p}{\partial} |
|
\newcommand{\pd}{\partial} |
|
|
|
\newtheoremstyle{theoremstyle} % name |
|
{\topsep} % Space above |
|
{\topsep} % Space below |
|
{} % Body font |
|
{} % Indent amount |
|
{\itshape} % Theorem head font |
|
{} % Punctuation after theorem head |
|
{.5em} % Space after theorem head |
|
{} % Theorem head spec (can be left empty, meaning ‘normal’) |
|
|
|
\theoremstyle{theoremstyle} |
|
\newtheorem{thm}{Theorem}[subsection] |
|
\newtheorem{lemma}[thm]{Lemma} |
|
\newtheorem{corr}{Corrolary}[thm] |
|
\newtheorem{claim}{Claim}[thm] |
|
|
|
\newtheorem{definition}{Definition}[subsection] |
|
\newtheorem{prop}{Proposition}[subsection] |
|
\newtheorem{ex}{Example}[subsection] |
|
|
|
\newtheoremstyle{ittheoremstyle} % name |
|
{\topsep} % Spaceabove |
|
{\topsep} % Space below |
|
{\itshape} % Body font |
|
{} % Indent amount |
|
{\scshape} % Theorem head font |
|
{.} % Punctuation after theorem head |
|
{.5em} % Space after theorem head |
|
{} % Theorem head spec (can be left empty, meaning ‘normal’) |
|
|
|
\theoremstyle{ittheoremstyle} |
|
|
|
|
|
\newcommand{\incfig}[1]{% |
|
\def\svgwidth{\columnwidth} |
|
\import{./img/}{#1.pdf_tex} |
|
} |
|
|
|
\newcommand{\incfigw}[2]{% |
|
\begin{center} |
|
\def\svgwidth{#2\columnwidth} |
|
\import{./img/}{#1.pdf_tex} |
|
\end{center} |
|
} |
|
|
|
|
|
\newcommand{\N}{\mathbb N} |
|
\newcommand{\R}{\mathbb R} |
|
\newcommand{\Z}{\mathbb Z} |
|
\newcommand{\C}{\mathbb C} |
|
\newcommand{\Q}{\mathbb Q} |
|
\newcommand{\contra}{\Rightarrow\Leftarrow} |
|
|
|
\setlength{\parskip}{0.5em} |
|
|
|
\setlength{\parindent}{0em} |
|
|
|
\newcommand{\nombreindice}{Contents By Week} |
|
\newlistof{lecture}{lecn}{\nombreindice} |
|
|
|
\newcounter{weekno} |
|
\newcounter{lecno} |
|
|
|
% definition of the commands used for the Spanish ToC; |
|
% \captce for chapters, \sectce for sections and |
|
% \ssectce for subsections |
|
\newcommand\captce[1]{% |
|
\addcontentsline{lecn}{chapter}{\protect\makebox[1.3em][l]{\thechapter}#1}} |
|
\newcommand\week[1]{\stepcounter{weekno}\setcounter{lecno}{0} |
|
\addcontentsline{lecn}{section}{\protect\makebox[4.8em][l]{Week \theweekno}#1}} |
|
\newcommand\lecture[1]{\stepcounter{lecno} |
|
\addcontentsline{lecn}{subsection}{\protect\makebox[6em][l]{Lecture \thelecno}#1} |
|
\marginnote{(Lecture \theweekno.\thelecno)} |
|
} |
|
|
|
\renewcommand{\cfttoctitlefont}{\normalfont\bfseries\Large} |
|
\renewcommand{\cftlecntitlefont}{\normalfont\bfseries\Large} |
|
|
|
\begin{document} |
|
|
|
|
|
|
|
\maketitle |
|
|
|
This document can be found at |
|
\href{https://alistairmichael.com/notes/applied.pdf}{alistairmichael.com/notes/applied.pdf} with source available at |
|
\href{https://git.topost.net/alistair/applied-notes}{git.topost.net/alistair/applied-notes}. |
|
It contains many typos and errors. |
|
|
|
\tableofcontents |
|
\listoflecture |
|
|
|
\pagebreak |
|
|
|
|
|
\week{} |
|
\chapter{Introduction} |
|
\lecture{} |
|
|
|
Applied mathematics is about creating and analyzing mathematical models of |
|
phenomena and processes in the real world. |
|
|
|
The textbook is M. Holmes \emph{Introduction to Foundations of Applied |
|
Mathematics, Springer.} |
|
|
|
|
|
\includegraphics[height=0.5\linewidth]{ab.pdf} |
|
|
|
|
|
\textbf{Modelling:} Creating models from words to equations, identifying key |
|
variables and processes. Models are as simple as possible. |
|
|
|
\emph{Every model is wrong, ... some are useful.} |
|
|
|
Limited insight and predictive power in highly complex systems or systems that |
|
are very sensitive to the initial conditions. |
|
|
|
\section{Topics} |
|
|
|
\begin{itemize} |
|
\item Dimensional analysis |
|
\item Perturbation methods |
|
\item Traffic flow models |
|
\item Continuum mechanics in 1D |
|
\item Elastic and viscoelastic materials |
|
\item Continuum mechanics in 3D |
|
\item Fluid flows |
|
\end{itemize} |
|
|
|
|
|
\chapter{Dimensional Analysis} |
|
|
|
\textbf{Dimensional Estimate:} Qualitative prediction of the order of magnitude |
|
of key quantities without using a model. |
|
|
|
\textbf{Non-dimensionalisation:} Reducing the numbner of effective indpendent |
|
parameters. |
|
|
|
So to create a dimensional estimate for the period of a pendulum $P$, we |
|
understand it may depend on time, distance ($L$), velocity, mass and gravity. |
|
|
|
$P = f(L, v_0, M, g)$ |
|
|
|
How can we combine these paramaters to get the correct dimension as a result? |
|
|
|
$$\frac L g = \frac {m} {m / s^2}$$ |
|
|
|
So as an estimate we can say $P = const \times \sqrt{ \frac L g }$. |
|
|
|
\lecture{} |
|
|
|
\begin{ex} |
|
Ex. Energy of the Bomb |
|
|
|
We want to determine the energy released by a bomb using dimensional analysis. |
|
|
|
We know the radius of the shockwave increases with time, at a rate which depends on the |
|
air density and the energy released. |
|
|
|
\textbf{Radius $R$:} $d$ |
|
|
|
\textbf{Energy:} $m \frac {d^2}{t^2}$ |
|
|
|
\textbf{Air Density ($\rho$):} $\frac {m}{d^3}$ |
|
|
|
So $\frac \rho E = \frac {t^2} {d^5} \implies t^2 \frac E \rho = d^5 $ |
|
|
|
So $R(t) = C \cdot \frac E \rho ^{\frac 15} \cdot t^{\frac 25}$ |
|
\end{ex} |
|
|
|
\subsubsection{Non-dimensionalization} |
|
|
|
\begin{enumerate} |
|
\item[1.] Introduce new units to each independent quantity to reduce the number |
|
of parameters. |
|
\item[2.] Identify relative magnitude of processes and identify dominant |
|
terms and small perturbations. |
|
\item[3.] Choose the non-dimensional form of the equation so that all terms |
|
are either order 1 or small (perturbation). |
|
\end{enumerate} |
|
|
|
\begin{ex} Pendulum |
|
|
|
$\phi$ is angle to overtical, $l$ is length of the pendulum. |
|
|
|
|
|
|
|
\begin{align*} |
|
ml \frac {d^2 \phi}{dt^2} &= -mg \sin(\phi) \\ |
|
\frac {d^2 \phi}{dt^2} &= -\frac gl \sin(\phi) \\ |
|
\end{align*} |
|
|
|
$\frac gl$ can be eliminited by choosing a new time unit. |
|
|
|
Call $t' = \frac t {t_0}$ the nondimensional time unit. |
|
|
|
$g$ = distance / time$^2$. $l$ = distance, $\phi = $ 1. $\sqrt{{\frac lg}} |
|
= $ time $\implies $ choose $t_0 = \sqrt{\frac lg} $ |
|
|
|
\marginnote{$\phi$ is a ratio of two distances it doesnt have a dimension.} |
|
$$\frac {d^2 \phi}{d{t'}^2} = -\sin (\phi)$$ |
|
|
|
Which has no parameters. |
|
\end{ex} |
|
|
|
\begin{ex} Damped Pendulum |
|
|
|
\begin{align*} |
|
ml \frac {d^2\phi}{dt^2} &= - mg\sin(\phi) - kl\frac {d\phi}{dt} \\ |
|
\frac {d^2\phi}{dt^2} &= - \frac gl\sin(\phi) - \frac km\frac {d\phi}{dt} \\ |
|
\end{align*} |
|
|
|
So there are two parameters $\frac km$ and $\frac gl$. |
|
|
|
So choose a new time unit $t = t' t_0$ to make one of the coefficients |
|
disappear. And make $t_0 = \sqrt{\frac lg} $ |
|
|
|
\begin{align*} |
|
\frac {d^2\phi}{dt^2} &= - \frac gl\sin(\phi) {t_0}^2 - \frac km\frac |
|
{d\phi}{dt} t_0 \\ |
|
\frac {d^2\phi}{dt^2} &= - \sin(\phi)- \frac km\frac |
|
{d\phi}{dt} \sqrt{\frac lg} \\ |
|
\frac {d^2\phi}{dt^2} &= - \sin(\phi)- \alpha\frac |
|
{d\phi}{dt} |
|
\end{align*} |
|
|
|
Natural time unit proportional to the period. |
|
|
|
The solution depends on a single parameter $\alpha = (\frac km)\sqrt{\frac l g}$ |
|
\end{ex} |
|
|
|
\begin{ex} |
|
Rabbits and Foxes. |
|
|
|
\begin{align*} |
|
\frac {dR}{dt} &= rR \left(1- \frac RK\right) - pF \frac R{R + S} \\ |
|
\frac {dF}{dt} &= \gamma p F \frac {R}{R + S} - \delta F \\ |
|
\end{align*} |
|
|
|
We want to reduce the parameters: it is a six-dimensional parameter space |
|
for a model of two species. |
|
|
|
\marginnote{So the name of the game is choosing units to absorb all the coefficients.} |
|
Choose natural units for the variables $R$, $F$ and $t$: which look like the |
|
most concrete ones. |
|
|
|
\begin{align*} |
|
R &= R_0 R' \\ |
|
F &= F_0 F' \\ |
|
t &= t_0 t' \\ |
|
\end{align*} |
|
|
|
Exercise. |
|
\end{ex} |
|
|
|
\lecture{} |
|
|
|
\subsubsection{Non-Dimensional Variables and Parameters} |
|
|
|
Where we have a variable $x$ and introduce $x = x_0 x'$ such that $x_0$ |
|
minimises the number of parameters in the model, we call $x'$ a non-dimensional |
|
variable. |
|
|
|
\begin{ex} Cancer Growth (PDE) |
|
|
|
$$\frac {\partial C}{\partial t} = rC \left(1 - \frac CK\right) + D \Delta C$$ |
|
|
|
\begin{itemize} |
|
\item First term is growth from cell division. |
|
\item Second term is diffusion due to random cell movements. |
|
\item $\Delta$ is the laplacian. $\sum \frac {\partial^2}{\partial {x_i}^2}$ |
|
\end{itemize} |
|
|
|
|
|
What are the dimensions of $r$ and $D$? |
|
|
|
\begin{itemize} |
|
\item $r = \frac 1 {time}$ |
|
\item $D = \frac {dist^2}{time}$ |
|
\end{itemize} |
|
|
|
We want to get something with the dimension $\frac {dist}{time}$. |
|
|
|
So $v \sim \sqrt{{rD}}$. |
|
|
|
So use $t_0 = \frac 1r$ $x_0 = \sqrt{{\frac Dr}}$ and $C_0 = K$. |
|
|
|
\begin{align*} |
|
C &= C'K \\ |
|
t &= \frac {t'}r \\ |
|
x &= x'\sqrt{\frac Dr} \\ |
|
\\ |
|
\frac {\partial C'}{\partial t'} &= C' (1 - C') + \Delta' C' |
|
\end{align*} |
|
|
|
\end{ex} |
|
|
|
\chapter{Perturbation Methods} |
|
|
|
Transform the original model into non-dimensional form and identify a small |
|
non-dimensional parameter $\varepsilon < 1$. |
|
|
|
Assume that in the special case $\varepsilon = 0$, that the model simplifies |
|
into a problem that is fully solvable analytically. |
|
|
|
This method gives good approximations for small epsilons. |
|
|
|
\begin{ex} |
|
$m \frac {d^2 y}{dt^2} = -mg $, $y(0) = 0$, $\frac {dy}{dt}|_{t = 0} = v_0$. |
|
|
|
Introduce a nondimensional variables for $t$ and height $y$. |
|
|
|
$y = y_0 t'$ and $t = t_0 t'$ so choose $t_0 = v_0 / g$ and $y_0 = |
|
{v_0}^2 / g$ which leads to the non-dimensional problem |
|
|
|
$$\frac {d^2 y'}{dt'^2} = -1,\quad\quad y'(0)=0,\quad\quad \frac |
|
{dy'}{dt'}|_{t' = 0} = 1$$ |
|
|
|
Its still fully solvable: $y'(t') = - \frac 12 {t'}^2 + t'$. If we want |
|
to generalise it, for example saying gravitational acceleration is no |
|
longer constant. Where $R$ is the radius of the Earth then |
|
|
|
$$\frac {d^2y}{dt^2} = - g \frac {R^2}{(R + y)^2}$$ |
|
|
|
Using the same non-dimensionalisation as before we get: |
|
|
|
$$\frac {d^2 y'}{dt'^2} = -\left(1 + \frac{{v_0}^2}{Rg}y'\right)^{-2},\quad y'(0)=0,\quad \frac |
|
{dy'}{dt'}|_{t' = 0} = 1$$ |
|
|
|
So the small non-dimensional parameter is $\varepsilon = \frac {{v_0}^2}{Rg}$ |
|
|
|
So can find approximate solutions when $\varepsilon \ll 1$ is close to $0$. |
|
|
|
|
|
\end{ex} |
|
|
|
|
|
|
|
\section{Perturbation Methods} |
|
|
|
\week{} |
|
\lecture{} |
|
|
|
\begin{ex} |
|
Solve $x^2 + 2\varepsilon x - 1 = 0$ with perturbation methods. |
|
|
|
The equation has a small non-dimensional parameter $\varepsilon$. |
|
|
|
It can be solved exactly: $x_{1,2} = -\varepsilon \pm |
|
\sqrt{{\varepsilon^2} + 1}$ but we will use it as an example. |
|
See it is easily solvable when $\varepsilon = 0$. |
|
|
|
When $\varepsilon \ll 1$ we can use a taylor series expansion (which is |
|
only valid around $0$ when $\varepsilon \ne 0$). |
|
|
|
{ |
|
$$x_{1,2}(\varepsilon) \approx \pm 1 - \varepsilon \pm \frac 12 \varepsilon^2 \mp \frac 18 \varepsilon^4 + \cdots$$ |
|
\captionof{figure}{Taylor Series Expansion}} |
|
|
|
Look for an approximation of the form: $x(\varepsilon) = c_0 + x_1 |
|
\varepsilon + x_2 \varepsilon^2$ |
|
|
|
Assume $\varepsilon \ll 1$ and $x_0, x_1, x_2,\dots$ are unknown |
|
coefficients of magnitude $\sim 1$. We want to be able to order the |
|
terms by the magnitude according only to the power of $\varepsilon$. |
|
|
|
In $(x_0 + x_1\varepsilon + \dots)^2 + 2\varepsilon(x_0 + x_1\varepsilon |
|
+ \dots) - 1 = 0$ it is clear that the $x_0$ terms dominate, and we can |
|
take out the dominant terms and group them by the powers of epsilon. |
|
|
|
The largest terms where $\varepsilon ^ 0$ |
|
|
|
$${x_0}^2 - 1 = 0 \quad \quad x_0 = \pm 1$$ |
|
|
|
The next largest terms, terms proportional to $\varepsilon^1$: |
|
|
|
$$2x_0 x_1 \varepsilon + 2\varepsilon x_0 = 0 \quad \quad x_1 = -1$$ |
|
|
|
Terms proportional to $\varepsilon^2$: |
|
|
|
$${x_1}^2 \varepsilon^2 + 2 x_0 x_2 \varepsilon^2 + 2x_1 \varepsilon^2 = |
|
0 \quad \quad x_2 = \pm \frac 12$$ |
|
|
|
We can neglect the rest assuming $\varepsilon^3$ is insignificantly |
|
small. |
|
|
|
So we have the approximate solution: |
|
|
|
\begin{align*} |
|
x(\varepsilon) &= x_0 + x_1 \varepsilon + x_2 \varepsilon^2 + \cdots \\ |
|
x(\varepsilon) &= \pm 1 - \varepsilon \pm \frac 12 \varepsilon^2 + \cdots |
|
\end{align*} |
|
\marginnote{Which matches the taylor series.} |
|
|
|
|
|
All we need to be able to solve these kinds of problems is to be able to |
|
solve for the coefficients along the way and use the approximation of |
|
$\varepsilon$ close to $0$: we do not need the explicit form that we |
|
have in this example. |
|
|
|
|
|
\end{ex} |
|
|
|
\textbf{Perturbation method: Solve $F(y; \varepsilon) = 0$ for $\varepsilon \ll 1$ assuming |
|
that $F(y; 0) = 0$ has an exact solution $y_0$ to find an approximate |
|
analytical solution.} |
|
|
|
\begin{definition}{Perturbation Method}{} |
|
|
|
Look for a solution in the form of an asymptotic expansion |
|
|
|
$$y(\varepsilon) = y_0 + y_1 f_1(\varepsilon) + y_2f_2(\varepsilon) + |
|
\cdots $$ |
|
|
|
Where $1 \gg f_1(\varepsilon) \gg f_2(\varepsilon) \gg \cdots$ when |
|
$\varepsilon \ll 1$ (for example a power series like in Example 2.0.1). |
|
|
|
\end{definition} |
|
|
|
\begin{ex} |
|
Find an asymptotic solution for $x^3 + \varepsilon x - 1 = 0$ when |
|
$\varepsilon \ll 1$ |
|
|
|
\marginnote{These problems are artificially simple but it is fine for the |
|
coefficients to contain parameters.} |
|
|
|
Special case: $\varepsilon = 0,\; x^3 - 1 = 0$ has the exact solution $x = |
|
1$. |
|
|
|
$(x_0 + x_1 \varepsilon + x_2\varepsilon^2 \cdots) ^3 + \varepsilon(x_0 |
|
+ x_1\varepsilon + x_2 \varepsilon^2 \cdots) - 1 = 0$ |
|
|
|
\begin{align*} |
|
e^0 &: x_0^3 - 1 = 0 &\implies x_0 = 1 \\ |
|
e^1 &: 3x_0^2x_1 + x_0 = 0, \; 3x_1 + 1 = 0 &\implies x_1 = \frac 13 \\ |
|
e^2 &: 3x_0^2x_2 + 3x_1^2x_0 + x_1 = 0,\; 3x_2 + \frac 13 - \frac 13 = |
|
0 &\implies x_2 = 0\\ |
|
e^3 &: {...} &\implies x_3 = \frac 1 {81} |
|
\end{align*} |
|
|
|
|
|
{ |
|
$$x = 1 - \frac 13 \varepsilon + \frac 1 {81}\varepsilon^3 + \cdots$$ |
|
\captionof{figure}{The Asymptotic Series Solution} |
|
} |
|
|
|
\end{ex} |
|
|
|
\begin{ex} |
|
Point mass in a non-constant gravitational field: |
|
|
|
\begin{equation} |
|
\frac {d^2y'}{dt'^2} = - \left(1 + \frac {{v_0}^2}{Rg}y'\right)^{-2}, y'(0) |
|
= 0, \frac {dy'}{dt'}(t = 0) = 1 |
|
\end{equation} |
|
|
|
Small non-dimensional parameter $\varepsilon = \frac {{v_0}^2}{Rg}$. It |
|
is not solvable for nonzero $\varepsilon$ but can be approximated when |
|
$\varepsilon \ll 1$. |
|
|
|
\marginnote{You could drop the primes after the change of variables to |
|
make the notation less confusing.} |
|
|
|
$$\varepsilon = 0 \implies \frac {d^2y'}{dt'^2} = -1 $$ |
|
|
|
$y'(t') = - \frac 1 2 t'^2 + t'$ |
|
|
|
$$\frac {d^2y'}{dt'^2}\left(1 + \frac {{v_0}^2}{Rg}y'\right)^{2} = -1 $$ |
|
|
|
\marginnote{Rearranging (1)} |
|
|
|
\begin{align*} |
|
\varepsilon^0:& \frac {d^2y_0}{dt'^2} = -1, y_0(0) = 0, \frac |
|
{dy_0}{dt'}(t'=0) = 1 && y_0(t') = - \frac 12 t'^2 + t' \\ |
|
\varepsilon^1:& \frac {d^2y_1}{dt'^2} + 2y_0 \frac {d^2y_0}{dt'^2} = |
|
0 &&\implies \frac{d^2 y_1}{2t'^2} = -t'^2 + 2t' \\ |
|
&\text{And solve with initial condition} |
|
&& \implies y_1(t') = - \frac 13 t'^3 + \frac 1 {12}t'^4 |
|
\end{align*} |
|
|
|
$$y'(t') = - \frac 12 t'^2 + t' - \varepsilon \frac {t'^3(4 - t')}{12} + |
|
\cdots$$ |
|
|
|
\end{ex} |
|
|
|
\subsection{Singular Perturbation} |
|
|
|
\lecture{} |
|
|
|
There are some cases where $\varepsilon = 0$ and small $\varepsilon \ne 0$ |
|
are actually qualitatively different (according to the model). |
|
|
|
A typical case is when setting $\varepsilon = 0$ results in a change of |
|
order of the equation (i.e. when the epsilon is on the highest order term). |
|
These are called \emph{singular perturbation} problems, and the standard |
|
method does not work. |
|
|
|
\begin{ex} $\varepsilon ^2 + x - 1 = 0, \; \varepsilon << 1$ |
|
|
|
If we use the regular perturbation method: |
|
|
|
\[ |
|
x(\varepsilon) = x_0 + x_1 \varepsilon + x_2 \varepsilon^2 + \dots |
|
\quad\quad\qquad (*) |
|
\] |
|
|
|
\marginnote{Substitute $x(\varepsilon)$ into the equation.} |
|
$\varepsilon(x_0 + x_1 \varepsilon + \dots) ^2 + (x_0 + x_1\varepsilon + |
|
\dots) - 1 = 0$ |
|
|
|
|
|
Largest $\varepsilon^0$: $x_0 - 1 = 0 \implies x_0 = 1$. |
|
|
|
\marginnote{Using solution $x_0 = 1$} |
|
$\varepsilon^1: $ $x_0^2 + x_1 = 0 \implies x_1 = -1$ |
|
|
|
So $x = 1 - \varepsilon + 2\varepsilon^2 + \dots$ |
|
|
|
What is the problem with this approach? This is a quadratic equation but |
|
we only got \textbf{one solution.} This is exactly solvable, so we know |
|
the solution should be $x_{1,2} = \frac 1 {2\varepsilon}(-1 \pm |
|
\sqrt{{1-4\varepsilon}})$. |
|
|
|
\begin{center} |
|
\begin{tikzpicture}{testp} |
|
\begin{axis}[xlabel = $\varepsilon$, ylabel=$x$, xtick={0,0.1}] |
|
\addplot[domain=0:0.1, samples=90,color=blue,]{ 1/(2 * x) * (-1 + sqrt(1 -4 * x))}; |
|
\addplot[domain=0:0.1, samples=90,color=blue]{ 1/(2 *x) * (-1 - sqrt(1 - 4* x))}; |
|
\end{axis} |
|
\end{tikzpicture} |
|
\end{center} |
|
|
|
Graphing this solution we can see we have large negative numbers in $x$ as |
|
$\varepsilon$ approaches zero so the assumption about the form of the |
|
solution $(*)$ is not correct: all out coefficients are not order $1$. |
|
|
|
When we set $\varepsilon = 0$ it reduces the order of the equation to |
|
$1$. |
|
\end{ex} |
|
|
|
\subsubsection{Solving Singular Perturbation Problems} |
|
|
|
Use a transformation of variables $x = e^\alpha z$, where |
|
$\alpha$ is determined based on the condition that the highest order |
|
term participates in the dominant part of the equation when $\varepsilon |
|
\to 0$. Then solve it as a regular perturbation problem. |
|
|
|
\begin{ex} Solving $\varepsilon x^2 + x - 1 = 0$. |
|
|
|
\begin{align*} |
|
x = \varepsilon^\alpha z &\implies \varepsilon^{2\alpha + 1}z^2 + |
|
\varepsilon^\alpha z - 1 = 0 |
|
\end{align*} |
|
|
|
We need to match terms of the same order, and there always needs to |
|
be at least two of them so they can cancel out to zero. So we need |
|
to decide which powers of epsilon match up. |
|
|
|
Either: |
|
|
|
\begin{itemize} |
|
\item $2 \alpha + 1 = 0$ |
|
\item $\alpha = 0$, or |
|
\item $2 \alpha + 1 = \alpha$ |
|
\end{itemize} |
|
|
|
|
|
Case 1: |
|
\begin{align*} |
|
2 \alpha + 1 &= 0 \implies \alpha = \frac 12 \\ |
|
\varepsilon^0 z^2 + \varepsilon ^{-\frac 12} - 1 &= 0\\ |
|
\varepsilon^0 z^2 - 1 &= 0 |
|
\end{align*} |
|
|
|
When $\varepsilon \to 0$ the first and last terms are of the same |
|
magnitude but the middle is larger: even though we chose the first |
|
and last to be the dominant terms. This doesnt work so we have to |
|
drop this case. |
|
|
|
Case 2: $\alpha = 0$ |
|
|
|
\begin{align*} |
|
\varepsilon z^2 + \varepsilon^0z - 1 &= 0 \\ |
|
\varepsilon z^2 + z - 1 &= 0 \\ |
|
\end{align*} |
|
|
|
This is not even a transformation and it is stilla singular |
|
perturbation problem. |
|
|
|
Case 3: $2\alpha + 1 = \alpha \implies \alpha = -1$ |
|
|
|
\begin{align*} |
|
\varepsilon^{-1} z^2 + \varepsilon^{-1} - 1 &= 0 \\ |
|
\implies z^2 + z - \varepsilon &= 0 \\ |
|
\end{align*} |
|
|
|
So now the third term $-1$ is small compared to the very large |
|
$\varepsilon^{-1}$ when $\varepsilon \to 0$. If we set $\varepsilon |
|
= 0$ we still have two roots. |
|
\begin{itemize} |
|
\item The highset order term is part of the dominant balance. |
|
\item We don't loose one of the solutions when $\varepsilon = |
|
0$. |
|
\end{itemize} |
|
|
|
|
|
Substitute $z = z_0 + z_1 \varepsilon + \dots$ to get |
|
\begin{align*} |
|
(z_0 + z_1\varepsilon + z_2\varepsilon^2+\dots)^2 + (z_0 + z_1 |
|
\varepsilon + z_2 \varepsilon^2) - \varepsilon &= 0 \\ |
|
z_0^2 + z_0 = 0 \implies z_0 &= 0, -1 \\ |
|
2z_0 z_1 + z_1 - 1 = 0 \implies z_1 &= 1,-1 \\ |
|
\dots |
|
\end{align*} |
|
|
|
Two solutions: $z = \varepsilon - \varepsilon^2 + \dots$ and $z = -1 |
|
- \varepsilon + \varepsilon^2 + \dots$. |
|
|
|
Now do the backwards transformation $x = e^{-1}z$: |
|
\begin{align*} |
|
x &= 1 - \varepsilon + \cdots \\ |
|
x &= -\varepsilon^{-1} - 1 + \varepsilon + \cdots |
|
\end{align*} |
|
\end{ex} |
|
|
|
\hrule |
|
\lecture{} |
|
|
|
Generally considering an $f(x;\varepsilon)$, $x = ?$ when $\varepsilon |
|
\ll 1$ epsilon allowing us to approximate a solution. |
|
\begin{itemize} |
|
\item Substitute a series approximation: $x = x_0 + x_1 \varepsilon |
|
+ x_2 \varepsilon^2 + \cdots$ |
|
\item Use larger numbers of terms and eventually truncate when there |
|
is enough accuracy. |
|
\end{itemize} |
|
|
|
Boundary conditions: |
|
|
|
Eg $\varepsilon y'' 2y' + 2y = 0$ want to find $y$ in $[0,1]$ with |
|
boundary conditions for example $y(0) = 0$ and $y(1) = 1$. |
|
|
|
Note $\varepsilon$ multiplies the term with the highset derivative. If |
|
$\varepsilon = 0$ then it becomes a first order ode and cannot satisfy |
|
both boundary conditions, so it is a singular perturbation problem. |
|
|
|
In this case we can directly solve it for reference: |
|
|
|
\begin{center} |
|
\begin{tikzpicture}{testp} |
|
\begin{axis}[width=0.9\linewidth,height=2in,xlabel = $x$, ylabel=$y$, xtick={0,0.5}] |
|
\addplot[domain=0:0.5, samples=200,color=blue]{-2.71964267074 * |
|
(exp(-18.994 * x) - exp(-1.0553 * x ) }; |
|
\addlegendentry{$\varepsilon= 0.1$} |
|
\addplot[domain=0:0.5, samples=200,color=red]{-2.73205 * |
|
(exp(-198.995 * x) - exp(-1.00505 * x ) }; |
|
\addlegendentry{$\varepsilon= 0.1$} |
|
\addplot[domain=0:0.5, samples=200,color=green]{-2.71964267074 * |
|
(exp(-1999 * x) - exp(-1.0005 * x ) }; |
|
\addlegendentry{$\varepsilon= 0.1$} |
|
\end{axis} |
|
\end{tikzpicture} |
|
\end{center} |
|
|
|
You can see there is a sharp transition point in the solution |
|
called `boundary layers' |
|
|
|
\subsubsection{Boundary Layers, Matched Asymptotic Expansions} |
|
|
|
For a singular perterbation in ODE boundary value problems, for example |
|
$\varepsilon y '' + 2y' + 2y = 0 $ where $y(0) = 0, y(1) = 1$, obviously |
|
this is singular since $\varepsilon = 0$ makes the highest order term |
|
disappear. Standard perterbation method is only good far from the sharp |
|
transition point shown above---the `boundary layers'. |
|
|
|
Matched asymptotic expansions makes an approximation that is valid in |
|
the entire domain. |
|
|
|
\begin{itemize} |
|
\item Use rescaling $x = \varepsilon^{\alpha} z$. We want to create |
|
a transformation that absorbs the singularity so we can solve |
|
the transformed equation using the normal perterbation method |
|
and then transform it back. |
|
|
|
|
|
So if we were to set $\varepsilon = 0$ the order should not |
|
change. |
|
\item Chain rule transforms derivatives: |
|
|
|
$Y'(z) = \frac{dY}{dx} \frac {dx}{dz} = y'(x) \varepsilon^a$ |
|
|
|
$Y''(z) = y''(x) \varepsilon^{2a} $ |
|
\item After change of variables to the rescaled equation do a |
|
regular perturbation so we still get a nonzero second derivative |
|
term. |
|
|
|
$\varepsilon y''(x) 2y'(x) + 2y(x) = 0 \to$ |
|
|
|
Rescaled equation: $\varepsilon^{1-2\alpha}Y''(z) + \varepsilon^{-\alpha}2Y'(z) + |
|
2Y(z) = 0$ |
|
|
|
\week{}\lecture{} |
|
\item Determine $\alpha$ |
|
|
|
Choose an alpha so the term with the highest derivative is part |
|
of the dominant balance. |
|
\begin{enumerate} |
|
\item $1 - 2 \alpha = -\alpha \implies \alpha = 1$ |
|
\item $1 - 2\alpha = 0 \implies \alpha = \frac 12$ |
|
\item $-\alpha = 0 \implies \alpha = 0$ (no rescaling) |
|
\end{enumerate} |
|
|
|
\begin{enumerate} |
|
\item $1 - 2 \alpha = -\alpha \implies \alpha = 1$ |
|
\begin{align*} |
|
\varepsilon ^{-1} Y''(z) + \varepsilon^{-1}2Y'(z) + |
|
2Y(z) &= 0 \implies \\ |
|
Y''(z) + 2Y'(zz) + \varepsilon2Y(z) &= 0 |
|
\end{align*} |
|
The two largest terms (with |
|
$\varepsilon^{-1}$) can cancel out and the other is a |
|
small correction |
|
|
|
So this one is okay. |
|
|
|
\item $1 - 2\alpha = 0 \implies \alpha = \frac 12$ |
|
|
|
$Y''(z ) + \varepsilon^{-\frac 12} 2Y'(z) + 2Y(z)$ |
|
|
|
The largest term (the middle term) does not have another |
|
term of the same order which can balance it out to reach |
|
a possible solution, so this cannot produce a solution. |
|
\end{enumerate} |
|
|
|
\item So we choose $\alpha = 1$, and can now sovle with |
|
perterbation. |
|
|
|
\emph{Boundary condition:} We only expect this to give a good |
|
solution near $x=$, not for the whole domain, so we can only use the boundary |
|
condition at $x = 0$. |
|
|
|
$\alpha = 1 \implies Y''(z) + 2Y'(z) + \varepsilon2Y(z) = |
|
0,\quad Y(0) = 0$ |
|
|
|
\begin{align*} |
|
Y(z) &= Y_0(z) + \varepsilon Y_1(z) + \cdots \\ |
|
0 &= (Y_0(z) + \varepsilon Y_1(z) + \dots)'' + 2(Y_0(z) + |
|
\varepsilon Y_1(z) + \dots)' + \varepsilon(Y_0(z) + |
|
\varepsilon Y_1(z) + \dots) \\ |
|
\end{align*} |
|
\begin{itemize} |
|
\item [$e^0$] We get $Y_0'' + 2Y_0'(z) = 0$, which is |
|
solvable: |
|
|
|
$Y_0(z) = e^{\lambda z} \implies \lambda^2 + 2\lambda = |
|
0 \implies \lambda_1 = 0, \lambda_2 = -2$ so |
|
|
|
$y_0(z) = c(1 - e^{-2z})$, with an unknown constant $c$. |
|
|
|
Transforming back: $x = \varepsilon z \implies$ |
|
|
|
\marginnote{The boundary layer solution, which is only |
|
valid around $x = 0$.} |
|
$$y_0(x) = c(1 - e^{\frac{-2x} \varepsilon})$$ |
|
\end{itemize} |
|
|
|
\item Matching the two solutions: |
|
|
|
Boundary layer (sometimes called inner solution), valid near $x = 0$. $y_{BL}(x) = c(1 - e^{\frac{-2x} |
|
\varepsilon})$. |
|
|
|
\marginnote{$y_{out}$ obtained before using normal perterbation} |
|
Outer solution: $y_{out} = e^{1 - x}$. Valid away from $x = 0$. |
|
|
|
We want to combine them to crate an approximation that is valid |
|
for the whole domain. |
|
|
|
Find the \textbf{common limit}: |
|
|
|
$y_{BL} (x, \varepsilon \to 0) = c$ |
|
|
|
$Y_{out}(x \to 0) = e$. |
|
|
|
This implies $e^1 = c$, $y_{BL} = e(1 - e^{\frac{-2x}{\varepsilon}})$ |
|
|
|
\item Combine to get composite solution: |
|
\begin{align*} |
|
y(x) &= y(x)_{OUT} + y(x)_{BL} - \text{common limit} \\ |
|
&= y(x) e^{1-x} + e(1 - e^{-2x/\varepsilon}) - e \\ |
|
&= e^{1-x} - e^{1 - \frac {2x}e} |
|
\end{align*} |
|
|
|
|
|
\end{itemize} |
|
|
|
By solving directly we know the exact solution for $\varepsilon = 0.1$ |
|
is |
|
|
|
$-2.874 \left(e^{-18.944 x}-e^{-1.056 x}\right)$ |
|
|
|
\begin{center} |
|
\begin{tikzpicture}{testp} |
|
\begin{axis}[width=0.9\linewidth,height=2in,xlabel = $x$, ylabel=$y$, xtick={0,0.5}] |
|
\addplot[domain=0:0.5, samples=200,color=red]{exp(1-x) }; |
|
\addlegendentry{$y_{out}(x)$} |
|
\addplot[domain=0:0.5, samples=200,color=blue]{ e * (1 - |
|
exp((-2 * x) / |
|
0.1))}; |
|
\addlegendentry{$y_{BL}(x)$} |
|
\addplot[domain=0:0.5, samples=200,color=green]{exp(1-x) - exp(1 - |
|
2 * x / 0.1)}; |
|
\addlegendentry{Combined} |
|
\addplot[domain=0:0.5, samples=200,color=orange]{-2.71964267074 * |
|
(exp(-18.994 * x) - exp(-1.0553 * x ) }; |
|
\addlegendentry{Exact} |
|
\end{axis} |
|
\end{tikzpicture} |
|
\end{center} |
|
|
|
|
|
\lecture{} |
|
|
|
Mathematica. |
|
|
|
\lecture{} |
|
|
|
\subsubsection{Methods of Multiple Scales} |
|
|
|
\begin{ex} |
|
Nonlinear pendulum without damping: |
|
|
|
$ml\theta'' = -mg \sin(\theta)$ |
|
|
|
In nondimensiona form: |
|
$\theta '' = \sin(\theta),\quad\theta(t) = ?$ |
|
|
|
To do perterbation we add an assumption, $\theta(0) = \varepsilon$, |
|
$\theta'(0) = 0$, and a small initial angle $\epsilon \ll 1$ |
|
|
|
Trying regular perterbation: $\theta(t) = \theta_0 (t) + \varepsilon |
|
\theta_1(t) + \cdots. $ |
|
Use taylor expansion of $\sin(\theta)$ around $0$ for small |
|
oscillations. $\sin(\theta) \approx \theta - \frac {\theta^3} 6 + |
|
\cdots $ So when oscillation is small you can use the taylor series to |
|
represent $\sin(\theta)$. |
|
|
|
Sub the expansion in the initial conditions: $\theta_0(0) + |
|
\varepsilon \theta_1(0) + ... = \varepsilon. $ so $\theta_0(0) = 0$, |
|
$\theta_1(0) = 1$. And |
|
|
|
$\theta_0'(0) + \varepsilon\theta_1'(0) + ... = 0 \implies |
|
\theta'_0(0) = 0, \theta_1'(0) = 0$. |
|
|
|
\begin{align*} |
|
\varepsilon^0 :&& \theta_0'' = -\theta_0' + \frac 16 |
|
(\theta_0^3) + \cdots = -\sin(\theta_0) \\ |
|
&&\theta_0(0) = 0\quad \theta_0'(0) = 0 \\ |
|
&&\theta_0(t) = 0 \\ |
|
\varepsilon^1:&& \theta_1'' = -\theta_1 \\ |
|
&& \theta_1(0) = 1, \theta_1'(0) = 0 \\ |
|
\theta_1(t) = \cos(t) \implies \theta(t) = |
|
\varepsilon \cos(t) |
|
\end{align*} |
|
|
|
The approximation is not uniform in time. |
|
|
|
Solving for the next power of $\varepsilon$... |
|
|
|
$$\theta(t) = \varepsilon\cos(t) - \frac 1 {16} \varepsilon^3 t \sin(t) + |
|
\cdots$$ |
|
|
|
For very large $t$ the second term may be larger than the first |
|
which is a problem since we expect the order to be based only on the |
|
power of $\varepsilon$. |
|
|
|
\begin{center} |
|
\begin{tikzpicture}{testp} |
|
\begin{axis}[width=0.9\linewidth,height=2in,xlabel = $t$, |
|
xtick={180,190,200}] |
|
\addplot[domain=0:50, samples=200,color=blue]{ exp(-0.1/(2*x)) * |
|
sin(1/2 * sqrt(4 - 0.1^2) * deg(x)) }; |
|
\addlegendentry{Exact} |
|
\addplot[domain=0:50, samples=200,color=red]{ (1 - 0.1*x/2) * |
|
sin(deg(x)) }; |
|
\addlegendentry{Approx} |
|
\end{axis} |
|
\end{tikzpicture} |
|
\end{center} |
|
|
|
There are two timescales in the problem the oscillation period and the |
|
slow drift in the phase of the oscillation which is not described by |
|
the first term so we can add parameters to the time to split it into |
|
multiple timescales. |
|
|
|
\end{ex} |
|
|
|
\week{} |
|
\lecture{} |
|
|
|
The ordering of the terms, their relative magnitude, may change in a |
|
way which violates our assumptions. This is what the method of multiple |
|
scales addresses. |
|
|
|
There are two characteristic timescales: oscillation period, |
|
the time of decay. |
|
|
|
\subsubsection{Method of multiple scales} |
|
|
|
\begin{itemize} |
|
\marginnote{(the equation is already non-dimensional)} |
|
\item Introduce two new nondimensional time variables to rescale. |
|
|
|
$t_1 = t$ and $t_2 = \varepsilon^\alpha t$ so $y(t) \implies Y(t_1,t_2)$ |
|
|
|
\item Use the chain rule: |
|
|
|
$$\frac d{dt} = \frac {\partial}{\partial t_1} |
|
+ \varepsilon^\alpha \frac \partial{\partial t_2}$$ |
|
|
|
$$\frac{d^2}{dt^2} = \frac{\partial^2}{\partial {t_1}^2} + |
|
2\varepsilon^\alpha \frac {\partial^2}{\partial t_1\partial t_2} + |
|
\varepsilon^{2\alpha} \frac {\partial^2}{\partial {t_2}^2}$$ |
|
|
|
\item Initial conditions: |
|
|
|
$Y(0,0) = 0, \quad \frac{\partial Y}{\partial t_1}(0,0) + |
|
\varepsilon^\alpha \frac {\partial Y}{\partial t_2}(0,0) = 1$ |
|
\item Substitute |
|
|
|
$$Y(t_1,t_2) = Y_0(t_1,t_2) + \varepsilon Y_1(t_1,t_2) + |
|
O(\varepsilon^2)$$ |
|
\item |
|
|
|
Largest:$\varepsilon^0: \frac{\partial^2 Y_0}{\partial t_1^2} + Y_0 = |
|
0$ |
|
|
|
$Y_0(0,0) = 0, \quad \frac{\partial Y_0} {\partial t_1}(0,0) = 1$ |
|
|
|
We know there is a general solution: |
|
|
|
$$Y_0(t_1, t_2) = A(t_2) \sin(t_1) + B(t_2)\cos(t_2)$$ |
|
|
|
Where $B(0) = 0$ and $A(0) = 1$ is the initial condition, by |
|
susbstituting the initial conditions in the general solution. |
|
|
|
We will choose $A(t)$ and $B(t)$ in such a way as to avoid |
|
resonancnes. |
|
|
|
|
|
In order for the friction terms to be balanced by the mixed |
|
derivatives choose $\alpha = 1$, (slow time: $t_2 = \varepsilon^1 |
|
t$) |
|
|
|
|
|
Terms $\varepsilon^1$: |
|
|
|
$\frac {\partial^2 Y_1} {\partial t_1^2} + Y_1 = - \frac {\partial |
|
y_0} {\partial t_1} - 2 \frac {\partial ^2 Y_0} {\partial t_1 |
|
\partial t_2}$ |
|
|
|
Substitute |
|
|
|
$\frac {\partial^2 Y_1} {\partial t_1^2} + Y_1 = - A(t_2)\cos(t_1) + |
|
B(t_2) \sin(t_1) - 2A'\cos(t_1) + 2B'(t_2) \sin(t_1)$ |
|
|
|
Group together forcing terms. |
|
|
|
$\frac {\partial^2 Y_1} {\partial t_1^2} + Y_1 = [- A(t_2) + 2A'(t_2)]\cos(t_1) + |
|
[B(t_2) \sin(t_1) + 2B'(t_2)] \sin(t_1)$ |
|
|
|
To avoid resonant forcing which leads to secular terms we need $A + 2A' |
|
= 0$ and $B + 2B' = 0$. So we need to solve these differential |
|
equations: |
|
|
|
$$A(t_2) = A(0) e^{\frac{-t_2}2} \quad \quad B(t_2) = B(0) |
|
e^{\frac{-t_2}2}$$ |
|
|
|
\item Transform back and substitute |
|
|
|
$y(t) = e^{\frac {-\varepsilon t}2 } \sin(t)$ |
|
|
|
\end{itemize} |
|
|
|
|
|
\chapter{Traffic Flow Models} |
|
\lecture{} |
|
|
|
\textbf{Agent based models:} rules for individual cars and how they |
|
interact can be stochastic. can be used for simulations, not suitable |
|
for mathematical analysis. |
|
|
|
\textbf{Conntinuum model:} car density function and velocity |
|
distribution etc. Treated as a continuous variable. |
|
|
|
\begin{definition}{}{} Density function (1D) |
|
|
|
$\rho(x,t) = \frac {\text{The number of cars in } [x - \Delta x, x + \Delta x]} { |
|
\text{distance } 2\Delta x}$ |
|
|
|
The $\Delta x$ needs to be large compared to the size of the cars but |
|
small compared to the region. |
|
\end{definition} |
|
|
|
\begin{definition}{}{} Flux $J(x,t)$ |
|
|
|
The number of cars that passing a position $x$ per unit time |
|
$\frac {[t-\Delta t, t + \Delta t]} {2\Delta t}$. |
|
|
|
For constant velocity of cars with distantce $d$ $J = \frac v {l + d} = |
|
v \rho$. |
|
\end{definition} |
|
|
|
\begin{definition}{}{} Balance Law |
|
|
|
Conservation law for cars along the region $[x - \Delta x, x + |
|
\Delta x]$ of the road. |
|
|
|
Initial cars + cars going in - cars going out. |
|
|
|
$$J(x - \Delta x, t) \Delta t - J(x + \Delta x , t) \Delta t = |
|
\frac{J(x - \Delta x, t) - J( x + \Delta x, t)}{2\Delta x}$$ |
|
|
|
Taking the limit $\Delta t \to 0$ and $\Delta x \to 0$ |
|
|
|
$$\frac {\partial \rho}{\partial t} = - \frac {\partial J}{\partial |
|
x}$$ |
|
|
|
The flux can also be written as $J = \rho v$ so we have the equation |
|
|
|
$$\frac {\partial \rho}{\partial t} + \frac {\partial (\rho |
|
v)}{\partial x} = 0$$ |
|
\end{definition} |
|
|
|
In order to obtain a closed system we need a relationship between the |
|
velocity and density of traffic. |
|
|
|
\begin{definition}{}{} |
|
Constitutive law |
|
|
|
$v(\rho)$ is known as the constitutive law |
|
|
|
It is natural to expect $v(\rho)$ is monotonously decreasing: as the |
|
density of cars increases the traffic slows down. And $v(0) = v_M$ when |
|
there are no other cars around traffic moves at the maximum allowable |
|
velocity. |
|
|
|
This can be based on observations from real traffic. |
|
\end{definition} |
|
|
|
The simplest constitutive law is greenshield's law assumes a linear |
|
function: |
|
|
|
$$v(\rho) = v_M\left(1 - \frac \rho{\rho_M}\right)$$ |
|
|
|
substituting into the balance law this gives |
|
|
|
$$\frac {\partial \rho}{\partial t} + c(\rho) \frac {\partial (\rho |
|
)}{\partial x} = 0 \quad\quad \text{where} $$ |
|
|
|
|
|
$$c(\rho) = v_M\left(1 - \frac {2\rho}{\rho_M}\right)$$ |
|
|
|
If you assume an almost uniform density and they propogate like waves |
|
you can get wave velocity, which is the velocity of propogation of |
|
perterbations in the density of traffic: |
|
|
|
$$c(\rho) = v(\rho) + \rho v'(\rho)'$$ |
|
|
|
$$\frac {\partial \rho} { \partial t} + c(\rho) \frac {\partial \rho}{\partial x} = |
|
0$$ |
|
|
|
\subsubsection{Solving} |
|
|
|
Consider the simplest case, which is consntant velocity: $v(\rho) = a$. |
|
But we don't assume density is uniform: $v$ is independent of the density, |
|
which is not realistic for traffic it |
|
is just a simpler maths problem; advection equation for fluid flow. |
|
|
|
The solution is $\rho(x,t) = f(x - at)$, can check this with |
|
substitution. |
|
|
|
$$\frac {\partial \rho}{\partial t} + a \frac {\partial \rho }{\partial x} |
|
= 0$$ |
|
|
|
Initial: $\rho (x,0) = f(x) $ |
|
|
|
Boundary: $\rho (x = 0, t) = g(t)$. |
|
|
|
Don't need boundary conditions on both sides of the domain, because it |
|
only goes in one direction. |
|
|
|
Solution : $$\rho(x,t) = \begin{cases}f(x - at) & x > at \\ g(t - \frac |
|
)xa & x < at\end{cases}$$ |
|
|
|
Nonconstant velocity: |
|
|
|
$$\frac {\partial \rho}{\partial t} + c (\rho) \frac {\partial |
|
\rho}{\partial x} = 0, \quad \rho(x,t=0) = f(x)$$ |
|
|
|
The density changes according to the velocity. |
|
|
|
\week{} |
|
|
|
\lecture{} |
|
|
|
\begin{ex} |
|
\end{ex} |
|
|
|
\begin{itemize} |
|
\item uniform positive velocity |
|
\item density linearly linked to velocity according to greenfields law |
|
\end{itemize} |
|
|
|
|
|
|
|
Perterbations propagate along the roads according to wave equation: |
|
areformation of hte balance equation: |
|
|
|
$$\frac {\partial \rho}{\partial t} + c (\rho) \frac {\partial \rho}{\partial x} |
|
= 0$$ |
|
|
|
Where $c(\rho) = J' (\rho) = v + \rho v'(\rho)$. |
|
|
|
\textbf{Special case with constatnt velocity} gives the travelling wave |
|
solution: $\rho(x,t) = f(x - ct)$ the density distribution is transported along hte x axis with |
|
constatnt velocity without any distortion. The shape is exactly the |
|
same. This is when the cars dont care about the density. |
|
|
|
\textbf{Method of characteristic trajectories}: $x(t) = x_0 + ct$. |
|
find th trajectory ending at $(x,t) \implies x_0 = x(t) - ct$; |
|
The density is conserved along the |
|
lines $x(t) = x_0 + c(\rho(x_0))t$. (A transformation into an arbitrary |
|
space such that there are lies where $\rho$ is constant to make it easy to |
|
solve.) |
|
|
|
|
|
The same method works for non-constant velocity, we know the waves propogate at |
|
constant speed (because of our assumption) but we don't know the slope. |
|
|
|
For any $(x,t)$ find an $x_0$ such that $x = x_0 + c(f(x_0))t$ then $\rho(x,t) |
|
= \rho(x_0, 0) = f(x_0)$. |
|
|
|
\emph{characteristics:} points $x_0$ that move to the point $(x,t)$. |
|
|
|
\begin{itemize} |
|
\item There may be regions with no characteristics |
|
\item there may be regions with multipl eintersecing characteristics |
|
\end{itemize} |
|
|
|
\section{Rankine-Hugoniot condition:} |
|
|
|
\begin{itemize} |
|
\item If there are discontinities in the initial conditions of density |
|
($\rho$), you can |
|
end up with diverging or overlapping characteristic lines and there can |
|
be no $\rho_l < \rho_r$ or no $\rho_l > \rho_r$ unique solutions for a particular region. |
|
\end{itemize} |
|
|
|
\begin{center} |
|
\begin{tikzpicture}{testp} |
|
\begin{axis}[xmin=0, xmax=10, ymin=0, ymax=10, ticks=none, xlabel=x, ylabel=t] |
|
\addplot[domain=0:7.5,samples=90,color=red]{ 2.5 * (x - 5)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 2 * x}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 2 * (x - 1)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 2 * (x - 2)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 2 * (x - 3)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 2 * (x - 4)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 2 * (x - 5)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 3 * (x - 5)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 3 * (x - 6)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 3 * (x - 7)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 3 * (x - 8)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 3 * (x - 9)}; |
|
\end{axis} |
|
\end{tikzpicture} |
|
\end{center} |
|
|
|
|
|
You can integrate the balance equation to get information about the |
|
discontinuity. Integrate $x$ in $[s(t) - \varepsilon, s(t) + \varepsilon]$, |
|
givin |
|
|
|
$$\int_{s - \varepsilon}^{s + \varepsilon} \frac {\partial \rho}{\partial t} dx |
|
+ J(\rho_r) - J(\rho_l) = 0$$ |
|
|
|
Can differentiate this again using \textbf{Leibniz rule}: the moving integrating bounds |
|
changes the value of the integral. |
|
|
|
If they are constant: |
|
|
|
$$\frac d {dt} \int_\alpha ^\beta f(x,t) dx = \int_\alpha^\beta \frac {\partial |
|
f} {\partial t} dx$$ |
|
|
|
Otherwise |
|
|
|
$$\cdots$$ |
|
|
|
Specifically for the non-constant functions $\alpha(t) = s - \varepsilon$ and $\beta(t) = s + |
|
\varepsilon$ |
|
|
|
$$\frac d {dt} \int_{x - \varepsilon}^{s + \varepsilon} \rho dx - s'(t) \rho_r + |
|
s'(t) + \rho_l + J(\rho_r) - J(\rho_l) = 0$$ |
|
|
|
|
|
$$\frac d {dt} \int_{x - \varepsilon}^{s + \varepsilon} \rho dx = \frac d {dt} |
|
(\epsilon \rho_l + \epsilon \rho _r) = 0$$ |
|
|
|
\marginnote{cool} |
|
$$s'(t) = \frac {J(\rho_r) - J(\rho l)} {\rho_r - \rho_L} \quad\quad\quad(*)$$ |
|
|
|
From this we can also show that the discontinuity moves with the velocity |
|
|
|
$$s'(t) = \frac 1 {\rho_r - \rho _l } \int_{\rho_l} ^{\rho_r} c(\rho) d\rho$$ |
|
|
|
|
|
\begin{itemize} |
|
\item This solution is not limited to greenshield's law it gives |
|
caharacteristic lines for intersections of the characteristic lines and |
|
they are related to the function $c(\rho)$ and the way these boundaries |
|
move through traffic (because of hte moving $s$ in the integral before). |
|
\end{itemize} |
|
|
|
Under \textbf{Greenshield's law} it becomes the average speed of the wave velocity: |
|
|
|
$$s'(t) = \frac 12 (c_R + c_L)$$ |
|
|
|
|
|
We want to solve these kinds of initial conditions cause the models often |
|
produce steps (singularities in the solution) or 'shocks' between fast and slow |
|
moving regions (because the faster cars `run into' the slow cars and teh |
|
slowdown region gets smaller until it is a jump). |
|
|
|
note that $c(\rho)$ is not strictly positive: perturbation waves can move |
|
backwards through the traffic, eg if a car suddenly stops on a highway. |
|
|
|
|
|
\lecture{} |
|
|
|
\textbf{Initial Conditions $\rho_l > \rho_r$} |
|
|
|
\begin{center} |
|
\begin{tikzpicture}{testp} |
|
\begin{axis}[xmin=0, xmax=10, ymin=0, ymax=10, ticks=none, xlabel=x, ylabel=t] |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 3 * x}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 3 * (x - 1)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 3 * (x - 2)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 3 * (x - 3)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 3 * (x - 4)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 3 * (x - 5)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 2 * (x - 5)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 2 * (x - 6)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 2 * (x - 7)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 2 * (x - 8)}; |
|
\addplot[domain=0:20,samples=90,color=black,very thick,dotted]{ 2 * (x - 9)}; |
|
\end{axis} |
|
\end{tikzpicture} |
|
\end{center} |
|
|
|
The initial jump is dissolved into a transition zone that gets wider |
|
over time: ``expansion fan''. In this case the Rankine-Huginot |
|
condition cant be used, since the jump discontinuity doesnt persist for |
|
greater $t$. |
|
|
|
\begin{ex} Left side maximum density and zero velocity, right side |
|
maximum velocity and zero density. |
|
|
|
\begin{itemize} |
|
\item Smooth out the initial step function with a linear interpolation |
|
|
|
$$f(x) = \begin{cases} \rho_M & x < -\varepsilon \\ |
|
\rho_M \frac{\varepsilon - x}{2\varepsilon} & -\varepsilon < x < |
|
\varepsilon \\ |
|
0 & \varepsilon < x |
|
\end{cases}$$ |
|
|
|
In the transition zone you have $-\varepsilon < x < \varepsilon$ |
|
and greenshields law to calculate wave velocity $c(f(x)) = V_M \frac x |
|
\varepsilon$. |
|
|
|
solve stuff |
|
\end{itemize} |
|
|
|
\end{ex} |
|
|
|
\week{} |
|
\lecture{} |
|
\chapter{Continuum Mechanics} |
|
|
|
Describing the movement and deformation dymaics of materials represented |
|
as a continuous medium (rather than representations as connected |
|
particles). For example elastic solids, flowing gases and liquids. |
|
|
|
|
|
\incfig{ref} |
|
|
|
|
|
\subsection{Coordinates} |
|
|
|
Consider two different coordinate systems: |
|
\begin{itemize} \item [] \textbf{Material coordinates:} label a pooint within the |
|
continous material and follow its movement along hte |
|
trajectories $X(A,t)$ where $A$ is the intial location: $X(A,0) |
|
= A$. (If you painted a dot on a piece of jelly). |
|
\item [] \textbf{Spatial coordinates:} Fixed coordinate system $x$ |
|
of an external observer. Does not follow the movement of the |
|
material. |
|
\end{itemize} |
|
|
|
\textbf{Material Coordinates:} |
|
|
|
\begin{itemize} |
|
\item Trajectory $X(A,t)$ where $X(A,0) = A$ |
|
|
|
\item Displacement: $U(A,t) = X(A,t) - A$ |
|
|
|
\item Velocity: $V(A,t) = \frac {\partial U(A,t)} {\partial t}$ |
|
\end{itemize} |
|
|
|
\textbf{Spatial Coordinates:} |
|
|
|
\begin{itemize} |
|
\item Displacement: $u(x,t)$ of the material $x$ at time $t$. |
|
To relate it to $U(A,t) $ find $A$ at $t = 0$ by solving $X(a,t) = x$ for |
|
$A = a(x,t)$ (a in the spatial coordinates). |
|
|
|
\item $u (x,t) = U(a(x,t) , t)$ |
|
\item $v(x,t) = \frac {\partial U(a(x,t),t)} {\partial t} = V(a(x,t),t)$ which is not equal to $\frac {\partial |
|
u}{\partial t}$ (there is an additional dependence on time, we |
|
evaluate by traversing the trajectory function $A$). |
|
\end{itemize} |
|
|
|
Example: |
|
|
|
Consider a material described by the trajectories $X(A,t) = A \frac{1 + |
|
2t}{1 + t}$. Initial domain: $0 < A < l_0$. |
|
|
|
The displacement from $A$ of a material point at time $t$ is |
|
|
|
$$U(A,t) = X (A,t) - A = A \frac t {1 + t}$$ |
|
|
|
\marginnote{Can check this is different from $\frac{\partial u}{\partial |
|
t}$. However it is the same as $\frac {\partial X}{\partial t}$ since $A$ |
|
is a constant.} |
|
|
|
The velocity: |
|
$$V(A,t) = \frac {\partial U(A,t)}{\partial t} = A \frac 1{(1 + t)^2}$$ |
|
|
|
The domain in material coordinates is $0 < A < l_0$, in spatial |
|
coordinates $X(0,t) = 0 < x < X(l_0,t) = l_0 \frac {1 + 2t}{1 + t}$ |
|
|
|
If final coordinate is $x$ at $t$ find the initial coordinate $A$ |
|
|
|
$$x = X(A,t) \implies A = x \frac {1 + t}{1 + 2t}$$ |
|
|
|
The displacement of a material point at $x$ at $t$ relative to its |
|
initial position: |
|
|
|
|
|
$$u(x,t) = U(x\frac{1 + 2t}{1 + t}, t) = x \frac t {1 + 2t}$$ |
|
|
|
$$v(x,t) = V(x \frac {1 + 2t}{1 + t}, t) = x \frac 1 {(1 + t)(1 + 2t)}$$ |
|
|
|
\subsection{Material Derivative} |
|
|
|
Consider a scalar quantity that describes a propertiy transported by the |
|
material. In spatial coordinates the distribution is described by |
|
$f(x,t) $ and in material coordinates is $F(A,t)$. They must satisfy the |
|
relationship $f(X(A,t),t) = F(A,t)$. |
|
|
|
\emph{What is the rate of change of this property in a material |
|
coordinate system (along the trajectories of movement of the material |
|
points).} |
|
|
|
\begin{align*} |
|
\frac {\partial F}{\partial t} &= \frac {\partial f}{\partial t} + |
|
\frac {\partial f} {\partial x} \frac {\partial X(A,t)} {\partial t} \\ |
|
&= |
|
\frac {\partial f}{\partial t} + V(A,t) \frac {\partial f}{\partial x} |
|
\\ &= |
|
\frac {\partial f}{\partial t} + v(x,t) \frac {\partial f} {\partial x} |
|
\end{align*} |
|
|
|
Material Derivative: $$\frac D {Dt} \equiv \frac \partial{\partial t} + |
|
v \frac \partial {\p x} \equiv \frac {\partial F}{\partial t}$$ |
|
|
|
\lecture{} |
|
|
|
\subsection{Continuity Equation} |
|
|
|
The law of conservation of mass |
|
|
|
$$\frac {\partial \rho}{\partial t} + \frac {\partial(v \rho)}{\partial |
|
x} = 0$$ |
|
|
|
This is not sufficient to solve, we need another equation to describe |
|
the velocity field $v(x,t)$. We use a continuum mechanics version of |
|
newtons law: conservation of momentum. |
|
|
|
\subsection{Momentum Equation} |
|
|
|
The motion of a point particle under forces: $ma = \sum F_i$. |
|
The rate of change of momentum: |
|
|
|
$$\frac d {dt} mv = F$$ |
|
|
|
Typoes of forces in continuum mechanics |
|
|
|
\begin{itemize} |
|
\item external surface forces |
|
\begin{itemize} |
|
\item Act along the boundary of the domain |
|
\end{itemize} |
|
\item external body forces |
|
\begin{itemize} |
|
\item Act on every material point inside the domain, for |
|
example gravity. |
|
\end{itemize} |
|
\item internal forces |
|
\begin{itemize} |
|
\item Stress ($\tau(x,t) = $ force / area) generated by different parts of |
|
the material acting on eachother. For example restoring |
|
force in elastic material. |
|
\end{itemize} |
|
\end{itemize} |
|
|
|
|
|
Consider a segment of material between $A_L$ and $A_R$ with trajectories |
|
$X(A_L,t) = \alpha (t) $ and $X(A_R,t) = \beta(t)$. |
|
|
|
The rate of change of the momentum of this segment is the sum of the |
|
forces acting on it, external body force and internal surface force. |
|
\marginnote{$S$ is the surface force. } |
|
|
|
$$\frac d {dt} \int_{\alpha(t)}^{\beta (t)} S \rho v \; dx = |
|
\int_{\alpha(t)}^{\beta(t)} S\rho f \; dx + S\tau(\beta,t) - S |
|
\tau(\alpha, t)$$ |
|
|
|
$$\frac d {dt}\int_{\alpha(t)}^{\beta(t)} \rho v \; dx = |
|
\int_{\alpha(t)}^{\beta(t)} \rho f \; dx + \tau(\beta,t) - \tau(\alpha, |
|
t)$$ |
|
|
|
$$\beta'(t)\rho(\beta9t))v(\beta(t)) - |
|
\alpha'(t)\rho(\alpha(t))v(\alpha(t)) = \int_{\alpha(t)}^{\beta(t)} \rho |
|
f \; dx +\int_{\alpha(t)}^{\beta(t)} \frac {\partial \tau}{\partial x}\; |
|
dx$$ |
|
|
|
$\cdots$ |
|
|
|
$$\frac {\partial(\rho v)}{\partial t} + \frac {\partial (\rho |
|
v^2)}{\partial x} = \rho f + \frac {\partial \tau}{\partial x}$$ |
|
|
|
$$\rho v_t + v \rho_t + \rho v_x + v \frac {\partial \rho v}{\partial |
|
x} = \rho f +_ \frac {\partial \tau} {\partial x}$$ |
|
|
|
$\cdots$ |
|
|
|
Then you get the momentum equation as: |
|
|
|
$$\rho \frac {Dv}{Dt} = \rho f + \frac {\partial \tau}{\partial x}$$ |
|
|
|
\begin{center} |
|
Material Derivative $=$ body force $+$ internal stress |
|
\end{center} |
|
|
|
You can think of this as a generalisation of Newton's law with internal |
|
forces in addition to the external forces. |
|
|
|
\subsubsection{In Material Coordinates} |
|
|
|
\begin{itemize} |
|
\item Density $R(A,t) = R_0(A) \left(1 + \frac {\partial U}{\p |
|
A}\right)^{-1} \quad (*)$ |
|
\item Velocity $V(A,t)$ |
|
\item Stress $T(A,t)$ |
|
\item Body force $F(A,t)$ |
|
\end{itemize} |
|
|
|
All considered along the trajectories $X(A,t)$. |
|
|
|
Using the transformation |
|
|
|
$$\frac {\partial \tau}{\partial x} = \frac{\partial T}{\partial |
|
A}\left(\frac {\p U}{\p A} + 1\right)^{-1}$$ |
|
|
|
We need to apply this to the stress, and momentum equation. |
|
|
|
\begin{align*} |
|
R(A,t) \frac {\p V}{\p t} &= RF + \frac {\p T}{\p A}\left(\frac {\p |
|
U}{\p A} + 1\right)^{-1} \\ |
|
R \left(\frac {\p |
|
U}{\p A} + 1\right)\frac {\p V}{\p t} &= R\left(\frac {\p |
|
U}{\p A} + 1\right)F + \frac {\p T}{\p A} |
|
\\ |
|
\marginnote{By the continuum equation in material coordinates $(*)$} |
|
\\ |
|
R_0 \frac {\p V}{\p t} &= R_0 F + \frac {\p T}{\p A} |
|
\end{align*} |
|
|
|
\lecture{} |
|
|
|
\subsection{Constsitutive Laws} |
|
|
|
\subsubsection{Constitutive Law for Elastic Materials} |
|
|
|
Here, stress ($T$) is generated as a result of the displacement of the |
|
material elemeents $U(A,t)$. However $T \ne T(U)$ since then a constant |
|
displacement does not generate zero stress. Stress is generated by the |
|
relative elongation of material elements: \emph{strain}. |
|
|
|
\subsubsection*{Strain in 1D} |
|
|
|
Consider a small section $A, A + \Delta$ at $t = 0$. At time $t$, |
|
|
|
\begin{align*} |
|
A &\to A + U(A, t) \\ |
|
A + \Delta &\to A + \Delta + U(A + \Delta,t) |
|
\end{align*} |
|
|
|
After substituting the relative elongation (local strain) becomes |
|
|
|
$$\frac {l(t) - l_0}{l_0} = \frac {U(A + \Delta, t) - U(A,t)}{\Delta}$$ |
|
|
|
Then since $\Delta \ll 1$ we can use ethe taylor series approximation |
|
for $U(A + \Delta, t)$, $U(A + \Delta, t) \approx U(A) + \frac {\p |
|
U}{\p A} \Delta + \cdots$ |
|
|
|
$$\frac {l(t) - l_0}{l_0} = \frac {U(A + \Delta, t) - U(A,t)}{\Delta} = |
|
\frac {\p U}{\p A}$$ |
|
|
|
Hence in elastic materials the stress is a function of strain |
|
|
|
$$T = f(\frac {\p U}{\p A})$$ |
|
|
|
$f(0) = 0$, but it is otherwise determined by the properties of the |
|
material. |
|
|
|
Assuming a \textbf{linear elastic material}: $T = E\frac {\p U}{\p A}$, where $E$ is a |
|
constant. |
|
\marginnote{We can often choose a region around which a linear |
|
approximation is a good solution} |
|
|
|
Then we get the momentum equation in this form; |
|
|
|
$$R_0 \frac {\p^2 U}{\p t^2} = R_0 F + E \frac {\p^2 U}{\p A^2}$$ |
|
|
|
\emph{Note:} we still don't know about the external forces $F$. |
|
|
|
\begin{ex} |
|
Consider an elastic bungie cord with length $l_0$: how long is it |
|
when stretched by its own weight hanging free in vertical position. |
|
|
|
Initial position $0 < A < l_0$ at $t = 0$. |
|
|
|
Boundary conditions: $U(0,t) = 0$ at the fixed end. |
|
|
|
$T(l_0, t) = 0$ at the free end, $\implies$ $E \frac {\p U}{\p A}(A = l_0, t) = 0$ |
|
|
|
$U(0) = 0$, $U'(l_0) = 0$. |
|
|
|
$$U(A) = - \frac {R_0 g}{2E} + C_1 A$$ |
|
|
|
$$U'(A) = - \frac {R_0 g}E A + C_1$$ |
|
|
|
$$C_1 = \frac {R_0 g}E l_0$$ |
|
|
|
$$U(A) = \frac {R_0 g}{2 E} A(2 l_0 - A)$$ |
|
|
|
final length is $l = l_0 + U(l_0) = l_0 + \frac {R_0 g} {2E} {l_0}^2$ |
|
|
|
We can also calculate the stress distribution along the cord: |
|
|
|
$$T(A) = E U'(A) = R_0 g(l_0 - A)$$ |
|
|
|
The density distribution (non-constant because it is now stretyched.) |
|
|
|
$$R(A) = \frac {R_0} {1 + \frac {R_0 g}E (l_0 - A)}$$ |
|
|
|
Note we are working in material coordinates, if we want them as a |
|
function of spatial coordinates you need to invert $x = A + U(A)$. |
|
|
|
\end{ex} |
|
|
|
|
|
\week{} |
|
\lecture{} |
|
|
|
\section{Review of Continuum Mechanics in 1D} |
|
|
|
\incfig{ref} |
|
|
|
|
|
\subsubsection*{Material coordinate} |
|
|
|
\begin{itemize} |
|
\item $t$ time |
|
\item $A$ cross section |
|
\item $X(A,t)$ a spatial position |
|
\item |
|
$U(A,t) = X(A,t) - X(A,0)$ a displacement at a time from the |
|
\item |
|
$V = \p_t U(A,t) $ velocity |
|
|
|
\item $R(A,t) = \rho(X(A,t), t)$ density |
|
\item $T(A,t) = \tau(X(A,t), t)$ stress |
|
\item $F(A,t) = f(X(A,t), t)$ body force |
|
\end{itemize} |
|
|
|
\marginnote{bad handwriting idk if these are right} |
|
\subsubsection*{Spatial coordinate} |
|
\begin{itemize} |
|
\item $t$ time |
|
\item $a(x,t) = X^{-1}$ |
|
\item $x$ spatial position |
|
\item $u(X,t) = U(a(x,t) ,t)$ displacement |
|
\item $v = D_t u = \frac {\partial_t u}{1 - u_X}$ |
|
|
|
\item $\rho(x,t)$ density ($\frac {mass}{volume}$) |
|
\item $\tau(x,t)$ stress ($\frac {force}{area}$) |
|
\item $f(x,t)$ body force ($\frac {force}{mass}$) |
|
\end{itemize} |
|
|
|
\subsubsection*{Governing equations} |
|
|
|
|
|
$$\p_t f + \p _x(v \rho) = v \quad{\text{continuity equation}}$$ |
|
|
|
$$\rho(\p_t v + v v_x) = \rho f + \p _X \tau = D_t v |
|
\quad{\text{momentum eqn}}$$ |
|
|
|
Have the unknowns $R, V,T$ and only two equations so need a |
|
constitutive equation for $\tau$ or $T$. This can be found using a creep |
|
experiment: stretch a material an amount and see how much force is |
|
requried, or apply an amount of force and see how much it streteches. |
|
|
|
$$\frac{\ell} {\ell_0} = \lambda \quad\quad\text{Extension ratio}$$ |
|
|
|
\begin{itemize} |
|
\item Typically extension grows with stress. |
|
\end{itemize} |
|
|
|
At steady state $\p_{t t} U = 0 \implies \p _A T = 0$ stress $T$ is |
|
constant along the material (from the momentum equation) so we cannot have $T = T(U)$ stress is a |
|
function of displacement. |
|
|
|
|
|
$$ |
|
\begin{array}{|l|l|l|} |
|
\hline \text { Name } & \text { Ratio } & \text { Definition } \\ |
|
\hline \text { Lagrangian Strain } & \left(\ell-\ell_{0}\right) / \ell_{0} & \epsilon=U_{A} \\ |
|
\hline \text { Eulerian Strain } & \left(\ell-\ell_{0}\right) / \ell & \epsilon_{e}=u_{x} \\ |
|
\hline \text { Green Strain } & \left(\ell^{2}-\ell_{0}^{2}\right) /\left(2 \ell_{0}^{2}\right) & \epsilon_{g}=U_{A}+\frac{1}{2} U_{A}^{2} \\ |
|
\hline \text { Almansi Strain } & \left(\ell^{2}-\ell_{0}^{2}\right) /\left(2 \ell^{2}\right) & \epsilon_{a}=u_{x}-\frac{1}{2} u_{x}^{2} \\ |
|
\hline \text { Midpoint Strain } & 2\left(\ell-\ell_{0}\right) /\left(\ell+\ell_{0}\right) & \epsilon_{m}=U_{A} /\left(1+\frac{1}{2} U_{A}\right) \\ |
|
\hline \text { Hencky Strain } & \ln \left(\ell / \ell_{0}\right) & \epsilon_{h}=\ln \left(1+U_{A}\right) \\ |
|
\hline |
|
\end{array} |
|
$$ |
|
|
|
|
|
\lecture{} |
|
Eg For lagrangian strain: $\epsilon \sim \frac {\ell - \ell_0}{\ell_0}$ |
|
|
|
\begin{align*} |
|
\frac{l - l_0}{l_0} &= \frac {X(A + \Delta A, t) - X(A - \Delta t, t) - |
|
2\Delta A - A \times A}{2 \Delta A} \\ |
|
&= \frac {X(A + \Delta A, t) - (A + \Delta A) - (X(A - |
|
\Delta t, t) - (A - \Delta A))}{2\Delta A} \\ |
|
& \quad\quad \text{Taylor expansion} \\ |
|
&= \frac {U(A,t) + \Delta A U_A(A,t) + O(\Delta A^2) - (U(A,t) |
|
- \Delta A(U_A(A,t) + O(\Delta A)^2)) } {2\Delta A}\\ |
|
&= U_A(A,t) + O(\Delta A) \\ |
|
&\implies \lim_{\Delta A \to 0} \frac {\ell - |
|
\ell_0}{\ell_0} = U_A(A,t) =: \epsilon \\ |
|
\end{align*} |
|
|
|
$\implies$ that for elastic materials |
|
|
|
$$T = T(\epsilon) = T(U_A) $$ |
|
|
|
Now the momentum equation is given by |
|
\marginnote{(non-linear) wave equation} |
|
|
|
$$R_0 U_{tt} = R_0 F + T'(U_A) U_{AA}$$ |
|
|
|
$$\lambda = \frac {\ell }{\ell_0} = \frac {\ell - \ell_0 + \ell_0} {\ell_0} \sim |
|
\epsilon + 1 $$ |
|
|
|
$\implies$ the reference state has $\lambda = 1 \Leftrightarrow \epsilon = 0$ |
|
|
|
For a small enough $\epsilon$ any material will exhibit a linear stress--strain |
|
relation. |
|
|
|
$$T(\epsilon) = E\epsilon - E\p_A U$$ |
|
|
|
Wher $E$ is young's modulus also known as the elastic modulus. |
|
|
|
Physical dimensions: |
|
|
|
\begin{itemize} |
|
\item[] $[T] = \frac {force}{area}$ |
|
\item[] $[\epsilon] = [\frac {\ell - \ell_0}{\ell_0}] = 1$ |
|
\item[] $[E] = \frac{force}{area}$ |
|
\end{itemize} |
|
|
|
\begin{center} |
|
\begin{tabular}{l|l|l} |
|
\hline Material & Young's modulus $(\mathrm{GPa})$ & Density $\left(\mathrm{kg} / \mathrm{m}^{3}\right)$ \\ |
|
\hline Diamond & 1000 & 3500 \\ |
|
\hline Stainless steel & 200 & 8030 \\ |
|
\hline Glass & 65 & 2600 \\ |
|
\hline White oak & 12 & 770 \\ |
|
\hline Beeswax & $0.2$ & 960 \\ |
|
\hline Rubber & $0.007$ & 1200 \\ |
|
\hline Silica aerogel & $0.001$ & 100 \\ |
|
\hline |
|
\end{tabular} |
|
\end{center} |
|
|
|
The momentum equation: $R_0 U_{tt} = R_0 F + EU_{AA}$ multiplied by $\frac |
|
1{R_0}$ gives $U{tt} = F + \frac E {R_0} U_{AA}$, which is the linear constant |
|
coefficient wave equation. So we can understand that the characteristic frequencies |
|
are proportional to wave speed $c = \sqrt{\frac E{R_0}}$. For example a high |
|
elastic modulus in the body of a brass instrument has a high wave speed so |
|
expresses those frequencies more. |
|
|
|
Another way of imagining this is if the linear elastic material is imagined as a |
|
chain of beads connected by Hookean springs (at reference length $l_0$), then $F = K(\ell - \ell_0)$ where |
|
$K$ is the spring constant. |
|
|
|
\incfig{linearelasticspring} |
|
|
|
|
|
\begin{align*} |
|
T &=E \frac {\ell - \ell_0} {\ell_0} \\ |
|
&\text{multiply by $\sigma$, cross sectional area} \\ |
|
\sigma T &= \frac { \sigma E}{ \ell_0 } (\ell - \ell_0) \\ |
|
&= K |
|
\end{align*} |
|
|
|
Spring constant is $\frac {\sigma E}{\ell_0}$ |
|
|
|
So we can use \textbf{the creep test} to find young's modulus. |
|
|
|
We need to solve the system of equations |
|
|
|
$$\begin{cases} 0 = EU_{AA} \\ 0 = U(0) \\ U(\ell_0) = \ell - \ell_0 \end{cases}$$ |
|
|
|
\begin{align*} |
|
U &= \alpha A = B \\ |
|
U(0) &= 0 = \beta \\ |
|
U(\ell_0) &= \ell - \ell_0 = \alpha\ell_0 \\ |
|
\implies U(A) = \frac {\ell - \ell_0} A |
|
\end{align*} |
|
|
|
Material density |
|
|
|
\begin{align*} |
|
R &= \frac {R_0}{1 + U_A} \\ |
|
\implies R &= \frac{R_0} {1 + \frac {\ell - \ell_0} {\ell_0}} = \frac {R_0} {\frac |
|
\ell {\ell_0}} = \frac {\ell_0}{\ell} R_0 |
|
\end{align*} |
|
|
|
\lecture{} |
|
|
|
\subsubsection{Bungee rope at steady state} |
|
|
|
\begin{itemize} |
|
\item [] |
|
$\ell$ is the length at rest while hanging: extended undner its own weight. |
|
\item[] |
|
$\ell_0$ is the length of the reference configuration: bungee rope under no |
|
force, at rest |
|
\end{itemize} |
|
|
|
$R_0 U_{t t} = R_0 F + EU_{ AA }$ |
|
|
|
$R_0 U_{t t} = 0$ |
|
|
|
$R_0 F$ is body force |
|
|
|
The body force is obviously gravity so $F = g = 9.81 m/{s^2}$ |
|
|
|
So the momentum equation is $U = R_0 g + EU_{A A}$ |
|
|
|
We dont nkow the value of $U$ at $A = \ell_0$ but we know the boundary |
|
condition: $T(\ell_0) = 0 \implies E U_A(\ell_0) = 0$, and $U(0) = 0$ |
|
|
|
Solve |
|
|
|
$$\begin{cases}0 = R_0 g + EU_{A A}\quad\quad(*)\\ U(0) = 0 \\ U_A(A = \ell_0) = 0 \end{cases}$$ |
|
|
|
THe easiest way is to integrate $(*)$ on $(A, \ell_0)$ to get $0 = R_0 g(\ell_0 |
|
- A) + E(U_A(\ell_0) - U_A(A))$ |
|
|
|
$\implies U_A = \frac {R_0}E g(\ell_0 - A) \quad \quad (* * )$ |
|
|
|
Ie, the tension on the bungee cord decreases linearly along the bungee cord |
|
since $U \propto T$. (This makes sence since as you go down the cord there is |
|
less weight beneath each cross section.) |
|
|
|
\begin{center} |
|
\begin{tikzpicture}{testp} |
|
\begin{axis}[ytick={0}, xtick={0.1},xlabel=$A$, |
|
xticklabels={$\ell_0$},yticklabels={$\frac {R_0}E |
|
g\ell_0$},] |
|
\addplot[domain=0:0.1, samples=90,color=blue,]{ -x / 4}; |
|
\end{axis} |
|
\end{tikzpicture} |
|
\end{center} |
|
|
|
Integrate $(* * ) $ on $(0, A)$ to get |
|
|
|
$U = \frac {R_0} E g (\ell_0 A - \frac {A^2}{2})$ |
|
|
|
$U(\ell_0) = \frac{R_0 g}{E} \frac {{\ell_0}^2} 2$ |
|
|
|
\begin{center} |
|
\begin{tikzpicture}{hy} |
|
\begin{axis}[ytick={1}, yticklabels={$\frac{R_0 g}{E} \frac |
|
{{\ell_0}^2} 2$}, xtick={0,2},xlabel=$A$, |
|
xticklabels={0, $2\ell_0$}, |
|
] |
|
\addplot[domain=0:2, samples=90,color=blue,]{ -(x -1)^2 + 1}; |
|
\end{axis} |
|
\end{tikzpicture} |
|
\end{center} |
|
|
|
So clearly $\ell$ grows like ${\ell_0}^2$ because of the cord's own weight. |
|
|
|
|
|
|
|
\section{Hyperelastic Materials} |
|
|
|
\begin{itemize} |
|
\item Nonlinear stress-strain relations. |
|
\item Energy balance |
|
\end{itemize} |
|
|
|
|
|
For example for rubber the stress strain raltion looks something like |
|
|
|
|
|
|
|
\begin{figure}[htpb] |
|
\centering |
|
\begin{tikzpicture}{testppop} |
|
\begin{axis}[ytick={},xtick={},xlabel=$\epsilon$, ylabel=$T$,] |
|
\addplot[domain=-1:2, samples=90,color=blue,]{ (x -0)^3 + 0}; |
|
\end{axis} |
|
\end{tikzpicture} |
|
\caption{Stress-Strain for rubber} |
|
\end{figure} |
|
|
|
\begin{align*} |
|
T(\epsilon) &= \psi'(\epsilon) \text{ where } \\ |
|
\psi(\epsilon) &:= \text{ Strain energy density} \\ |
|
\text{For example} &T(\epsilon) = K\epsilon^3 \\ |
|
T(\epsilon) &= K\epsilon^3 \implies \psi(\epsilon) = K \frac {\epsilon^4}4 \\ |
|
T(\epsilon) &= E\arctan(\epsilon) \\ |
|
T(\epsilon) &= E \epsilon + K \epsilon^3 \\ |
|
&\implies \psi(\epsilon) = E \frac {\epsilon^2} 2 + K \frac |
|
{\epsilon^4} 4 |
|
\end{align*} |
|
|
|
Implies the momenntum equation |
|
|
|
\begin{align*} |
|
R_0 U{z z} &= R_0 F + \p_A \psi'(U_A) \\ |
|
&= R_0 F + \psi'' (U_A) U_{A A} |
|
\end{align*} |
|
|
|
Typically need that $\psi'' > 0$ |
|
so that the problem is well posed: $T(\epsilon) = \psi'(\epsilon) $ is monotone |
|
in $\epsilon$. |
|
|
|
\subsubsection*{Energy balance} |
|
|
|
For a $A_\ell $ and $A_r$. (note $\sigma$ is the cross sectional area). |
|
|
|
Multiply momentum equatino by $\sigma U_t$. |
|
|
|
\begin{align*} |
|
\sigma R_0 U_t U_{tt} &= \sigma R_0 F V + \sigma \p_A \psi'(U_A) U_z \\ |
|
\text{integrate over} A_\ell \text{ to } A_r \\ |
|
\sigma R_0 \int_{A_\ell}{A_r} U_t U_{T t} \: dA &= \sigma R_0 |
|
\int_{A_\ell}^{A_R} FV \: dA + \sigma \int_{A_\ell} ^{A_r} \p_A \psi'(U_A) |
|
U_z \: dA \\ |
|
&\text{Now integrate by parts} |
|
\\ |
|
\sigma R_0 \int_{A_\ell}^{A_r} \frac d {dt} \frac {U_t^2} 2 \: da &= \sigma |
|
R_0 \int_{A_\ell} ^{A_r} FV \:dA + r[\psi'(U_A) U_t]_{A_\ell}^{A_R} - \sigma |
|
\int_{A_l}^{A_r} \psi'(U_A) U_{At} \: da \\ |
|
... \\ |
|
\implies |
|
\end{align*} |
|
$$\frac d{dt} \int_{A_\ell}^{A_r} \left[\sigma R_0 \frac {V^2}{2} + \sigma |
|
\psi(U_A)\right] \: dA = \sigma \left(R_0 \int_{A_\ell}^{A_r} FV \: da + |
|
TV|_{A_r} - TV |_{A_\ell}\right) $$ |
|
|
|
Where we can see |
|
|
|
\begin{itemize} |
|
\item |
|
$\sigma R_0 \frac {V^2}{2}$ : kinetic energy |
|
\item |
|
$ \sigma \psi(U_A)$ : strain energy |
|
\item |
|
$R_0 \int_{A_\ell}^{A_r} FV \: da $ : Power delivered by body forces |
|
\item |
|
$TV|_{A_r} - TV |_{A_\ell} $ : Power delivered by force acting on the tips |
|
|
|
\end{itemize} |
|
|
|
\section{Time-dependent elastic models} |
|
|
|
For a linearly elastic material: |
|
|
|
$$\begin{cases} R(A,t) = \frac {R_0} {1 + U_A} \quad(i)\\ |
|
R_0 U_{t t} = R_0 F + EU_{AA}\quad(ii)\end{cases}$$ |
|
|
|
These are not coupled. You should solve $(ii)$ then $(i)$ |
|
|
|
From math2100: The wave equation on $\R$ |
|
|
|
$$(*) \begin{cases} |
|
U_{t t} - c^2 U_{A A} = 0\\ |
|
U(t = 0, A) = f(A)\\ |
|
\p_t U(t = 0, A) = \rho(A) |
|
\end{cases}$$ |
|
|
|
You can get the solution using the d'Alembert solution to the wave equation. |
|
This superimposes two waves, one that moves right and one that moves left at |
|
the same speed. |
|
|
|
$$U(t,a) = \frac 12 (f(A - ct) + f(A + ct)) + \frac 1 {2c} \int_{A - ct}^{A + |
|
ct} g(\hat A) \: d\hat A$$ |
|
|
|
Remark: This formula also gives a solution to $(*)$ on an interval with derichlet zero |
|
boundary conditions ($U(t, A = 0) = 0 = U(t, A = \ell)$: when there is no |
|
|
|
displacement at the endpoints of $(0,\ell)$). |
|
|
|
To get this initial condition you can extend $f,g$ as odd $2\ell$ periodic |
|
functions. |
|
|
|
\begin{figure}[htpb] |
|
\centering |
|
\incfig{exten} |
|
\caption{Extended function that is periodic on $2\ell$} |
|
\end{figure} |
|
|
|
\begin{ex} |
|
|
|
With a slinky: |
|
|
|
$A = 0$ and $A = 10 = \ell_0$ |
|
|
|
Intially move the loop of $A = 4$ to the left by $\Delta A = 2$, then |
|
release it. |
|
|
|
\incfigw{exspr}{0.5} |
|
|
|
How do you write $U(t = 0, A) = U_0(A)$? |
|
|
|
\end{ex} |
|
|
|
\week{} |
|
\lecture{} |
|
|
|
We know |
|
\begin{align*} |
|
U_0(A = 0) &= 0 \\ |
|
U_0(A = 10) &= 0 \\ |
|
U_0(A = 4) &= -2 \\ |
|
\end{align*} |
|
|
|
At $t = 0$ the slinky satisfies $EU_{A A} = 0 \implies U_A$ is a constant. |
|
|
|
\begin{align*} |
|
U_A &= \frac {-2 - 0}{4 - 0} = -\frac 1 2 \implies U = -\frac 12 A \\ |
|
U_A &= \frac {0 - (-2)}{10 - 4} = \frac 13 \implies U = \frac 13 (A - 10) \\ |
|
U_0 &=\begin{cases} -\frac 12 A & 0 \le A \le 4 \\ \frac 13 (A - 10) & 4 < A |
|
\le 10\end{cases} \\ |
|
\end{align*} |
|
The time-dependent problem is given by |
|
|
|
\begin{align*} |
|
U_{t t} &= \frac E{R_0} U_{A A} \\ |
|
U(0) &= 0 = U(10) \\ |
|
U(t = 0, A) &=\begin{cases} -\frac 12 A & 0 \le A \le 4 \\ \frac 13 (A - 10) |
|
& 4 < A \le 10 \end{cases} \\ |
|
\p_t U(t = 0, A) &= 0 |
|
\end{align*} |
|
|
|
$\implies$ the solution is given by D'Alembert's formula where $\rho=0$ and $f$ |
|
is the $2\ell$ periodic odd extension of $U_0$. |
|
|
|
|
|
\begin{ex} |
|
Finite-length bar under external forcing (resonance) |
|
|
|
$U_{t t} = \frac E {R_0} U_{A A} + F(A,t) $ on the interval $(0,\ell)$. |
|
|
|
$U_A(t = 0, A = 0) = 0 = U_A(t = 0, A = \ell) $ |
|
|
|
$F(A,t) = a(t) \cos(\kappa_n A) $ where $a(t) = \sin(\omega t)$ |
|
|
|
$\omega$ is angular frequency of the forcing. |
|
|
|
$\kappa_n = \frac {n \pi}\ell$ where $n$ is the wave number of the specific |
|
mode of the forcing term. |
|
|
|
\begin{center} |
|
\begin{tikzpicture}{notes-figure9} |
|
\begin{axis}[xtick={3.141},xticklabels={$\ell$}, ytick={},] |
|
\addplot[domain=0:3.141 , samples=90, color=green]{ cos(deg(x))}; |
|
\addlegendentry{$n = 1$} |
|
\addplot[domain=0:3.141 , samples=90, color=black]{ cos(2 * deg(x))}; |
|
\addlegendentry{$n = 2$} |
|
\addplot[domain=0:3.141 , samples=90, color=blue]{ cos(3 * deg(x))}; |
|
\addlegendentry{$n = 3$} |
|
\end{axis} |
|
\end{tikzpicture} |
|
\end{center} |
|
|
|
|
|
We can solve with Laplace transform to turn the pde into an ODE, transform |
|
the entire system of equations to get rid of the time derivatives |
|
|
|
|
|
$$ |
|
F(s)=L(f(t))=\int_{0}^{\infty} e^{-s t} f(t) d t |
|
$$ |
|
|
|
and the opposite |
|
|
|
$$ |
|
L(t f(t))=-\frac{d F(s)}{d s} |
|
$$ |
|
|
|
\begin{align*} |
|
(*) \begin{cases} |
|
s^2 \hat U - c^2 \hat U_{A A} &= \hat F (A,s) = \hat a(s) \cos(\kappa_n A) \\ |
|
&\text{Boundary conditions:} \\ |
|
\hat U_A(s, A=0) &= 0 = \hat U_A(s, A = \ell) |
|
\end{cases} |
|
\end{align*} |
|
|
|
To solve $(*)$ first consider the homogeneous equation: $\hat U_{hom} - \frac |
|
{c^2}{s^2} \hat U_{A A} = 0$ |
|
|
|
$$\hat U_{hom} = \alpha \exp(\frac cs A) + \beta \exp(-\frac cs A)$$ |
|
|
|
Find a particular solution : |
|
|
|
$$\hat U _{part} = \gamma \sin(\kappa_n A) + \delta \cos(\kappa_n A)$$ |
|
|
|
|
|
\begin{align*} |
|
\implies \p_{A A} U_{part} &= -\kappa_n^2 (\gamma \sin(\kappa_n A) + \delta |
|
\cos(\kappa_n A)) \\ |
|
s^2 \hat U_{part} - c^2 \p_{A A} \hat U_{part} &= (s^2 + c^2 |
|
\kappa_n^2)(\sin(\kappa_n a) + \delta \cos(\kappa_n A)) \\ |
|
&= \hat a(s) \cos(\kappa_n A) \\ |
|
\implies \gamma = 0 \text{ and} \delta = \frac {\hat a(s)}{ s^2 + c^2 |
|
\kappa_n^2} |
|
\end{align*} |
|
|
|
$\implies $ the general solution of $(*)$ is |
|
|
|
$$ |
|
\hat U = \frac {\hat a(s)}{ s^2 + c^2 |
|
\kappa_n^2} \cos(\kappa_n A) + \alpha \exp(\frac cs A) + \beta \exp(-\frac cs A) |
|
$$ |
|
|
|
Boundary conditions: |
|
|
|
\begin{align*} |
|
\hat U_A &= -\kappa_n \frac {\hat a}{s^2 + \kappa_n^2 c^2} \sin(\kappa_n A ) |
|
+ \frac cs (\alpha \exp(\frac cs A) - \beta \exp(-\frac cs A)) \\ |
|
U = \hat U_A(A = 0) &= \frac cs (\alpha - \beta) \implies \alpha = \beta \\ |
|
U = \hat U_A(A = \ell) &= \frac cs (\alpha \exp(\frac cs \ell) - \beta |
|
\exp(-\frac cs \ell))\alpha \\ |
|
& \text{since }\ell \ne 0, \; \alpha \exp(\frac cs \ell) - \beta |
|
\exp(-\frac cs \ell) |
|
\implies \alpha = 0 \\ |
|
&\implies \hat U = \frac {\hat a} {s^2 \kappa_n^2 |
|
c^2} \cos (\kappa_n A) |
|
\end{align*} |
|
|
|
\lecture{} |
|
|
|
Now we compute the inverse laplace transform using the convolution theorem. |
|
|
|
\begin{align*} |
|
\hat U &= \frac {\hat a} {s^2 \kappa_n^2 c^2} \cos (\kappa_n A) \\ |
|
&\text{We have to use the convolution theorem} \\ |
|
L^{-1}(\hat U) &= \cos(\kappa_n A) L^{-1}(\frac 1 {s^2 + \kappa_n^2 |
|
c^2}) \times u(t) |
|
\quad\quad \text{(note $\cos(\kappa_n A)$ does not depend on $s$)} \\ |
|
L^{-1}(\frac 1 {s^2 + \kappa_n^2 |
|
c^2}) \times u(t) &= \frac 1 {k_n c} \int _0^t \sin(\omega(t - \hat t)) |
|
\sin(\kappa_n c \hat t) \: :d\hat t \\ |
|
&= \frac 1 {\kappa_n c} \begin{cases} \frac 1 {\omega^2 |
|
\kappa_n^2 c^2}(\omega \sin(\kappa_n ct) - \kappa_n c |
|
\sin(\omega t)) & \omega \ne \kappa_n c\\ |
|
-\frac 12 |
|
\mathcolorbox{yellow}{t} |
|
\cos(\kappa_n ct) + \frac 1 {2\kappa_n c} |
|
\sin(\kappa_n ct) & \omega = \kappa_n c |
|
\end{cases} \\ |
|
\end{align*} |
|
|
|
$\implies \omega = \kappa_n c$ is the resonance cases; the amplitude of $U$ |
|
grows linearly in time (due to the highlighted \mathcolorbox{yellow}{$t$}). The |
|
reason is $\kappa_n c$ are the eigenmodes of the structure. |
|
|
|
|
|
\marginnote{Remember $R_0$ is density, $E$ is elastic modulus.} |
|
|
|
Since the eigenfrequenceies $\kappa_n c$ are proportional to $c = \sqrt{\frac E |
|
{R_0}}$ one can probe the material propertiy $\sqrt{\frac E {R_0}}$ through an |
|
acoustic signal (frequancies of resonances in the material will be $\kappa_n |
|
\sqrt{\frac E {R_0}}$). |
|
|
|
\end{ex} |
|
|
|
\subsection{Viscoelastic Materials} |
|
|
|
\begin{itemize} |
|
\item Introducing damping to the stress-strain relation. |
|
\end{itemize} |
|
|
|
\begin{ex} |
|
|
|
For example, considering the slinky again. It is unrealistic for it to |
|
oscillate forever, instead some sort of damping, eg air-external friction, |
|
internal friction (resistance to deformation) will slow it down. |
|
|
|
If you have a spring hanging freely with a weight and it is in equilibrium, you |
|
have the spring force $F_s = -K_s u = mg$, and a damping force $F_d-c\dot u = 0$. |
|
When it moves the damping force creates a resistance. Resitive force is proportional |
|
to the velocity: $F_d = -c \dot u$. |
|
|
|
Newtond's Law: $m \ddot u = F_d + F_s$ gives the damped oscillation equation |
|
|
|
$$m \ddot u + c\dot u + K u = 0$$ |
|
|
|
These systems provide a framework to model viscoelasticity. |
|
\end{ex} |
|
|
|
\subsubsection{Kelvin-Voigt element} |
|
|
|
You have a (horizontally) a spring and dashpot in parallel. |
|
|
|
$$F = - (F_s + F_d)$$ |
|
|
|
So displacements are equal: $u_s = u_d = u$. |
|
|
|
Spring |
|
constant $K$. |
|
|
|
$$F_s = -K u$$ |
|
$$F_d = -c \dot u$$ |
|
|
|
To create a constitutive law for stress, identify $F$ by stress $T$ and $u$ b |
|
strain $\epsilon$. Also replace $\dot F$ by $\partial_t T$ (stress rate), and |
|
$\dot u$ by $\p_t \epsilon$ (strain rate). |
|
|
|
$$T = E(\epsilon + \tau_1 \p_t \epsilon ) $$ |
|
|
|
$E$ elastic modulus |
|
|
|
$\tau_1$ dissipation time of the strain. |
|
|
|
\subsubsection{Maxwell element} |
|
|
|
Spring and dashpot in series. |
|
|
|
Displacements sum to $u$: $u = u_s + u_d$ |
|
|
|
Forces are equal: $F = - F_s = -F_d$ |
|
|
|
Constitutive law: |
|
|
|
$F = -F_s= ku_s \implies u_s = \frac FK \implies \dot u_s = |
|
\frac {\dot F }K $ |
|
|
|
$F = -F_d = c \dot u_d \implies \dot u_d = \frac FC$ |
|
|
|
So |
|
$\dot u = \dot u _s + \dot u_d = \frac {\dot F}k + \frac F c$ |
|
|
|
$(* * ) \quad\quad c \dot u = \frac ck \dot F + F$ |
|
|
|
Which gives the constituive law |
|
|
|
$$T + \tau_0 \p_t T = E\tau_1 \p_t \epsilon$$ |
|
|
|
|
|
\lecture{} |
|
|
|
\subsubsection{Standard Linear Model} |
|
|
|
\incfigw{dsprlin}{0.7} |
|
|
|
To model the spring and despot in series use the resulf for the Maxwell element |
|
$(* * ) $. |
|
|
|
$$c_1 u' = \frac {c_1}{k_1} F_1' + F_1 = \frac {c_1}{k_1} (F' - k_0 u') + F - |
|
k_0 u$$ |
|
|
|
$$\implies k_0 u + u(c_1 + \frac {k_0 }{k_1} c_1) = F + \frac {c_1}{k_1} F'$$ |
|
|
|
Rewrite this as |
|
|
|
$$F + a_1 F' = a_2 u + a_3 u'$$ |
|
|
|
Where $a_1 = \frac {c_1} {k_1}$, $a_2 = k_0$ and $a_3 = c_1 + \frac |
|
{k_0}{k_1}c_1$ so that $a_3 = c_1 + \frac {k_0}{k_1} c_1 > \frac {k_0}{k_1} c_1 |
|
= a_1 a_2$ |
|
|
|
This gives the constitutive law |
|
|
|
$$T + \tau_0 T_t = E(\epsilon + \tau_1 \p_t \epsilon)$$ |
|
|
|
\begin{align*} |
|
\tau_0 &= a_1 = \frac {c_1}{k_1} \\ |
|
E &= a_2 = k_0 \\ |
|
E T_1 &= a_3 |
|
\end{align*} |
|
|
|
And as we showed |
|
|
|
$ET_1 = a_3 > a_1 a_2 = ET_0$ |
|
|
|
$E = a_2 = k_0$ |
|
|
|
Which $\implies \tau_1 > \tau_0$, where $\tau_0$ and $\tau_1$ are dissipation |
|
times, $E$ the elastic modulus. |
|
|
|
|
|
$$\begin{cases} |
|
R = \frac {R_0}{1 + U_A} \\ |
|
R_0 U_{t t} = R_0 F + T_A \\ |
|
T + \tau_0 T_t = E(U_A + \tau_1 U_{At}) |
|
\end{cases}$$ |
|
|
|
\subsubsection{Summary} |
|
|
|
\emph{Standard linear model} |
|
$$T + \tau_0 \p_t T = E(U_A + \tau_1 \p_t \epsilon)$$ |
|
|
|
\emph{Kelvin-Voigt model} |
|
$$T = E(\epsilon + \tau_1 \p_t \epsilon ) $$ |
|
|
|
Note the above two models involve the strain (in the right hand side), so they will eventually go back to |
|
the reference configuration. These are viscoelastic solids. |
|
|
|
\emph{Maxwell Model} |
|
$$T + \tau_0 \p_t T = E\tau_1 \p_t \epsilon$$ |
|
|
|
This is a visco-elastic fluid since it only includes the strain rate, it will |
|
not return to the reference configuration. It is elastic on fast time scales, |
|
but on a long timescale you have a creeping motion that deforms the fluid. |
|
|
|
\subsubsection{Relaxation Function} |
|
|
|
This unifies the above three models, by adding constant $C$ to the standard linear model. |
|
|
|
$$T + \tau_0 \p_t T = E(C\epsilon \tau_1 \p_t \epsilon)$$ |
|
|
|
\begin{itemize} |
|
\item[] $\tau_0 =0, C = 1$ Kelvin-Voigt model |
|
\item[] $C = 0$ Maxwell model |
|
\item[] $C = 1$ Standard linear model |
|
\end{itemize} |
|
|
|
Solving |
|
$$\begin{cases} T + \tau_0 \p_t T = h \\ |
|
T(t=0) = 0 |
|
\end{cases}$$ |
|
|
|
where $h = E(C\epsilon + \tau_1 \p_t \epsilon)$, (and assuming $\epsilon(0) = 0$) |
|
|
|
|
|
$$T = \int_0^t E(C + (\frac {\tau_1}{\tau_0} - c)e^{\frac {-(t - |
|
\hat t)}{\tau_0}})\p_t \epsilon{\hat t}\: d\hat t$$ |
|
|
|
Where the relaxation function |
|
|
|
$$G(t - \hat t) = E(C + (\frac {\tau_1}{\tau_0} - c)e^{\frac {-(t - |
|
\hat t)}{\tau_0}})$$ |
|
|
|
\begin{itemize} |
|
\item Standard linear model |
|
$$G(t) = E(1 + (\frac {\tau_1}{\tau_0} - 1) e^{\frac {-t}{\tau_0}})$$ |
|
\item Maxwell element |
|
|
|
$$G(t) = E \frac {\tau_1}{\tau_0} e^{\frac {-t}{\tau_0}}$$ |
|
|
|
\item Kelvin-Voigt |
|
|
|
Note cannot simpoly sat $\tau_0$ to zero but it is attainable using a |
|
different method. |
|
|
|
$$T = E(\epsilon + \tau_1 \p_t \epsilon) = \int_0^1 E(\p_t \epsilon(\hat |
|
t) + \tau_1 \rho_0 (t - \hat t) \p _t \epsilon(\hat t) ) \: d \hat t |
|
\implies$$ |
|
|
|
$$G(t) = E(1 + \tau_1 \delta_0(t))$$ |
|
\end{itemize} |
|
|
|
A general formulation for constitutive relation for stress ($T$) is |
|
|
|
$$T = \int_0^t G(t - \hat t) \epsilon_t (\hat t) \: d\hat t = G * \p_t \epsilon $$ |
|
|
|
Here $G \ge 0$ can be choses quite arbitrarily (want it to be integrable, decay |
|
to zero). |
|
|
|
This is the ``visco-elastic law of relaxation type''. |
|
|
|
\week{} |
|
\lecture{} |
|
|
|
Instead you may also write the stress as an integral against the strain and |
|
integrate by parts in $T_2$ using $\epsilon (0) = 0$: |
|
|
|
$$T_2 = \int_0^t E \frac {T_1}{T_0} E^{-\frac {t - \hat t}{\tau_0}} \p_t |
|
\epsilon(\hat t) \: dt = |
|
E\frac {T_1}{T_0} \epsilon(t) - \int_0^t E\frac {T_1}{T_2} e^{-t - \hat t} |
|
{\tau_0} \epsilon(\hat t) \: d \hat t$$ |
|
|
|
For example for the standard linear model ($C = 1$) this implies |
|
|
|
$$T = T_1 + T_2 = \underbrace{E \frac {T_1}{T_0} |
|
\epsilon(t)}_{\begin{minipage}[adjusting]{5em}instantaneous elastic stress \end{minipage}} |
|
- \underbrace{\frac E {\tau_0} \int_0^t \left( \frac {\tau_1}{\tau_0} - 1\right ) e^{\frac |
|
{t - \hat t}{\tau_0}} \epsilon(\hat t) \: d \hat |
|
t}_{\begin{minipage}{5em}damping\end{minipage}} |
|
$$ |
|
|
|
{The computations clearly get complicated so often you can't solve |
|
viscoelastic models explicitly very easiliy} |
|
|
|
\begin{ex} |
|
Maxwell-viscoelastic bar occupying $0 \le A < \infty$ with body force $F = |
|
0$. |
|
|
|
Momentum equation: |
|
|
|
$$(*) \qquad\quad R_{0} U_{t t} - T_A = \p_A (G \star \epsilon_t) = \int_0^t G(t - \hat t) |
|
U_{AA_t}$$ |
|
|
|
Where $G(t) = E \frac {\tau_1}{\tau_0} e^{\frac {-t} {\tau_0}}$. |
|
|
|
We have the ``far field condition'' on the right, the displacement vanishes at |
|
infinity: |
|
|
|
$$\lim_{A \to \infty} U(A,t) = 0$$ |
|
|
|
$$U(t=0, A) = 0 = \p_t U(t = 0, A)$$ |
|
|
|
How to enforce $U_A(t , A=0) = F(t) $ experimentally? |
|
|
|
Apply a force to the left endpoint: $T(A = 0, t) = \int_0^t G(t - \hat t) F_t(\hat t) \: d\hat t $ |
|
|
|
Use laplace transform on $(*)$, $\hat U = L(U)$ |
|
|
|
$$\begin{cases} |
|
R_0 s^2 \hat U = L(G) s \hat U_{A A} $ where $L(G) = \frac {\hat E}{s + |
|
\frac 1 {\tau_0}} \\ |
|
\hat U_A (A = 0, s) = \hat F(s) \\ |
|
\lim_{A \to \infty} \hat U(A, s) = 0 |
|
\end{cases}$$ |
|
|
|
|
|
\begin{align*} |
|
R_0 s \hat U = \frac {\hat E}{s + \frac 1 {\tau_0} } \hat U_{A A} |
|
&\implies \frac {R_0 s}{\hat E}(s + \frac 1 {\tau_0}) \hat U = \hat U |
|
_{A A} \\ |
|
&\implies \hat U = \alpha (s) e ^{\omega (s) A} + \beta(s) e^{-\omega(s) A} |
|
\qquad {\text{ $($ODE for }\hat U)}\\ |
|
\text{where } \omega(s) &=\sqrt{\frac {R_0} s {\hat E} (s + \frac 1 |
|
{\tau_0})} \\ |
|
\end{align*} |
|
|
|
\marginnote{laplace transform only works when the real part of $s$ is large |
|
enough compared to the imaginary part.} |
|
|
|
Note that $Re(\omega) > 0$ if $Re(s) $ is large enough. |
|
Therefore $\alpha(s) = $ to guarantee that $\lim_{A \to \infty} \hat U(A) = |
|
0$. |
|
|
|
Find the other constant: |
|
|
|
\begin{align*} |
|
\implies & \hat U(A,s) = \beta(s) e^{-\omega(s) A} \\ |
|
\implies & \hat U_A(A,s) = -\omega(s) \beta(s) e^{-\omega(s) A} \\ |
|
{\text{2nd boundary cond. }} \hat F(s) &= \hat U_A(A = 0, s) \\ |
|
&=-\omega9s) \beta(s) \\ |
|
\implies \beta(s) &= \frac {-\hat F(s)}{\omega(s)} \\ |
|
\implies \hat U(A, s) &= -\frac 1 {-\omega (s)} \hat F(s) e^{-\omega(s) A} \\ |
|
\end{align*} |
|
|
|
Now using the convolution theorem to inverse laplace transform: |
|
|
|
\begin{align*} |
|
U(A,t) &= F(t) \star \underbrace{L^{-1}\left(-\frac 1 {-\omega (s)} e^{-\omega(s) A} |
|
\right)}_{(* *)} \\ |
|
(* *) &= \sqrt{\frac {\hat E}{R_0}} |
|
e^{\frac {-t}{\tau_0}}I_0 \left(\frac {-1}{2\tau_0} \sqrt{t^2 - \frac {R_0}{\hat E} |
|
A^2}\right) \times H\left(t - \sqrt{\frac {R_0}{\hat E}} A\right)\\ |
|
I_0 & {\text{ is the Bessel function of the 1st kind }} \\ |
|
H & \text{ is the heaviside function } \\ |
|
U &= \int_0^t Q(A,\hat t) H(\hat t - \sqrt{\frac {R_0}{\hat E}} A) F (t - \hat t) |
|
\hat H \\ |
|
&=\int_{\sqrt{\frac {R_0}{\hat E}}A}^t Q(A, \hat t) F(t - \hat t) \: d\hat t \times H\left(t - \sqrt |
|
{\frac{R_0}{\hat E}}A\right) \\ |
|
\implies & U \text{ is zero wherever } t < \sqrt{\frac {R_0}A} , \\ |
|
\implies & A > t\sqrt {\frac {\hat E} {R_0}} = t \sqrt{\frac {E \tau_1}{\tau_0} |
|
\frac 1 {R_0}} \\ |
|
\text{For } F(t) &= \sin(t) \\ |
|
\end{align*} |
|
|
|
The speed of propagation in the linearly elastic material $c = \sqrt{\frac E {R_0}}$ is smaller |
|
than the propagation in the viscoelastic model $\sqrt{\frac E{R_0} \frac {\tau_1}{\tau_0} }$ |
|
$(\tau_1 > \tau_0)$. And for $\tau_0 \to 0$ (no memory, purely viscous) the speed of propagation |
|
become infinite just like in the diffusion equation. |
|
\end{ex} |
|
|
|
\lecture{} |
|
|
|
\subsubsection{Dynamic (complex) modulus} |
|
|
|
\begin{itemize} |
|
\item Dynamic mechanical analysis |
|
\item ie ``The {Vibration test}'' |
|
\end{itemize} |
|
|
|
A bar of infinite length connected by a rod to a wheel that forces the bar periodically. We are |
|
aiming to determine the properties of the bar, so do not assume anything about its elastic |
|
properties. |
|
|
|
|
|
$U(A = 0, t) = a \sim(\omega t)$ apply a sinusoidal displacement at the tip. |
|
|
|
$\lim_{A \to \infty} U(A,t) = 0$. |
|
|
|
$T + \tau_0 \p_t T = E(\epsilon + \tau_1 n\p _t \epsilon) = E(U_a + \tau_1 U_{At}) \qquad (*)$ |
|
|
|
$R_0 U_{t t} = T_a \qquad(* *)$ |
|
|
|
it is easier to write everything in complex coordinates: |
|
|
|
$U(A = 0, t) = e^{i \omega t}$, $(\sin(\omega t))$ is the complex part |
|
|
|
$U(A,t) = \bar U(A) e^{i \omega t}$ |
|
|
|
$T(A,t) = \bar T (A) e^{i \omega t}$ |
|
|
|
Substitute in $(*)$ |
|
|
|
\begin{align*} |
|
\bar T(A) e^{i \omega t} (1 + i\omega \tau_0)) &= E \bar U_A (A) e^{i \omega t} (1 + \tau_1 i |
|
\omega) \\ |
|
\bar T(A) &= E \frac {1 + i \omega \tau_1}{1 + i \omega \tau_0} \bar U_A(A) = E^* \bar U_A(A) |
|
\end{align*} |
|
|
|
Which is the complex stress-strain relation under oscillatory conditions (because it looks very |
|
similar to the linear law for it but in terms of complex numbers). |
|
|
|
\begin{align*} |
|
E^* &= E \frac {1 + i\omega \tau_1}{1 + i \omega \tau_0} = E \frac {(1 + i \omega \tau_1)(1 - i |
|
\omega \tau_0)}{(1 + i \omega \tau_0)(1 - i \omega \tau_0)} \\ |
|
&= E \frac {1 + \omega^2 \tau_0 \tau_1}{1 + \omega^2 {\tau_0} ^2} + i E \frac {\omega (\tau_1 |
|
- \tau_0)}{1 + \omega^2 {\tau_0}^2} \\ |
|
&= E' + i E'' \\ |
|
E': & \text{ storage modulus (`the elastic part') } \\ |
|
E'': & \text{ loss modulus (`the viscous part') } \\ |
|
& \text{from the momentum equation } (* * ) : R_0 U_{t t} = T_A \\ |
|
R_0 \bar U(A) (-\omega^2) e^{i \omega t} &= \bar T_A(A) e^{i \omega t} |
|
= E* \bar U_{A A} (A) e^{i \omega t} \\ |
|
\frac {-\omega^2 R_0}{E^*} \bar U(A) &= \bar U_{A A} (A) \\ |
|
\text{ make } \kappa &= \sqrt{\frac {\omega^2 R_0} {E^*}} \\ |
|
&\text{The general solution is: } |
|
\bar U(A) = \alpha e^{i \kappa A} + \beta e^{i \kappa A} \\ |
|
\kappa^2 &= \frac {\omega^2 R_0} E \left(\frac {1 + \omega^2 \tau_0 \tau_1}{1 + \omega^2 |
|
\tau_1^2} + i \frac {\omega(\tau_0 - \tau_1)} {1 + \omega^2 \tau_1^2}\right) \\ |
|
Im(\kappa) < 0 &\implies Re (i\kappa) > 0\\ |
|
& \implies \alpha = 0 \text{ to satisfy } \lim_{A \to \infty} \bar U(A) = 0 \\ |
|
U(A = 0, t) &= \alpha e^{i \omega t} \\ |
|
& \implies \bar U(0) e^{i \omega t} = \beta e^{-i \kappa 0} e^{i \omega t} |
|
= \beta e^{i \omega t} \\ |
|
& \implies a = \beta \\ |
|
&\implies \bar U(A) = a e^{i \kappa A} |
|
\end{align*} |
|
|
|
Now want the stress at the endpoint $A = 0$ |
|
|
|
\begin{align*} |
|
T(A = 0, t) &= \bar T(A = 0) e^{i \omega t} = E^* \bar U_A(A = 0) e^{i \omega t} \\ |
|
&= E^*(-i) \kappa a e^{i \omega t}\\ |
|
\end{align*} |
|
|
|
\begin{itemize} |
|
\item You increase the freqency forcing and for every frequency you measure the phase |
|
difference. |
|
\end{itemize} |
|
|
|
|
|
\lecture{} |
|
|
|
What is the displacement $U$ and stress $T$ at the left end point. |
|
|
|
$U(A = 0, t) = \bar U(A = 0) e^{i \omega t} a e^{i \kappa 0}e^{i \omega t} = a e^{i \omega _t}$ |
|
|
|
$T(A = 0, t) = \bar T(A = 0) e^{i \omega t} = E^* (-i) \kappa a e^{i \omega t}$ |
|
|
|
$T$ In polar coordinates: |
|
\begin{align*} |
|
-i &= e^{-\frac \pi 2}, \\ |
|
\kappa a e^{i \omega t} &= \kappa a e^{i \omega t} \quad\text{and, } \\ |
|
E^* &= re^{i \delta} = E^* \sqrt{R_0 \frac {\omega^2}{E^*}} = \sqrt{R_0 \omega^2 |
|
E^*} \\ |
|
r^2 e^{2 i \delta} &= R_0 \omega^2 E^* = R_0 \omega^2(E' + iE'') \\ |
|
& \text{take the argument } 2\delta = \arctan \frac {E''}{E'} \\ |
|
\implies \delta &= \frac 12 \arctan {\frac {\omega(\tau_1 - \tau_0)}{1 + |
|
\omega^2 \tau_0 \tau_1}} \\ |
|
\implies T(A = 0, t) &= a r \exp(i (\omega t - \frac \pi 2 + \delta)) |
|
\end{align*} |
|
|
|
Phase difference: $\delta - \frac \pi 2$ |
|
|
|
\subsubsection*{Interpretation} |
|
|
|
\begin{enumerate} |
|
\item |
|
If you have an elastic material, then $\tau_1 = \tau_0 = 0$, so $\delta = 0$ and |
|
the phase difference is $-\frac \pi 2$ independently of the frequency $\omega$ |
|
\item |
|
if you have $\tau_0 = 0$ (the Kelvin-Voigt model of a viscoelastic solid) $\delta = \frac 12 |
|
\arctan(\omega \tau_1)$, which is a monotone function in $\omega$ |
|
\item |
|
In general $\delta$ has a maximum of $\omega = \frac 1 {\sqrt{\tau_0 |
|
\tau_1}}$ |
|
\begin{align*} |
|
g &= \frac {\omega (\tau_1 - \tau_0)}{1 + \omega^2 \tau_0 \tau_1} \\ |
|
g' &= \frac {\sqrt{\tau_1 - \tau_0}}{1 + \omega_2 \tau_0 \tau_1} - \frac |
|
{\omega(\tau_1 - \tau_0)}{(1 + \omega^2 \tau_0 \tau_1)^2 } 2 \omega |
|
\tau_0 \tau_1 = 0\\ |
|
0 &= 1 + \omega^2 \tau_0 \tau_1 - 2\omega^2 \tau_0 \tau_1 \\ |
|
\omega &= \frac 1 {\sqrt {\tau_0 \tau_1}} |
|
\end{align*} |
|
\end{enumerate} |
|
|
|
|
|
\pagebreak |
|
\chapter{Continuum Mechanics in 3D} |
|
|
|
\begin{description} |
|
\item[Material Coordinates] (lagrangian) |
|
|
|
\item $A \in \R^3$ the material coordinate/position at $t = 0$ |
|
\item $X(A,t) \in \R^3$, $X(A,0) = A$ |
|
\item $U(A,t) = X(A,t) - A \in \R^3$ displacement |
|
\item $V(A,t) = \p_t U(A,t) $ velocity |
|
\item[Spatial/Eulerian coordinates] $x\in \R^3$ |
|
\item $u(x,t) = U(a(x,t) , t) \Leftrightarrow u(X(A,t) , t) = U(A,t)$ |
|
\item $v(x,t) = V(a(x,t) , t) \Leftrightarrow v(X(A,t) , t) = V(A,t)$ |
|
\end{description} |
|
|
|
Usually we use material coordinates for solids and spatial coordinates for |
|
fluids. |
|
|
|
\subsection{Deformation Gradient} |
|
|
|
Consider $A_0 \in \R^3$ and $A = A_0 + \Delta A \in \R^3$ |
|
|
|
\marginnote{Use the taylor expansion only keeping the first derivative term} |
|
$x = X(A_0 + \Delta A, t) = \underbrace{X(A_0, t)}_{X_0} + \underbrace{\nabla_A |
|
X(A_0,t)}_{F(A_0, t)} \Delta A + o(\Delta A)$ |
|
|
|
\begin{definition}{}{} Deformation gradient |
|
$$\nabla _A X(A_0,t) = F(A_0, t) $$ |
|
|
|
$x - x_0 \sim F(A_0, t) (A - A_0)$ |
|
|
|
eg at $t = 0$ $X(A,t = 0) = A)$, and $F = \nabla_A X = \nabla_A A = \mathbb I$ |
|
|
|
\end{definition} |
|
|
|
\subsubsection{Impenetrability} |
|
|
|
The assumption is you can always go from material coordinates to spatial |
|
coordinates and back: so the $A$ must always have an inverse, so $\det F \ne 0$. |
|
|
|
$$J = \det F = \det \underbrace{\nabla_A X}_{\text{Jacobi matrix}} \ne 0$$ |
|
|
|
This is essentially the impenetrability assumption from 1D continuum mechanics in |
|
three dimensions. |
|
|
|
At $t = 0$, $X(A, t= 0) = A \implies F = \mathbb I \implies J = 1$ |
|
|
|
Since $J(t = 0) = 1$, $J \ne 0$ always implies $J(t) > 0$ |
|
|
|
\week{} |
|
\lecture |
|
|
|
\begin{ex} Uniform dilation |
|
|
|
Motion according to $A \mapsto X(A,t) = \alpha(t) A$ where |
|
$\alpha(t) \ge 0, \alpha(0) = 1$. |
|
|
|
Consider a sphere centred at $\hat 0 = \begin{pmatrix} 0\\0\\0 |
|
\end{pmatrix} $ with radius $r = 1$ at $t = 0$ and at $t = 1$, |
|
$\alpha(1) = 2$, (the circle is being stretched according to its |
|
radius linearly increasing). |
|
|
|
Material Coordinates |
|
|
|
$U(A,t) = X(A, t) - A = (\alpha(t) - 1)A $ |
|
|
|
$V(A,t) = \p_t U = \alpha' (t) A$ |
|
|
|
Spatial coordinates |
|
|
|
Invert $x$: $A = \frac 1 {\alpha (t)} x$ |
|
|
|
$u(x,t) = U(\frac 1 {\alpha(t)}x, t) = \alpha(t) - 1 \frac 1\ |
|
{\alpha(t)} x = (1 - \frac 1 {\alpha(t)}) x$ |
|
|
|
$v(x,t) = v(\frac 1{\alpha(t)} x, t) = \alpha ' (t) \frac 1 |
|
{\alphta(t)} x$ |
|
|
|
$\implies \p_t u \ne v$ like in one dimension. |
|
|
|
$F = \nabla_A X = \nabla _A(\alpha(t) A) = \alpha(t) \nabla_A |
|
\begin{pmatrix} A_1\\A_2\\A_3 \end{pmatrix} = \alpha(t) |
|
\begin{pmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{pmatrix}$ |
|
|
|
To guarantee $J > 0$, $J = \det F = \det (\alpha(t) \mathbb I) = \alpha(t)^3$ |
|
|
|
$\implies$ we need $\alpha(t) > 0$. |
|
\end{ex} |
|
|
|
\begin{ex} |
|
Simple shear |
|
|
|
Motion according to: |
|
|
|
\marginnote{There is a shear along the $y$ axis} |
|
$$ |
|
X(A,t) = \begin{pmatrix} A_1 + \alpha(t) A_2 \\ A_2 \\ A_3 \end{pmatrix} |
|
\quad\text{Where } \alpha(0) = 0 \text{ so } X(A,0) = A |
|
$$ |
|
|
|
\incfig{shear} |
|
Material coordinate: |
|
|
|
$U(A,t) = \begin{pmatrix} \alpha(t) A_2 \\ 0 \\ 0 \end{pmatrix} |
|
\qquad \qquad |
|
V(A,t) = \begin{pmatrix} \alpha'(t) A_2 \\ 0 \\0 \end{pmatrix}$ |
|
|
|
Spatial coordinate: |
|
|
|
Invert $X = \begin{pmatrix} A_1 + \alpha(t) A_2 \\ A_2 \\ A_3 \end{pmatrix} |
|
= \begin{pmatrix} x_1 \\ x_2 \\ x_3 \end{pmatrix} $. |
|
|
|
$A_1 = x_1 - \alpha(t) A_2 = x_1 - \alpha(t) x_2$ |
|
|
|
$A_2 = x_2$ and $A_3 = x_3$ So |
|
|
|
$$A = \begin{pmatrix} x_1 - \alpha(t) x_2 \\ x_2 \\ x_3 \end{pmatrix} $$ |
|
|
|
So |
|
|
|
$$u(x,t) = \begin{pmatrix} \alpha(t) x_2 \\ 0 \\0 \end{pmatrix} $$ |
|
|
|
$$v(x,t) = V\left( \begin{pmatrix} x_1 - \alpha(t) x_2 \\ x_2 \\ x_3 \end{pmatrix} |
|
, t\right) = \begin{pmatrix} \alpha ' (t)x_2 \\ 0\\ 0 \end{pmatrix} $$ |
|
|
|
Deformation gradient |
|
|
|
$F = \nabla_A X = \nabla _a \begin{pmatrix} A_1 + \alpha(t) A_2 \\ A-2 \\ A_3 |
|
\end{pmatrix} = \begin{pmatrix} 1&\alpha(t)&0 \\ 0&1&0 \\ 0&0&1 |
|
\end{pmatrix} $ |
|
|
|
$\det F = 1 > 0$ which does not imply any additional constraints. |
|
|
|
|
|
\end{ex} |
|
|
|
|
|
\subsection{Material Derivative} |
|
|
|
Let $F = (F_1,F_2,F_3)$ respecitvely to $f$ be a quantity in material respective to spatial |
|
coordinates then we ahve |
|
|
|
\marginnote{$X = \begin{pmatrix} X_1\\X_2\\X_3 \end{pmatrix} $} |
|
|
|
$$F(A,t) = f(X(A,t), t) $$ |
|
|
|
\begin{align*} |
|
\frac {\p F }{\p t} &= \frac d {dt} f (X(A,t) , t) \\ |
|
&= \p_t f + \p _{X_1} f |
|
\frac {\p{X_1}}{\p t} + \p _{X_2} f \frac {\p X_2} {\p t} + \p _{X_3} f \frac |
|
{\p X_3}{\p t} \\ |
|
&= \p _t f + \nabla f \cdot \p _ t X \\ |
|
&= \p_t f + \nabla f \p _t U \\ |
|
&= p_t f + \nabla f \cdot v \\ |
|
&= \p_t f + (v \cdot \nabla) f = (\p_t + v \cdot \nabla) f \\ |
|
\end{align*} |
|
|
|
{$(v \cdot \nabla)$ is notation for a differential operator and does |
|
not make sense on its own, rather it is notation for the previous statement} |
|
|
|
\begin{align*} |
|
(v \cdot \nabla) f = \nabla f \cdot v &= |
|
\begin{pmatrix} |
|
f'_{X_1} & f'_{X_2} & f'_{X_3} \\ |
|
f''_{X_1} & f''_{X_2} & f''_{X_3} \\ |
|
f'''_{X_1} & f'''_{X_2} & f'''_{X_3} \\ |
|
\end{pmatrix} \cdot \begin{pmatrix} v_1 \\ v_2 \\ v_3 \end{pmatrix} |
|
\\ |
|
&= \begin{pmatrix} |
|
f'_{X_1} v_1 + f'_{X_2} v_2 + f'_{X_3} v_3 \\ |
|
f''_{X_1}v_1 + f''_{X_2}v_2 + f''_{X_3}v_3 \\ |
|
f'''_{X_1}v_1 + f'''_{X_2}v_2 + f'''_{X_3}v_3 \\ |
|
\end{pmatrix} \\ |
|
&= \begin{pmatrix} \nabla f^1 \cdot v \\ \nabla f^2 \cdot v \\ \nabla |
|
f^3 \cdot v\end{pmatrix} |
|
\end{align*} |
|
|
|
$\implies$ the material derivative: |
|
|
|
$$D_t f = \p_t f + \nabla f \cdot v$$ |
|
|
|
The rate of change change from the perspective of fixed world coordinates, |
|
so how both the thing moves as a result of its own movement and it being |
|
transported by the material. |
|
|
|
\begin{ex} |
|
Uniform dilation |
|
|
|
$u(x,t) = 1 - \frac 1 {\alpha(t)}x$ |
|
|
|
$v(x,t) = \frac {\alpha'} {\alpha} x$ |
|
|
|
|
|
\begin{align*} |
|
D_t u &= \p _t u + \nabla u \cdot v = \frac {\alpha^2} {\alpha^2} x + (1 - |
|
\frac {\alpha (t)}) \mathbb I \cdot \frac {\alpha ' }{\alpha } x \\ |
|
&= \frac {\alpha^2} {\alpha^2}x + (1 - \frac 1 {\alpha(t)}) \cdot |
|
\frac {\alpha '}{\alpha} x \\ |
|
&= \frac {\alpha'}{\alpha} x \\ |
|
\end{align*} |
|
|
|
Which is just the velocity field. |
|
|
|
\end{ex} |
|
|
|
\pagebreak |
|
\appendix |
|
|
|
\chapter{Laplace Tranform Review}\label{laplace} |
|
|
|
\marginnote{Taken from MATH2100 Workbook} |
|
|
|
$$ |
|
F(s)=L(f(t))=\int_{0}^{\infty} e^{-s t} f(t) d t |
|
$$ |
|
Powers |
|
$$ |
|
L\left(t^{n}\right)=\frac{n !}{s^{n+1}} \quad L\left(t^{a}\right)=\frac{\Gamma(a+1)}{s^{a+1}} |
|
$$ |
|
Sine and Cosine |
|
$$ |
|
L(\cos (\alpha t))=\frac{s}{s^{2}+\alpha^{2}} \quad \text { and } \quad L(\sin (\alpha t))=\frac{\alpha}{s^{2}+\alpha^{2}} |
|
$$ |
|
Dirac Delta function |
|
$$ |
|
L(\delta(t-k))=e^{-k s} |
|
$$ |
|
Linearity |
|
$$ |
|
L(a g(t)+b h(t))=a G(s)+b H(s) |
|
$$ |
|
First Shifting Theorem |
|
$$ |
|
L\left(e^{a t} f(t)\right)=F(s-a) |
|
$$ |
|
Second Shifting Theorem |
|
$$ |
|
L(f(t-k) u(t-k))=e^{-k s} F(s) |
|
$$ |
|
Convolution Theorem |
|
$$ |
|
L\left(\int_{0}^{t} f(\tau) g(t-\tau) d \tau\right)=F(s) G(s) |
|
$$ |
|
Differentials |
|
$$ |
|
\begin{gathered} |
|
L\left(f^{(n)}(t)\right)=s^{n} F(s)-s^{n-1} f(0)-s^{n-2} \dot{f}(0)-\ldots-s f^{n-2}(0)-f^{n-1}(0) \\ |
|
L(t f(t))=-\frac{d F(s)}{d s} |
|
\end{gathered} |
|
$$ |
|
Periodic Functions $f(t+p)=f(t)$ |
|
$$ |
|
L(f(t))=\frac{1}{1-e^{-s p}} \int_{0}^{p} e^{-s t} f(t) d t |
|
$$ |
|
|
|
|
|
|
|
|
|
|
|
\end{document}
|
|
|