Commit 7f948a3b authored by conmccoid's avatar conmccoid
Browse files

Initial commit on CRM summer school 2021 notes, responses to DD26 proceedings paper reviewers

parent a5923f5d
\newcommand{\ddt}{\frac{\partial}{\partial t}}
\title{Conference Notes}
\chapter{Stationary Iterative Methods, Felix Kwok}
\section{Multiphysics problems, Monday May 31st}
Single physics problems have specialized discretizations and solvers to take advantage of the specifics of the problems, ie. fast Poisson solvers for heat transfer, ray tracing for wave propagation., upwinding for advection, etc.
Obviously, this doesn't work if you have multiple physics acting at once, ie. a fast Poisson solver isn't great if there's wave propagation as well as heat transfer.
\subsection{Multiphysics Example Problem: Flow in porous media}
Goal: track the evolution of underground fluids, as well as chemical concentration.
Physics involved: fluid dynamics, diffusion, chemical reactions, capillary action.
\subsubsection{Darcy flow model:}
single phase flow, for one fluid only,
\vec{v} = -K(\vec{x}) \Delta p
Conservation of mass is important: the amount accumulated must equal the mass in minus the mass out.
\ddt \left ( \int_V \rho dV \right ) = & -\int_{\partial V} \vec{v} \cdot \vec{n} + \tilde{q}\\
= & -\int_V \nabla \cdot \vec{v} + \tilde{q} \\
\implies \ddt( \phi \rho) = & \nabla \cdot \left ( K(\vec{x}) \rho \nabla p \right) + q
where $\tilde{q}$ is a possible source term and $V$ is the volume in question.
\subsubsection{Two phase flow}
Consider two immiscible phases, ie. oil and water.
Define the saturation of each, $S_w$ and $S_o$, such that $S_w + S_o = 1$.
Each has its own flow,
\vec{v}_w = -K(\vec{x}) \lambda_w(S_w) \nabla p, \quad \vec{v}_o = -K(\vec{x}) \lambda_o(S_w) \nabla p
with each flow depending on their respective saturation.
Let $S=S_w$ and $S_o=1-S$.
\ddt (\phi \rho_w S) = & \nabla \cdot (K(\vec{x}) \rho_w \lambda_w(S) \nabla p) + q_w, \\
\ddt (\phi \rho_o (1-S)) = & \nabla \cdot (K(\vec{x}) \rho_o \lambda_o(S) \nabla p) + q_o.
If everything is incompressible then $\theta$, $\rho_o$ and $\rho_w$ are constant.
Adding the equations leads to
0 = \nabla \cdot (K(\vec{x})(\lambda_w(S) + \lambda_o(S)) \nabla p) + \frac{q_w}{\rho_w} + \frac{q_o}{\rho_o}.
This is an elliptic equation of the form
-\nabla \cdot (K_T(\vec{x},S) \nabla p) = \tilde{q}.
Define the total fluid velocity as $\vec{v}_T = -K(\vec{x}) (\lambda_w (S) + \lambda_o(S)) \nabla p$.
This measures the velocity of the total fluid particle, composed of both oil and water.
This often varies slowly with time.
In the special case in 1D without sources $\vec{v}_T$ is constant.
The velocity of water can now be expressed as a function of total velocity:
\vec{v}_w = f(S) \vec{v}_T := \frac{\lambda_w(S)}{\lambda_w(S) + \lambda_o (S)} \vec{v}_T.
This leads to
\ddt (\phi S) + \nabla \cdot (f(S) \vec{v}_T) = \frac{q_w}{\rho_w}.
In the special case mentioned above this is a hyperbolic conservation law.
The two problems are coupled, since $K_T$ depends on $S$ and $\vec{v}_T$ on $p$.
For fixed $S$ pressure satisfies a linear elliptic equation.
For fixed $\vec{v}_T$ saturation satisfies a nonlinear hyperbolic conservation law.
Therefore, different discretizations and solvers are needed.
\subsubsection{Elliptic problems}
Take our pressure equation with some boundary conditions (Dirichlet, Neumann or Robin).
We need to discretize the PDE with finite differences, finite volume or finite element methods.
We can have pressure at one set of grid points and $K$ between these grid points, staggering them.
For this elliptic problem:
\item the operator is isotropic (no pref. direction)
\item the operator is local (sparse matrix)
\item the operator is unbounded (ill-conditioned matrices)
\item the solution is smooth if $K$ is continuous
\item there is infinite speed of propagation (solution at any given point depends on data everywhere)
Solution methods:
\item Gaussian elimination (direct, good for simple, small problems)
\item Iterative methods (hey, title drop)
\subsubsection{Hyperbolic problems}
Our saturation variable satisfies the transport equation.
When there is no source it is an example of a hyperbolic conservation law.
The simplest example is the advection equation, when $f(S) = S$ and $\phi=1$:
\ddt S + \nabla \cdot (S \vec{v}_T) = 0.
This has solution $S(\vec{x},t) = g(\vec{x} - t \cdot \vec{v}_T)$, where $g(\vec{x}) = S(\vec{x},0)$.
Ergo, the solution values `transport' with speed and direction $\vec{v}_T$ (ray tracing).
Compare and constrast: elliptic problem had infinite speed of propagation, information affected everywhere.
Hyperbolic problem had specific speed of propagation, information affects a single location.
There is preferential directions.
Information does not disperse.
Solution is exactly as smooth as $g$, which may not be.
In applications $f(S)$ is nonlinear, which can cause discontinuities in solutions, ie. shock waves.
Mass conservation is necessary in the solutions.
\subsection{Multiphysics problem 2: cooling by fluid injection}
A heated metal object is cooled by a liquid flowing around it at prescribed velocity and temperature.
As the fluid cools the metal the metal heats the fluid, changing its dynamics.
\item[$\Omega_s$] part of the domain containing solid
\item[$\Omega_f$] part containing fluid
\item[$\Gamma$] interface
\item[$K_s$] heat conductivity, $s$ for solid, $f$ for fluid
\item[$T$] temperature
\item[$C_s$] specific heat capacity
\item[$\rho_s$] density
\item[$\vec{v}$] fluid velocity (unknown)
\ddt (\rho_s C_s T) = & \nabla \cdot (K_s \nabla T) + q \\
\ddt (\rho_f C_f T) + \nabla \cdot (\vec{v} \rho_f C_f T) = & \nabla \cdot (K_f \nabla T)
There are also mass conservation, which leads to momentum equations.
The fluid velocity is determined by Navier-Stokes.
Momentum conservation uses the Cauchy stress tensor.
Change in momentum causes velocity change.
All heat leaving the solid enters the fluid.
Ultimately this gives five coupled equations.
This gives different physics in each region:
\item Navier-Stokes on $\Omega_f$;
\item diffusion on $\Omega_s$;
\item advection-diffusion on $\Omega_f$.
Once NS is discretized we get a saddle-point problem:
\begin{bmatrix} A & B \\ B^\top & 0 \end{bmatrix} \begin{bmatrix} \vec{v} \\ \vec{p} \end{bmatrix} = \begin{bmatrix} \vec{f} \\ \vec{g} \end{bmatrix}
where $\vec{v}$ is velocity and $\vec{p}$ is pressure, which are the components of a solution to NS.
\section{Iterative Methods}
Suppose we look for the solution to a nonlinear system $F(\vec{x})$ where $F: D \subset \bbr^n \to \bbr^n$ is Frechet-differentiable.
There are no direct methods for solving this problem for arbitrary $F$.
We use iterative methods, constructing $\set{\vec{x}^i}_i \in \bbn$ starting from an initial guess $\vec{x}^0$ and trending towards the solution $\vec{x}^*$.
ex. Newton's method, fixed point methods, Jacobi and Gauss-Seidel.
Block stationary methods can be applied to the cooling by fluid injection problem described above.
The momentum and pressure can be split into blocks, and also the heat equation into the two domains.
\section{CG and Preconditioning, Tuesday June 1st}
This slideshow is mostly a review of conjugate gradient methods.
Solving $Ax=f$ is same as minimizing $x^\top A x - 2 x^\top f$.
\chapter{Domain Decomposition, Victorita Dolean}
Original Schwarz method: solve the problem on a subdomain then take the solution at the interface of a second subdomain as an initial condition to solve the problem there.
Repeat, passing interface conditions back and forth until convergence.
Jacobi Schwarz: same as above, but process is done in parallel.
Rather than wait for the first subdomain to finish, the second subdomain is solved using information from the previous step.
Restricted additive Schwarz:
u^{n+1} = \sum_{i=1}^2 E_i (\chi_i u_i^{n+1})
where $\chi_i$ is the respective partition of unity and $E_i$ is the respective extension operator.
Additive Schwarz is the same without the partitions of unity.
\section{Two-Level DDM}
\begin{defn}[Strong scalability]
How the solution time varies with the number of processors for a fixed total problem size.
\begin{defn}[Weak scalability]
How the solution time varies with the number of processors for a fixed problem size per processor.
At each stage of Schwarz, add a coarse space correction by finding the parts of the problem that are slowing the procedure down and solving it separately over the ``coarse'' space.
In essence, the problem is split into two, one solved over a coarse space and another solved over the more standard subdomains.
Look up fictitious space lemma, Nepomnyaschikh 1991.
ffddm is an implementation of parallel solvers in FreeFEM, namely overlapping Schwarz domain decomposition methods.
\href{}{Link to a tutorial.}
\ No newline at end of file
% Master Preamble
\newcommand{\dxdy}[2]{\frac{d #1}{d #2}}
\newcommand{\dxdyk}[3]{\frac{d^{#3} #1}{d {#2}^{#3}}}
\newcommand{\pdxdy}[2]{\frac{\partial #1}{\partial #2}}
\newcommand{\liminfty}[1]{\lim_{#1 \to \infty}}
\newcommand{\limab}[2]{\lim_{#1 \to #2}}
\newcommand{\abs}[1]{\left \vert #1 \right \vert}
\newcommand{\norm}[1]{\left \Vert #1 \right \Vert}
\newcommand{\order}[1]{\mathcal{O} \left ( #1 \right )}
\newcommand{\set}[1]{\left \{ #1 \right \}}
\newcommand{\Set}[2]{\left \{ #1 \ \middle \vert \ #2 \right \}}
\newcommand{\vmat}[1]{\begin{vmatrix} #1 \end{vmatrix}}
\ No newline at end of file
The authors wish to thank the reviewers for their positive reports. The abstract of this paper follows. The guidelines of DDM state this abstract should only appear on the web and not within the manuscript.
"Newton-Raphson preconditioned by Schwarz methods does not have sufficient convergence criteria. We explore an alternating Schwarz method accelerated by Newton-Raphson to find an example where the underlying Schwarz method converges but the Newton-Raphson acceleration fails. Alternating Schwarz is posed as a fixed point iteration to make use of theory for generic root-finding methods. An algorithm is proposed combining several aspects of this theory and others to guarantee convergence."
To the point made by the second reviewer, the situation is indeed significantly more complicated in higher dimensions. Whereas in 1D cycles can only exist by alternating between two sides of the fixed point, in nD these cycles can reposition themselves anywhere on the hypersphere.
\ No newline at end of file
Supports Markdown
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