Skip to content
Snippets Groups Projects
Commit 68a6a189 authored by Justin Pittman's avatar Justin Pittman
Browse files

Moving Numerical Files

parent 1e82fa5b
No related branches found
No related tags found
1 merge request!86To update materials for Numerical Methods weeks in Q1
Pipeline #207357 passed
Stage: test
Hide whitespace changes
Inline Side-by-side
Showing
with 318 additions and 0 deletions
book/numerical_analysis/ExcessMaterial.ipynb 0 → 100644
+ 293
− 0
View file @ 68a6a189
%% Cell type:code id: tags:
``` python
<details><summary><strong>Example 1: Induced Seismicity</strong></summary>
---
Subsurface activities (e.g. extraction of hydrocarbons or cyclic injection) can induce strains in the subsurface. These strains have accumulated energy that can be (partially) released by fault reactivation, possibly leading to seismic events. And the question posed is: how can we mitigate the induced seismicity?
<img src="figs/example1.png" alt= “” width="350px">
```{figure} /figs/example1.png
---
width: 350px
name: Regres
---
Cross-section of fault in subsurface.
```
***Analytical Solution for the Fault Plane (x=0)***
We need to apply Green functions that generate displacement dependencies, i.e. we only care about displacement along the fault!
$$G_x(0,y) = \frac{y-a}{4}\ln\left[(y-a)^2\right] - \frac{y+b}{4}\ln\left[(y+b)^2\right] + \frac{y+a}{4}\ln\left[(y+a)^2\right] - \frac{y-b}{4}\ln\left[(y-b)^2\right]$$
$$G_y(0,y) = \frac{\pi}{4}(|y+a|-|y-a|-|y-b|+|y+b|)$$
And the function describing the shear stress:
$$G_{xy}(0,y) = \frac{1}{2}\ln\left[\frac{(y-a)^2(y+a)^2}{(y-b)^2(y+b)^2}\right]$$
***Numerical Solution (using Finite Volume Methods)***
For the numerical solution, you need to apply some discretization techniques. It could be finite element, finite volume or finite difference (methods that will be covered later on in MUDE). Here we will use the finite volume method to get to a numerical solution for the problem.
$$\pmb{g}(\pmb{\omega}) = \underbrace{V(\pmb{a}^{n+1}-\pmb{a}^n+\Delta t\pmb{q}^{n+1})}_\text{accumulation} + \overbrace{\Delta t \sum_{\delta}|\delta|\pmb{f}_{\delta}^{n+1}}^\text{fluxes} = \pmb{0}$$
where we have $\pmb{q}$ representing the source/sink, $\pmb{f}_{\delta}$ the flux over $\delta$, given by
$$\pmb{f}_{\delta} = - \begin{pmatrix}\color{orange}{\left(\mathbb{C}:\dfrac{\nabla \pmb{u} + \nabla \pmb{u}^T}{2} - \mathbb{B}p - \mathbb{A}T\right)\cdot\pmb{n}} \\ \dfrac{\rho}{\mu}\color{blue}{\mathbb{K}\pmb{n}\cdot(\nabla p} - \rho\color{blue}{\pmb{g}}) - \rho\color{blue}{\mathbb{B}\pmb{n}\cdot\frac{\partial\pmb{u}}{\partial t}}\end{pmatrix}$$
with the orange and blue terms saying respect to momentum and fluid fluxes, respectively. Finally, we describe $\pmb{a} = \{0,\rho\bar{\phi}\}$, where
$$\bar{\phi} = \phi_0 + \frac{(b-\phi_0)(1-b)}{K_d}(p-p_0)$$
***Comparison of Results (Analytical vs Numerical)***
<img src="figs/example1_2.png" alt= “” width="500px">
---
</details>
```
%% Output
Input In [5]
<details><summary><strong>Example 1: Induced Seismicity</strong></summary>
^
SyntaxError: invalid syntax
%% Cell type:code id: tags:
``` python
<details><summary><strong>Example 1: Induced Seismicity</strong></summary>
---
Subsurface activities (e.g. extraction of hydrocarbons or cyclic injection) can induce strains in the subsurface. These strains have accumulated energy that can be (partially) released by fault reactivation, possibly leading to seismic events. And the question posed is: how can we mitigate the induced seismicity?
<img src="figs/example1.png" alt= “” width="350px">
```{figure} /figs/example1.png
---
width: 350px
name: Regres
---
Cross-section of fault in subsurface.
```
***Analytical Solution for the Fault Plane (x=0)***
We need to apply Green functions that generate displacement dependencies, i.e. we only care about displacement along the fault!
$$G_x(0,y) = \frac{y-a}{4}\ln\left[(y-a)^2\right] - \frac{y+b}{4}\ln\left[(y+b)^2\right] + \frac{y+a}{4}\ln\left[(y+a)^2\right] - \frac{y-b}{4}\ln\left[(y-b)^2\right]$$
$$G_y(0,y) = \frac{\pi}{4}(|y+a|-|y-a|-|y-b|+|y+b|)$$
And the function describing the shear stress:
$$G_{xy}(0,y) = \frac{1}{2}\ln\left[\frac{(y-a)^2(y+a)^2}{(y-b)^2(y+b)^2}\right]$$
***Numerical Solution (using Finite Volume Methods)***
For the numerical solution, you need to apply some discretization techniques. It could be finite element, finite volume or finite difference (methods that will be covered later on in MUDE). Here we will use the finite volume method to get to a numerical solution for the problem.
$$\pmb{g}(\pmb{\omega}) = \underbrace{V(\pmb{a}^{n+1}-\pmb{a}^n+\Delta t\pmb{q}^{n+1})}_\text{accumulation} + \overbrace{\Delta t \sum_{\delta}|\delta|\pmb{f}_{\delta}^{n+1}}^\text{fluxes} = \pmb{0}$$
where we have $\pmb{q}$ representing the source/sink, $\pmb{f}_{\delta}$ the flux over $\delta$, given by
$$\pmb{f}_{\delta} = - \begin{pmatrix}\color{orange}{\left(\mathbb{C}:\dfrac{\nabla \pmb{u} + \nabla \pmb{u}^T}{2} - \mathbb{B}p - \mathbb{A}T\right)\cdot\pmb{n}} \\ \dfrac{\rho}{\mu}\color{blue}{\mathbb{K}\pmb{n}\cdot(\nabla p} - \rho\color{blue}{\pmb{g}}) - \rho\color{blue}{\mathbb{B}\pmb{n}\cdot\frac{\partial\pmb{u}}{\partial t}}\end{pmatrix}$$
with the orange and blue terms saying respect to momentum and fluid fluxes, respectively. Finally, we describe $\pmb{a} = \{0,\rho\bar{\phi}\}$, where
$$\bar{\phi} = \phi_0 + \frac{(b-\phi_0)(1-b)}{K_d}(p-p_0)$$
***Comparison of Results (Analytical vs Numerical)***
<img src="figs/example1_2.png" alt= “” width="500px">
---
</details>
```
%% Cell type:markdown id: tags:
<details><summary><strong>Example 2: Floating Wind Turbine</strong></summary>
---
A floating wind turbine is a system subject to several loads, namely aerodynamic (wind), hydrodynamic (wave, currents and tides) or Mooring systems. For a system like this, how can we determine the structural response?
<img src="figs/example_2_1.png" alt= “” width="400px">
For this system we will assume some simplifications that will make our life easier. Starting by its geometry, we are going to consider either the floater and the rotor-nacelle as lumped masses and tower as flexible. Regarding the loads, the wind will me modelled as a fixed rotor and the hydrodynamic effects will be simulated using linear combinations of sinusoidal waves. Finally, we will consider that both the aerodynamic and the hydrodynamic damping are linear.
**Analytical Solution**
Considering the assumptions and simplifications made above, we have a system described by 4 degrees of freedom in hands:
* Floater horizontal motion
* Floater vertical motion
* Floater rotation
* First tower fore-aft modal deflection
and the equation of motion (that you do not need to know how to derive!) is given by the following ODE:
$$m\ddot{\xi}(t) + c\dot{\xi}(t) + k\xi(t) = F(t)$$
in which $\xi$ represents the degree of freedom under analysis and $m$, $c$ and $k$ are coefficients to describe the mass, damping and stiffness (respectively) of the degree of freedom being studied. Finally, $F$ represents the external loads.
To solve this, we usually go to the frequency domain, where we suggest:
$$\xi(t) = \mathbb{R}\{\hat{\xi}(\omega)e^{i\omega t}\}$$
with this, we can write:
$$[-\omega^2m + i\omega c + k]\hat{\xi}(\omega) = \hat{F}(\omega) \Leftrightarrow \hat{\xi}(\omega) = \frac{\hat{F}(\omega)}{[-\omega^2m + i\omega c + k]} \equiv H(\omega)\hat{F}(\omega)$$
For a frequency $\omega$, the motion of the wind turbine is equal to a linear transfer function $H(\omega)$ times the loading force. This transfer function will only depend on the frequency, mass, damping, and stiffness of the system.
**Numerical Solution (using the Finite Difference Method)**
Once again, and as shown in the previous example, an approximation method should be used to apply and find a numerical solution. In this case we will use the Finite Difference Method. The discrete equations for air-water urge as:
$$\underbrace{\frac{1}{J}\frac{\pmb{U^*}-\pmb{U^n}}{\Delta t}}_\text{Discrete time derivative at intermediate step} = \overbrace{\pmb{P}(p^n,\phi^n)}^\text{Pressure gradient} + \underbrace{\frac{1}{2}(\pmb{F}(\pmb{U^*}, \pmb{u^*}, \phi^{n+1})}_\text{Momentum flux at intermediate step} + \overbrace{\pmb{F}(\pmb{U^n},\pmb{u^n},\phi^n))}^\text{Momentum flux at step $n$}$$
**Comparison of Results (Analytical vs Numerical)**
For a certain set of input conditions looking like
<img src="figs/example_2_2.png" alt= “” width="500px">
The final results obtained via analytical and numerical methods plotted look like
<img src="figs/example_2_3.png" alt= “” width="450px">
---
</details>
%% Cell type:markdown id: tags:
<details><summary><strong>Example 3: Dynamic Behaviour of a Tower</strong></summary>
---
Let us imagine a different example, where now the main goal is to predict the dynamic behavior of the tower, specially for the first five natural frequencies.
**Analytical Solution**
The simplest way to think about this problem and formulate it analytically would be to think of it as a cantilever beam, fixed/clamped on one side and free to move at the top. This behavior can be described according to the following equation:
$$f_n = 2\pi C_n\sqrt{\frac{EI}{AL\rho}} = \pi^2(2n-1)\sqrt{\frac{EI}{AL\rho}}$$
where we need to use as input the following parameters (either by getting them or estimating them):
* I - moment of inertia of the cross-section
* A - area of the cross-section
* $\rho$ - density of the material
* E - Young's Modulus of the material
* L - length of the beam, or height of the structure
**Numerical Solution (using the Finite Element Method)**
This time we use the Finite Element Method to simulate the approximately 17 000 degrees of freedom existing in our system. Doing so, we have made some relevant, but reasonable assumptions:
* Elastic and ortothropic behavior of the material, stone + mortar
* Flexible connection with the neighboring church (we are consider a really specific example where the tower under study is connected to a church)
* Rigid connection with the soil
Solving everything considering these assumptions, we get to the following results (where the theoretical natural frequency expected is found on top of the label of each image and the one obtained by our numerical method is showed under everything):
<img src="figs/example_3.png" alt= “” width="500px">
**Comparison of Results (Analytical vs Numerical)**
The results show an amazing accuracy, even though some approximations were taken and some differences were found for larger frequencies under this prediction, as we can see in the figure above.
---
</details>
%% Cell type:markdown id: tags:
<details><summary><strong>Example 4: Nonlinear Response of Composites</strong></summary>
---
For this example, we want to describe the nonlinear response of a layer in a fiber reinforced laminate. The response can be measured for basic load cases (e.g. uniaxial stress).
**Analytical Solution**
To describe this problem analytically, we explore the orthotropic plasticity of our system. Therefore, the stress states can be described by an evolving yield function given by
$$\pmb{\sigma^T}\pmb{P}\pmb{\sigma}+\pmb{\sigma^T}\pmb{p}-1\leq0$$
where
$$\pmb{P} = \begin{bmatrix}0 & 0 & 0 & 0 & 0 & 0 \\ 0 & \dfrac{1}{2}\alpha_1+2\alpha_{32} & -\dfrac{1}{2}\alpha_1+2\alpha_{32} & 0 & 0 & 0 \\ 0 & -\dfrac{1}{2}\alpha_1+2\alpha_{32} & \dfrac{1}{2}\alpha_1+2\alpha_{32} & 0 & 0 & 0 \\ 0 & 0 & 0 & \dfrac{1}{2}\alpha_1 & 0 & 0 \\ 0 & 0 & 0 & 0 & \dfrac{1}{2}\alpha_2 & 0 \\ 0 & 0 & 0 & 0 & 0 & \dfrac{1}{2}\alpha_2\end{bmatrix}$$
and
$$\pmb{p} = \begin{bmatrix}0 & \alpha_3 & \alpha_3 & 0 & 0 & 0\end{bmatrix}$$
The evolution of all the $\alpha$ parameters is based on experimental stress/strain curves.
**Numerical Solution (using the Finite Element Method)**
For the numerical solution, we will use a micromechanical finite element model to simulate our composite response. But what does this mean? That we will consider a geometry with random fiber distribution, periodic boundary conditions, and simple plasticity model for the polymer matrix. Therefore, our system will be described by
$$\begin{cases}\pmb{\nabla\cdot\sigma} = 0 \\ \epsilon=\pmb{\nabla^su}\end{cases}$$
where $\pmb{\sigma} = \pmb{D}(\epsilon-\epsilon^p)$
**Comparison of Results (Analytical vs Numerical)**
The results of both the numerical and the analytical analysis are shown below, where we can clearly see that these approaches have a really good agreement on some stress states (left). However, we can see that the analytical model is not able to captue all the interactions that can occur (right), as we find different possible modes with the numerical analysis that are all superpositioned for the analytical analysis done.
<img src="figs/example_4.png" alt= “” width="800px">
---
</details>
book/numerical_analysis/README.md 0 → 100644
+ 25
− 0
View file @ 68a6a189
##Note to self for integration of topics.
###Programming Assign
1.5
1.6
Dhruv: For the coding assignment for week 1.6, I believe the students should be brought up to speed with calculating the mean, root mean square (e2), obtaining the absolute value, maximum and minimum value inside an array (e_infinity) and creating a log log plot with a detailed legend. I believe that should be enough? They should be able to use these to write the small portion of code for calculating these errors and filling in the table?
###Wed Tutorial
1.5
1.6
Last year's archive is FD: 1D & 2D Diffusion using FD in python. Everything is provided and they can run and see the results. It's a read & understand assignment.
###Project
1.5
1.6
____________________________________________________________
book/numerical_analysis/figs/01s.png 0 → 100644
+ 0
− 0
View file @ 68a6a189
book/numerical_analysis/figs/01s.png

40.8 KiB

book/numerical_analysis/figs/comparison1.png 0 → 100644
+ 0
− 0
View file @ 68a6a189
book/numerical_analysis/figs/comparison1.png

44.3 KiB

book/numerical_analysis/figs/comparison2.png 0 → 100644
+ 0
− 0
View file @ 68a6a189
book/numerical_analysis/figs/comparison2.png

55.9 KiB

book/numerical_analysis/figs/comparison3.png 0 → 100644
+ 0
− 0
View file @ 68a6a189
book/numerical_analysis/figs/comparison3.png

61.8 KiB

book/numerical_analysis/figs/example1.png 0 → 100644
+ 0
− 0
View file @ 68a6a189
book/numerical_analysis/figs/example1.png

236 KiB

book/numerical_analysis/figs/example1_2.png 0 → 100644
+ 0
− 0
View file @ 68a6a189
book/numerical_analysis/figs/example1_2.png

132 KiB

book/numerical_analysis/figs/example_2_1.png 0 → 100644
+ 0
− 0
View file @ 68a6a189
book/numerical_analysis/figs/example_2_1.png

90.3 KiB

book/numerical_analysis/figs/example_2_2.png 0 → 100644
+ 0
− 0
View file @ 68a6a189
book/numerical_analysis/figs/example_2_2.png

171 KiB

book/numerical_analysis/figs/example_2_3.png 0 → 100644
+ 0
− 0
View file @ 68a6a189
book/numerical_analysis/figs/example_2_3.png

301 KiB

book/numerical_analysis/figs/example_3.png 0 → 100644
+ 0
− 0
View file @ 68a6a189
book/numerical_analysis/figs/example_3.png

367 KiB

0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment