SciCode / 1 / CG

Background:
The conjugate gradient method finds a unique minimizer of the quadratic form
\begin{equation}
f(\mathbf{x})=\frac{1}{2} \mathbf{x}^{\top} \mathbf{A} \mathbf{x}-\mathbf{x}^{\top} \mathbf{A} \mathbf{x}, \quad \mathbf{x} \in \mathbf{R}^n .
\end{equation}
The unique minimizer is evident due to the symmetry and positive definiteness of its Hessian matrix of second derivatives, and the fact that the minimizer, satisfying $\nabla f(x)=\mathbf{A x}-\mathbf{b}=0$, solves the initial problem.

This implies choosing the initial basis vector p₀ as the negation of the gradient of f at x=x₀. The gradient of ff is $Ax−b$. Beginning with an initial guess x₀, this implies setting $p_0=b−Ax_0$. The remaining basis vectors will be conjugate to the gradient, hence the name "conjugate gradient method". Note that p₀ is also the residual generated by this initial algorithm step.

The conjugation constraint is similar to an orthonormality constraint, which allows us to view the algorithm as an instance of Gram-Schmidt orthonormalization. This leads to the following expression:
$$
\mathbf{p}_k=\mathbf{r}_k-\sum_{i<k} \frac{\mathbf{p}_i^{\top} \mathbf{A} \mathbf{p}_k}{\mathbf{p}_i^{\top} \mathbf{A} \mathbf{p}_i} \mathbf{p}_i
$$
The next optimal location is therefore given by
x_(k+1)=x_(k)+α_(k)p_(k)
with
$$
\alpha_k=\frac{\mathbf{p}_k^{\top}\left(\mathbf{b}-\mathbf{A} \mathbf{x}_k\right)}{\mathbf{p}_k^{\top} \mathbf{A} \mathbf{p}_k}=\frac{\mathbf{p}_k^{\top} \mathbf{r}_k}{\mathbf{p}_k^{\top} \mathbf{A} \mathbf{p}_k},
$$

Background:
The conjugate gradient method finds a unique minimizer of the quadratic form
\begin{equation}
f(\mathbf{x})=\frac{1}{2} \mathbf{x}^{\top} \mathbf{A} \mathbf{x}-\mathbf{x}^{\top} \mathbf{A} \mathbf{x}, \quad \mathbf{x} \in \mathbf{R}^n .
\end{equation}
The unique minimizer is evident due to the symmetry and positive definiteness of its Hessian matrix of second derivatives, and the fact that the minimizer, satisfying $\nabla f(x)=\mathbf{A x}-\mathbf{b}=0$, solves the initial problem.

This implies choosing the initial basis vector $p_0$ as the negation of the gradient of $f$ at $x=x_0$. The gradient of ff is $Ax−b$. Beginning with an initial guess $x_0$, this implies setting $p_0=b−Ax_0$. The remaining basis vectors will be conjugate to the gradient, hence the name "conjugate gradient method". Note that $p_0$ is also the residual generated by this initial algorithm step.

The conjugation constraint is similar to an orthonormality constraint, which allows us to view the algorithm as an instance of Gram-Schmidt orthonormalization. This leads to the following expression:
$$
\mathbf{p}_k=\mathbf{r}_k-\sum_{i<k} \frac{\mathbf{p}_i^{\top} \mathbf{A} \mathbf{p}_k}{\mathbf{p}_i^{\top} \mathbf{A} \mathbf{p}_i} \mathbf{p}_i
$$
The next optimal location is therefore given by
$$
\mathbf{x}_{k+1}=\mathbf{x}_k+\alpha_k \mathbf{p}_k
$$
with
$$
\alpha_k=\frac{\mathbf{p}_k^{\top}\left(\mathbf{b}-\mathbf{A} \mathbf{x}_k\right)}{\mathbf{p}_k^{\top} \mathbf{A} \mathbf{p}_k}=\frac{\mathbf{p}_k^{\top} \mathbf{r}_k}{\mathbf{p}_k^{\top} \mathbf{A} \mathbf{p}_k},
$$

Create a function to solve the linear system Ax=b using the conjugate gradient method. This function takes a matrix A and a vector b as inputs.

Create a function to solve the linear system $\mathbf{A} \mathbf{x} = \mathbf{b}$ using the conjugate gradient method. This function takes a matrix $\mathbf{A}$ and a vector $\mathbf{b}$ as inputs.