hw2: done 1.1 and 1.2

This commit is contained in:
Claudio Maggioni (maggicl) 2021-04-10 12:39:47 +02:00
parent 430c222fb7
commit ecc4b8afff
3 changed files with 48 additions and 12 deletions

View File

@ -40,28 +40,55 @@
boundary conditions. As you can see, the matrix $A$ is not symmetric. Does an energy function of
the problem exist? Consider $N = 4$ and show your answer, explaining why it can or cannot exist.}
Answer is a energy function does not exist. Since A is not symmetric
(even if it is pd), the minimizer used for the c.g. method
(i.e. $\frac12 x^T A x - b^T x$ won't work
since $x^T A x$ might be negative and thus the minimizer does not point to
the solution of $Ax = b$ necessairly
The implementation of a function that generates $A$ w.r.t. the parameter $N$ can be found in the MATLAB script
\texttt{main.m} under Section 1.1.
The matrix $A$ and the vector $b$ appear in the following form:
\[A = \begin{bmatrix}
1 & 0 & 0 & 0 & \ldots & 0 \\
-1 & 2 & -1 & 0 & \ldots & 0 \\
0 & -1 & 2 & -1 & \ldots & 0 \\
\vdots & \vdots & \vdots & \vdots & \ddots & \vdots \\
0 & 0 & 0 & 0 & \ldots & 1 \\
\end{bmatrix}\;\;\;b = \begin{bmatrix}0\\h^2\\h^2\\\vdots\\0\end{bmatrix}\]
For $N = 4$, we can attempt to build a minimizer to the solve the system $Ax = b$. In order to find an $x$ such that $Ax = b$, we could define a minimizer $\phi(x)$ such that:
\[\phi(x) = \|b - Ax\|^2\]
Here $\phi(x) = 0$ would mean that $x$ is an exact solution of the system $Ax = b$. We can then attempt to write such minimizer for $N = 4$:
\[\phi(x) = \|b - Ax\|^2 = \left|\begin{bmatrix}0 - x_1\\\frac19 -x_1 + 2x_2 -x_3\\
\frac19 -x_2 +2x_3 -x_4\\0 - x_4\end{bmatrix}\right|^2 =\]\[= x_1^2 + \left(\frac19 - x_1 + 2x_2 - x_3\right)^2 + \left(\frac19 - x_2 + 2x_3 - x_4\right)^2 + x_4^2\]
\[\Delta \phi(x) = \begin{bmatrix}4x_1 - 4x_2 + 2x_3 -\frac29\\
-4x_1 +10x_2 -8x_3 + 2x_4 +\frac29\\
2x_1 -8x_2 +10x_3 + -4x_4 +\frac29\\
2x_2 - 4x_3 + 4x_4 -\frac29\end{bmatrix}\;\;\;\Delta^2 \phi(x) = \begin{bmatrix}4&-4&2&0\\-4&10&-8&2\\2&-8&10&-4\\0&2&-4&4\end{bmatrix}\]
As it can be seen from the Hessian calculation, the Hessian is positive definite forall $x$s. This means, by the sufficient condition of minimizers, that we can find a minimizer by solving $\Delta \phi(x) = 0$ (i.e. finding stationary points in the hyperspace defined by $\phi(x)$. Solving that, we find:
\[x = \begin{bmatrix}0\\\frac19\\\frac19\\0\end{bmatrix}\;\;\]
which is indeed the minimizer and solution of $Ax = b$. Therefore, $\phi(x)$ is a valid energy function and we can say an energy function for the problem exists.
\subsection{Once the new matrix has been derived, write the energy function related to the new problem
and the corresponding gradient and Hessian.}
we already enforce x(1) = x(n) = 0, since b(1) = b(n) = 0 and thus
A(1, :) * x = b(0) = 0 and same for n can be solved only for x(1) = x(n)
= 0size(A, 1)
As by the definitions of $A$ and $b$ given in the assignment, we already enforce $x_1 = x_n = 0$, since the first and the last term of the matrix-vector product $Ax$ are $x_1$ and $x_n$ respectively and the system of equations represented by $Ax = b$ would indeed include the equalities $x_1 = b_1 = 0$ and $x_n = b_n = 0$. Therefore, we can simply alter the matrix $A$ without any need to perform any other transformation in the system.
The objective is therefore $\phi(x) = (1/2)x^T\overline{A}x - b^x$ with a and b
defined above, gradient is = $\overline{A}x - b$, hessian is $= \overline{A}$
We therefore define $\overline{A}$ as a copy of $A$ where $\overline{A}_{2,1} = \overline{A}_{n-1, n} = 0$.
The objective function then becomes $\phi(x) = \frac12 x^T \overline{A} x - b^T x$.
And since the objective is a standard quadratic form, the gradient is $\overline{A}x - b$ while the Hessian is $\overline{A}$.
\subsection{Write the Conjugate Gradient algorithm in the pdf and implement it Matlab code in a function
called \texttt{CGSolve}.}
See page 112 (133 for pdf) for the algorithm implementation
The solution of this task can be found in Section 1.3 of the script \texttt{main.m}.
The solution of this task can be found in Section 1.3 of the script \texttt{main.m} under the function \texttt{CGSolve}.
\subsection{Solve the Poisson problem.}

View File

@ -4,10 +4,19 @@
% Note: exercises are not in the right order due to matlab constraints of
% functions inside of scripts.
A = [1 0 0 0;
-1 2 -1 0;
0 -1 2 1;
0 0 0 1];
b = [0; 1/9; 1/9; 0];
A \ b
clear
clc
close all
plots = 1;
plots = 0;
%% 1.4 - Solution for 1D Poisson for N=1000 using CG