# Conjugate gradient method

From https://en.wikipedia.org/wiki/Conjugate_gradient_method:

In mathematics, the conjugate gradient method is an algorithm for the numerical solution of particular systems of linear equations, namely those whose matrix is symmetric and positive-definite.

#### Description
Let $A\mathbf x = \mathbf b$ be a square system of $n$ linear equations.
<br>
<br>
<br>
#### The Algorithm:
<br>
$\begin{align}
& \mathbf{r}_0 := \mathbf{b} - \mathbf{A x}_0 \\
& \hbox{if } \mathbf{r}_{0} \text{ is sufficiently small, then return } \mathbf{x}_{0} \text{ as the result}\\
& \mathbf{p}_0 := \mathbf{r}_0 \\
& k := 0 \\
& \text{repeat} \\
& \qquad \alpha_k := \frac{\mathbf{r}_k^\mathsf{T} \mathbf{r}_k}{\mathbf{p}_k^\mathsf{T} \mathbf{A p}_k}  \\
& \qquad \mathbf{x}_{k+1} := \mathbf{x}_k + \alpha_k \mathbf{p}_k \\
& \qquad \mathbf{r}_{k+1} := \mathbf{r}_k - \alpha_k \mathbf{A p}_k \\
& \qquad \hbox{if } \mathbf{r}_{k+1} \text{ is sufficiently small, then exit loop} \\
& \qquad \beta_k := \frac{\mathbf{r}_{k+1}^\mathsf{T} \mathbf{r}_{k+1}}{\mathbf{r}_k^\mathsf{T} \mathbf{r}_k} \\
& \qquad \mathbf{p}_{k+1} := \mathbf{r}_{k+1} + \beta_k \mathbf{p}_k \\
& \qquad k := k + 1 \\
& \text{end repeat} \\
& \text{return } \mathbf{x}_{k+1} \text{ as the result}
\end{align}$
<br>
<br>

#### Stop Condition:
$$ \lVert X^{(k+1)} - X^{(k)} \rVert_2 \le 10^{-4}$$



```python
import numpy as np
```


```python
def conjgrad(A, b, initial_guess):
    x = initial_guess
    #Symbol of matrix multiplication in numpy is @
    r = b - A@x
    p = r
    
    while(1):
        Ap = A@p
        alpha = (r.T@r)/(p.T@Ap)
        alpha = np.asscalar(alpha)
        x_old = x
        x = x + alpha*p
        x_new = x
        #using norm2:
        if np.linalg.norm(x_new-x_old) < 10**(-4):
            break
        r_old = r
        r = r - alpha*Ap
        beta = (r.T@r)/(r_old.T@r_old)        
        beta = np.asscalar(beta)
        p = r + beta*p
    return x

```


```python
A = np.matrix([[2.0,5.0],
               [5.0,7.0]])

b = np.matrix([[11.0],[13.0]])

initialGuess = np.matrix([[1.0],[1.0]])

sol = conjgrad(A,b,initialGuess)

print ('A:')
print(A)

print ('\nb:')
print(b)

print('\nSolution:')
print(sol)
```

    A:
    [[2. 5.]
     [5. 7.]]
    
    b:
    [[11.]
     [13.]]
    
    Solution:
    [[-1.09090909]
     [ 2.63636364]]



```python

```

2675e05a6d53c8d9c6237cc293615eb67c84b840
ipynb
Jupyter Notebook
.ipynb_checkpoints/Conjugate_Gradient_Method-checkpoint.ipynb
enazari/iterative-methods-for-solvling-linear-systems-of-equations
1529c135923e616b7893b602688a748e037ddffd
MIT
2020-10-02T19:33:07.000Z
2021-01-03T21:04:06.000Z
.ipynb_checkpoints/Conjugate_Gradient_Method-checkpoint.ipynb
enazari/iterative-methods-for-solvling-linear-systems-of-equations
1529c135923e616b7893b602688a748e037ddffd
MIT
.ipynb_checkpoints/Conjugate_Gradient_Method-checkpoint.ipynb
enazari/iterative-methods-for-solvling-linear-systems-of-equations
1529c135923e616b7893b602688a748e037ddffd
MIT
2021-03-01T17:27:50.000Z
2021-03-01T17:27:50.000Z
Qwen/Qwen-72B
1. YES
2. YES