multiphysicsProblemSet2aProblem1.tex

%
% Copyright © 2014 Peeter Joot.  All Rights Reserved.
% Licenced as described in the file LICENSE under the root directory of this GIT repository.
%
\index{conjugate gradient}
\index{Gershgorin circle theorem}
\makeproblem{Conjugate gradient and Gershgorin circles.}{multiphysics:problemSet2a:1}{
\makesubproblem{}{multiphysics:problemSet2a:1a}
Consider an arbitrary circuit made by positive resistors and independent DC current sources. Prove, mathematically, that its nodal matrix is always symmetric and positive semi-definite.
\makesubproblem{}{multiphysics:problemSet2a:1b}
Consider an electrical network made by:
\begin{itemize}
\item \( N \times N \) square grid of resistors of value \( R \), where \( N \) is the number of resistors per edge. The grid nodes are numbered from \( 1 \) to \( (N +1)^2 \) .  Node 1 is a corner node
\item resistor \( R_{\textrm{g}} \) between each node of the grid and the reference node.
\item a non-ideal voltage source connected between node 1 and ground.  The voltage source has value \( V_{\textrm{s}} \) and internal (series) resistance \( R_{\textrm{s}} \)
\item three current sources between three randomly-selected nodes of the grid and ground. The source current must flow from the grid to the reference node as shown in \cref{fig:problemSet2a:problemSet2aFig1} (and not vice versa!)  Choose the current source values randomly between 10 mA and 100 mA.
\imageFigure{../figures/ece1254-multiphysics/problemSet2aFig1}{Subcircuit.}{fig:problemSet2a:problemSet2aFig1}{0.2}
\end{itemize}
Write a Matlab routine that generates a SPICE compatible netlist for this system (you can reuse the code you developed for the first problem set). Generate the modified nodal analysis equations
\begin{equation}\label{eqn:multiphysicsProblemSet2a:20}
\BG \Bx = \Bb,
\end{equation}
for a grid with \( N = 40, R = 0.1 \Omega, R_{\textrm{g}} = 1M \Omega, V_{\textrm{s}} = 2V, R_{\textrm{s}} = 0.1 \Omega\).
Then, write in Matlab your own routine for the conjugate gradient
method. Give to the user the possibility of specifying a preconditioning
matrix P. The routine shall stop iterations when the residual norm
satisfies
\begin{equation}\label{eqn:multiphysicsProblemSet2a:40}
\frac{\Norm{ \BG \Bx - \Bb }_2}{
\Norm{\Bb}_2} < \epsilon,
\end{equation}
where \( \epsilon \) is a threshold specified by the user. Use the conjugate gradient method to solve modified nodal analysis \cref{eqn:multiphysicsProblemSet2a:20} of the grid.

Does the conjugate gradient method converge?
\makesubproblem{}{multiphysics:problemSet2a:1c}
Discuss if the conjugate gradient method can be applied to the modified nodal analysis equations of the grid. Support your answer with some numerical results.
\makesubproblem{}{multiphysics:problemSet2a:1d}
Suggest a transformation of the grid that will lead to a circuit equivalent to the original one, but for which conjugate gradient can be used. We call this new circuit ``2D grid''. You will have to use it for all the questions that follow.
\makesubproblem{}{multiphysics:problemSet2a:1e}
Solve the ``2D grid'' circuit using three methods:
\begin{itemize}
\item your own LU decomposition
\item the conjugate gradient method with \( \epsilon = 10^{-3} \)
\item the conjugate gradient method with \( \epsilon = 10^{-3} \) and a tri-diagonal preconditioner \( \BP \)
\end{itemize}
for increasing N. Use \( R = 0.1 \Omega \). Plot the CPU time taken by the three
methods vs N, and the number of iterations taken by the two iterative
methods. Explore a range of N compatible with your PC speed and
memory, but make sure to reach fairly large values of N.
\makesubproblem{}{multiphysics:problemSet2a:1f}
Does preconditioning reduce the number of iterations required by conjugate gradient? Does preconditioning reduce CPU time?
\makesubproblem{}{multiphysics:problemSet2a:1g}
Try to make your code for the conjugate gradient method as efficient as possible. Show the improvements that you have obtained.
\makesubproblem{}{multiphysics:problemSet2a:1h}
Generate nodal analysis \cref{eqn:multiphysicsProblemSet2a:20} for ``2D grid'' with
\( N = 20\). Plot the eigenvalues of \( \BG \) and the Gershgorin circles in three
cases: \( R = 0.1 \Omega, R = 1 \Omega, R = 10 \Omega \). Verify Gershgorin Circle theorem
in the three cases.
\makesubproblem{}{multiphysics:problemSet2a:1i}
 Repeat the previous point for the preconditioned nodal matrix, and compare the circles obtained in the two cases.
\makesubproblem{}{multiphysics:problemSet2a:1j}
Let \( N = 20\), and solve ``2D grid'' with conjugate gradient with
and without preconditioning (use the tri-diagonal preconditioner). Plot
the normalized norm of the residual
\begin{equation}\label{eqn:multiphysicsProblemSet2a:60}
\frac{\Norm{ \BG \Bx - \Bb }_2}{
\Norm{\Bb}_2},
\end{equation}
versus the iteration counter for three cases: \( R = 0.1\Omega, R = 1\Omega, R =
10\Omega\). Discuss your findings in the light of Gershgorin circles.
} % makeproblem
%
\makeanswer{multiphysics:problemSet2a:1}{
\withproblemsetsParagraph{
\makeSubAnswer{}{multiphysics:problemSet2a:1a}
%
\paragraph{Symmetry}
%
A circuit with only resistors and constant current sources has the MNA form
%
\begin{equation}\label{eqn:multiphysicsProblemSet2a:80}
\BG \BV = \BI_{\textrm{S}}.
\end{equation}
%
The current sources, appearing only in the vector \( \BI_{\textrm{S}} \), do not effect the symmetry of the \textAndIndex{nodal matrix} \( \BG \).  To show that \( \BG \)  is symmetric, recall that
%for resistor node specification corresponding to \cref{fig:lecture2:lecture2Fig2}
this was the sum of stamp matrices of the form
%\imageFigure{../figures/ece1254-multiphysics/lecture2Fig2}{Resistor node specification.}{fig:lecture2:lecture2Fig2}{0.3}
%
\begin{equation}\label{eqn:multiphysicsProblemSet2a:100}
\kbordermatrix{
    & n_1      & n_2 \\
n_1 & \inv{R}  & -\inv{R} \\
n_2 & -\inv{R} & \inv{R}
}.
\end{equation}
%
To apply this to an actual system, more precision is required, since subsets of the stamps may not apply to \( \BG \) for resistors connected to the reference node.  Suppose that all \( N \) resistors can be enumerated, each connecting all pairs of nodes \( i, j \), where \( i < j \) with a resistance \( R_{ij} \).  The \( r,s\) element of nodal matrix is then the sum over all the stamps, as follows
%
\begin{equation}\label{eqn:multiphysicsProblemSet2a:120}
\begin{aligned}
\lr{\BG}_{r,s}
&=
\sum_{0 < i < j \le N} \inv{R_{ij}} \lr{
    \delta_{ri} \delta_{si}
   +\delta_{rj} \delta_{sj}
   -\delta_{ri} \delta_{sj}
   -\delta_{rj} \delta_{si}
}
+
\sum_{i = 0, 0 < j \le N} \inv{R_{0j}}
   \delta_{rj} \delta_{sj}
\\ &= \sum_{0 < i < j \le N} \inv{R_{ij}} \lr{
   \delta_{ri} \lr{ \delta_{si} - \delta_{sj} }
   +\delta_{rj} \lr{ \delta_{sj} - \delta_{si} }
}
+
\sum_{i = 0, 0 < j \le N} \inv{R_{0j}}
   \delta_{rj} \delta_{sj}
\\ &= \sum_{0 < i < j \le N} \inv{R_{ij}}
   \lr{ \delta_{ri} - \delta_{rj}} \lr{ \delta_{si} - \delta_{sj} }
+
\sum_{i = 0, 0 < j \le N} \inv{R_{0j}}
   \delta_{rj} \delta_{sj}.
\end{aligned}
\end{equation}
%
Here the sum is over all the sets of nodes \( i, j \), with an allowance for multiplicities in \( R_{ij} \) between any pair of nodes \( i, j \) (i.e. a parallel resistors between any pair of nodes).
%
Transposition requires swapping the indexes
%
\begin{equation}\label{eqn:multiphysicsProblemSet2a:140}
\lr{\BG^\T}_{r,s}
= \sum_{0 < i < j \le N} \inv{R_{ij}}
   \lr{ \delta_{si} - \delta_{sj}} \lr{ \delta_{ri} - \delta_{rj} }
+
\sum_{i = 0, 0 < j \le N} \inv{R_{0j}}
   \delta_{sj} \delta_{rj},
\end{equation}
%
which is clearly equivalent to \cref{eqn:multiphysicsProblemSet2a:120}, proving that \( \BG \) is symmetric.
%
\paragraph{Positive semi-definite}
\index{positive semi-definite}
%One way to prove that the nodal matrix is positive semi-definite, would be to show that all the eigenvalues \( \lambda_i \ge 0 \).
%Although the eigenvalues are all real (symmetric matrix), non-negative eigenvalues cannot be ruled out without some thought.  With real eigenvalues Gershgorin's theorem takes on a simple form, but using that to demonstrate the lower bound it provides is non-negative turns out to be non-trivial (and perhaps not possible).
%
%It turns out to be more profitable to ignore the nature of the eigenvalues and
The positive semi-definite nature of the nodal matrix can be demonstrated directly from the defining property
%
\begin{equation}\label{eqn:multiphysicsProblemSet2a:260}
\Bx^\T \BG \Bx \ge 0, \quad \forall \Bx.
\end{equation}
%
The delta filtering and the nice factorization of the products in the \( \BG \) summation facilitate this
%
%It's possible that Gershgorin's could be used to provetheorem to prove that the eigenvalues are all greater than or equal to zero.
%%%
%%%For row \( r \) of this matrix, that theorem states
%%%
%%%\begin{equation}\label{eqn:multiphysicsProblemSet2a:160}
%%%\Abs{\lambda - G_{rr} }
%%%\le \sum_{s \ne r} \Abs{ G_{rs} }.
%%%\end{equation}
%%%
%%%Because the eigenvalues are real, this can be written as
%%%
%%%\begin{equation}\label{eqn:multiphysicsProblemSet2a:180}
%%%G_{rr} - \sum_{s \ne r} \Abs{ G_{rs} } \le \lambda \le G_{rr} + \sum_{s \ne r} \Abs{ G_{rs} },
%%%\end{equation}
%%%
%%%so if we want to show that \( \lambda \ge 0 \), this is equivalent to showing
%%%
%%%\begin{equation}\label{eqn:multiphysicsProblemSet2a:200}
%%%G_{rr} \ge \sum_{s \ne r} \Abs{ G_{rs} },
%%%\end{equation}
%%%
%%%or
%%%\begin{equation}\label{eqn:multiphysicsProblemSet2a:220}
%%%\sum_{0 < i < j \le N} \inv{R_{ij}}
%%%\biglr{
%%%   \lr{ \delta_{ri} - \delta_{rj}}^2
%%%   -
%%%   \sum_{r \ne s} \Abs{\lr{ \delta_{ri} - \delta_{rj}} \lr{ \delta_{si} - \delta_{sj} } }
%%%   } \ge 0.
%%%\end{equation}
%%%
%%%This will clearly be the case if the inequality holds for each \( i, j \), or
%%%
%%%\begin{equation}\label{eqn:multiphysicsProblemSet2a:240}
%%%   \lr{ \delta_{ri} - \delta_{rj}}^2
%%%   \ge
%%%   \sum_{r \ne s} \Abs{\lr{ \delta_{ri} - \delta_{rj}} \lr{ \delta_{si} - \delta_{sj} }} .
%%%\end{equation}
%
\begin{equation}\label{eqn:multiphysicsProblemSet2a:280}
\begin{aligned}
\Bx^\T \BG \Bx
&=
\Bx^\T \sum_s \sum_{0 < i < j \le N} \inv{R_{ij}}
   \lr{ \delta_{ri} - \delta_{rj}} \lr{ \delta_{si} - \delta_{sj} } x_s
+
\Bx^\T \sum_s \sum_{i = 0, 0 < j \le N} \inv{R_{0j}}
   \delta_{rj} \delta_{sj} x_s
\\ &=
\sum_r x_r \sum_{0 < i < j \le N} \inv{R_{ij}}
   \lr{ \delta_{ri} - \delta_{rj}} \lr{ x_i - x_j }
+
\sum_r x_r \sum_{i = 0, 0 < j \le N} \inv{R_{0j}}
   \delta_{rj} x_j
\\ &=
\sum_{0 < i < j \le N} \inv{R_{ij}}
   \lr{ x_i - x_j }^2
+
\sum_{i = 0, 0 < j \le N} \inv{R_{0j}}
   x_j^2
\ge 0.
\end{aligned}
\end{equation}
%
With all the sums involving only resistors (positive) and squares of real numbers (voltages), the nodal matrix has been shown to be semi-positive definite as desired.
%
%FIXME: remove this on integration into the notes compilation.
%\input{../ece1254/conjugateGradientAlgorithm.tex}
%
\makeSubAnswer{}{multiphysics:problemSet2a:1b}
%
\paragraph{Generation and solution of the netlist equations}
%
Code for the netlist generation, the nodal equation generation were implemented respectively in
%
\begin{itemize}
\item \nbcite{ps2a:genResistorGridNl.m}{ps2a/genResistorGridNl.m}
\item \nbcite{ps2a:NodalAnalysis.m}{ps2a/NodalAnalysis.m}
\end{itemize}
%
These were called from \matlabFuncPath{usenetlistProblemBD(1, 1e-10)}{ps2a:usenetlistProblemBD.m}.
%\nbcite{ps2a:usenetlistProblemBD.m}{ps2a/usenetlistProblemBD.m}
%
The CG iteration for this matrix did not converge.
%It is possible that better initial iteration vectors \( \Bx^{(0)} \) could have allowed for convergence.
% FIXME: ... that would be testable by seeding the iteration with something near the G\b value.
%
\makeSubAnswer{}{multiphysics:problemSet2a:1c}
%
The CG iterations for the circuit above did not converge, oscillating, and eventually appearing to slowly diverge.
A sampling of that lack of convergence is illustrated in \cref{tab:multiphysicsProblemSet2a:300}, data collected using the CG implementation \nbcite{ps2a:cgQuarteroniPrecond.m}{ps2a/cgQuarteroniPrecond.m} called with \( \epsilon = 10^{-10} \) as noted above.
%
\captionedTable{Sampling of relative error at various iteration counts.}{tab:multiphysicsProblemSet2a:300}{
\begin{tabular}{|l|l|}
\hline
Iteration & Relative error \\ \hline
0 & \( 5.2015 \times 10^3 \) \\ \hline
200 & \( 0.0527 \times 10^3 \) \\ \hline
400 & \( 0.3300 \times 10^3 \) \\ \hline
600 & \( 0.7012 \times 10^3 \) \\ \hline
800 & \( 1.6846 \times 10^3 \) \\ \hline
1000 & \( 3.4659 \times 10^3 \) \\ \hline
1200 & \( 5.2740 \times 10^3 \) \\ \hline
1400 & \( 6.6699 \times 10^3 \) \\ \hline
1600 & \( 7.8395 \times 10^3 \) \\ \hline
\end{tabular}
}
%
This lack of convergence is not surprising since the circuit has voltage sources, which results in a non-symmetric nodal matrix.  \citep{shewchuk1994introduction} suggests resetting the residual periodically to avoid cumulative numerical error, and it seems possible that the divergence observed here is also exasperated by such cumulative numerical error.
%
\makeSubAnswer{}{multiphysics:problemSet2a:1d}
%
Application of Norton's theorem yields an equivalent current source value of
%
\begin{equation}\label{eqn:multiphysicsProblemSet2a:21}
I_{\textrm{s}}
=
\frac{V_{\textrm{s}}}{R_{\textrm{S}}},
\end{equation}
%
placed in parallel with a resistance of \( R_{\textrm{S}} \).  An exploration of Norton equivalents can be found in \cref{chap:simpleNortonEquivalents}.
%
A test of this transformation can be found in \nbcite{ps2a:compareBD.m}{ps2a/compareBD.m} which compares the voltages at each of the interior node locations with and without this source transformation.  The norm of the difference in voltage vectors produced by these respective netlist solutions (\( 4.6 \times 10^{-12} \)) is effectively zero.
%
\makeSubAnswer{}{multiphysics:problemSet2a:1e}
%
The time required to generate the nodal equations on my machine was reasonable for \( N \le 196 \).  For \( N = 256 \) Matlab ran out of memory generating the nodal equations.
The nodal equation variables \( G, b \) for \( N \in \{ 8, 16, 32, 64, 128, 196 \} \), were \matlabFunc{save()}d to .mat files for subsequent use 
in \nbcite{ps2a:nodalEquationsPartEH.m}{ps2a/nodalEquationsPartEH.m}
using a call to
\matlabFuncPath{nodalEquationsPartEH(1)}{ps2a:nodalEquationsPartEH.m}
.  Timing data was collected and plotted in
\nbcite{ps2a:partEtimings.m}{ps2a/partEtimings.m}.
%
LU timing data was collected for only \( N \le 32 \), and is plotted in \cref{fig:luTimings:luTimingsFig1}.
%
\mathImageFigure{../figures/ece1254-multiphysics/luTimingsFig1.pdf}{LU solution times.}{fig:luTimings:luTimingsFig1}{0.3}{ps2a:partEtimings.m}
%
While I was able to generate the nodal equations for \( N = 196 \), attempting to solve that set of nodal equations with the methods specified made my system unresponsive (with Matlab eventually using 9Gb of memory before I killed it).  Timing data was collected for \( N \le 128 \), and is plotted with and without tridiagonal preconditioning in \cref{fig:timingsCG:timingsCGFig2}.
%
\mathImageFigure{../figures/ece1254-multiphysics/timingsCGFig2.pdf}{CG solution times.}{fig:timingsCG:timingsCGFig2}{0.3}{ps2a:partEtimings.m}
%
The iteration counts with and without tridiagonal preconditioning are plotted in \cref{fig:iterationsCG:iterationsCGFig3}.
%
\mathImageFigure{../figures/ece1254-multiphysics/iterationsCGFig3.pdf}{CG iteration counts for convergence.}{fig:iterationsCG:iterationsCGFig3}{0.3}{ps2a:partEtimings.m}
%
\makeSubAnswer{}{multiphysics:problemSet2a:1f}
%
Preconditioning is seen to reduce the number of iterations required to meet the desired relative error.  However, this increased the time required to solve the system in the initial implementation of the preconditioned algorithm.
%
\makeSubAnswer{}{multiphysics:problemSet2a:1g}
%
Preconditioning is seen to have a runtime cost.  Without it there are also additional optimization available, as detailed above under ``CG without preconditioning''.  The runtime benefits of those optimizations turn out to be fairly small (at least compared to the additional preconditioning cost) and are plotted in \cref{fig:timingsCGwithOpt:timingsCGwithOptFig4}.
%
\mathImageFigure{../figures/ece1254-multiphysics/timingsCGwithOptFig4.pdf}{Benefits of CG optimized for no preconditioning.}{fig:timingsCGwithOpt:timingsCGwithOptFig4}{0.3}{ps2a:partEtimings.m}
%
As implemented, all the repeated calculations are already avoided, with temporarily used to hold any results used multiple times.  The only clear opportunity for gain is by precalculating LU factors for the preconditioning matrix instead of using \matlabText{G \(\backslash\) r} in each iteration.
%
Using \matlabFunc{[ L, U, P ] = lu(A)} command, we can replace \matlabText{A \(\backslash\) b } with \matlabText{U \(\backslash\)(L\(\backslash\)(P * b))}.
%
With this optimization the time required for CG with the preconditioner application begins to look closer to linear, but this optimization is not enough to drop the time required down to that of CG without the preconditioner.  This is plotted in \cref{fig:timingsCGwithOptLU:timingsCGwithOptLUFig5}.
%
\mathImageFigure{../figures/ece1254-multiphysics/timingsCGwithOptLUFig5.pdf}{LU optimization for the CG preconditioner.}{fig:timingsCGwithOptLU:timingsCGwithOptLUFig5}{0.3}{ps2a:partEtimings.m}
%
\makeSubAnswer{}{multiphysics:problemSet2a:1h}
%
See 
\matlabFuncPath{nodalEquationsPartEH(0)}{ps2a:nodalEquationsPartEH.m}
for the matlab code that generated the nodal equations.
The Gershgorin circles and eigenvalues were plotted with \matlabFuncPath{partH(1, 0)}{ps2a:partH.m}.
Plots for \( R = 0.1, 1, 10 \) can be found respectively in \cref{fig:gershgorinPlotR0_1:gershgorinPlotR0_1Fig6}, \cref{fig:gershgorinPlotR1:gershgorinPlotR1Fig7}, and \cref{fig:gershgorinPlotR10:gershgorinPlotR10Fig8}.
%
\mathImageFigure{../figures/ece1254-multiphysics/gershgorinPlotR0_1Fig6.pdf}{Gershgorin circles and eigenvalues for ``2D grid'' \( N = 20, R = 0.1 \Omega\).}{fig:gershgorinPlotR0_1:gershgorinPlotR0_1Fig6}{0.3}{ps2a:partH.m}
\mathImageFigure{../figures/ece1254-multiphysics/gershgorinPlotR1Fig7.pdf}{Gershgorin circles and eigenvalues for ``2D grid'' \( N = 20, R = 1 \Omega\).}{fig:gershgorinPlotR1:gershgorinPlotR1Fig7}{0.3}{ps2a:partH.m}
\mathImageFigure{../figures/ece1254-multiphysics/gershgorinPlotR10Fig8.pdf}{Gershgorin circles and eigenvalues for ``2D grid'' \( N = 20, R = 10 \Omega\).}{fig:gershgorinPlotR10:gershgorinPlotR10Fig8}{0.3}{ps2a:partH.m}
%
No eigenvalues are found outside of some surrounding Gershgorin circle.
%
\makeSubAnswer{}{multiphysics:problemSet2a:1i}
%
The Gershgorin circles and eigenvalues with tridiagonally preconditioned \( G \) were plotted with \matlabFuncPath{partH(1, 1)}{ps2a:partH.m}.
Plots for \( R = 0.1, 1, 10 \) can be found respectively in
\cref{fig:gershgorinPlotPreconditionedR0_1:gershgorinPlotPreconditionedR0_1Fig9},
\cref{fig:gershgorinPlotPreconditionedR1:gershgorinPlotPreconditionedR1Fig10},
and
\cref{fig:gershgorinPlotPreconditionedR10:gershgorinPlotPreconditionedR10Fig11}
.
%
\mathImageFigure{../figures/ece1254-multiphysics/gershgorinPlotPreconditionedR0_1Fig9.pdf}{Gershgorin circles and eigenvalues for tridiagonally preconditioned ``2D grid'' \( N = 20, R = 10 \Omega\).}{fig:gershgorinPlotPreconditionedR0_1:gershgorinPlotPreconditionedR0_1Fig9}{0.3}{ps2a:partH.m}
\mathImageFigure{../figures/ece1254-multiphysics/gershgorinPlotPreconditionedR1Fig10.pdf}{Gershgorin circles and eigenvalues for tridiagonally preconditioned ``2D grid'' \( N = 20, R = 10 \Omega\).}{fig:gershgorinPlotPreconditionedR1:gershgorinPlotPreconditionedR1Fig10}{0.3}{ps2a:partH.m}
\mathImageFigure{../figures/ece1254-multiphysics/gershgorinPlotPreconditionedR10Fig11.pdf}{Gershgorin circles and eigenvalues for tridiagonally preconditioned ``2D grid'' \( N = 20, R = 10 \Omega\).}{fig:gershgorinPlotPreconditionedR10:gershgorinPlotPreconditionedR10Fig11}{0.3}{ps2a:partH.m}
%
The imaginary component of the eigenvalues observed in these plots show that the transformation \( G \rightarrow P^{-1/2} G P^{-1/2} \) was not a symmetric-matrix preserving transformation.  It was effective in rescaling the eigenvalues of the problem, especially for \( R = 0.1 \), which helps speed convergence despite not being symmetric.
%non-symmetric (voltage source) nodal matrix system was not solvable using CG, I wonder if it is it possible that the destruction of the symmetric nature of the matrix made the iteration converge more slowly?
Note that the tridiagonal preconditioner used in this case was nothing more fancy than the original matrix \( G \) with all non-tridiagonal terms zeroed.
%
\makeSubAnswer{}{multiphysics:problemSet2a:1j}
%
Residual vs iteration plots for \( N = 20 \) were generated with \matlabFuncPath{partH(0, 1)}{ps2a:partH.m}.  For \( R = 0.1 \Omega, 1 \Omega \), and \( 10 \Omega \) respectively these are, respectively \cref{fig:residualsByIterationR0_1:residualsByIterationR0_1Fig6}, \cref{fig:residualsByIterationR1_0:residualsByIterationR1_0Fig7}, \cref{fig:residualsByIterationR10_0:residualsByIterationR10_0Fig8}.
%
\mathImageFigure{../figures/ece1254-multiphysics/residualsByIterationR0_1Fig6.pdf}{Residual vs iteration \( R = 0.1 \Omega \).}{fig:residualsByIterationR0_1:residualsByIterationR0_1Fig6}{0.3}{ps2a:partH.m}
\mathImageFigure{../figures/ece1254-multiphysics/residualsByIterationR1_0Fig7.pdf}{Residual vs iteration \( R = 1 \Omega \).}{fig:residualsByIterationR1_0:residualsByIterationR1_0Fig7}{0.3}{ps2a:partH.m}
\mathImageFigure{../figures/ece1254-multiphysics/residualsByIterationR10_0Fig8.pdf}{Residual vs iteration \( R = 10 \Omega \).}{fig:residualsByIterationR10_0:residualsByIterationR10_0Fig8}{0.3}{ps2a:partH.m}
%
%\paragraph{Grading remarks:} For the \( R = 0.1 \Omega \): ``A log y axis would have been useful here''.
%
Without preconditioning the Gershgorin circles for \( R = 0.1 \Omega \) were the most spread.  Preconditioning reduces the extent of these circles, also reducing the number of iterations required.  All of these had similar post conditioning Gershgorin circles, and post conditioning eigenvalue distributions.  Despite that, it is clear that more than just the spread of the Gershgorin circles is required to account for the rate of decrease of the residual, since the iteration counts seen in the plots for \( R = 1 \Omega, 10 \Omega\) still fall off much faster.
%
\paragraph{Grading remark:} ``Correct. But what we should ultimately look at is the conditioning number. In this case, from the Gershgorin circles, we cannot estimate a lower bound for the min eigenvalue, and hence we cannot estimate an upper bound on the conditioning number.''
}
}