Numerical Solution of least action principle

Question

I am trying to numerically find the path of least action between two points (ignoring the time step normalizing factor): $$S=\sum_i (x_{i+1}-x_{i})^2/2-V(x_i)$$ I don't have the potential in explicit form, but only in a look-up table. The initial $x_0$ and final point $x_N$ are fixed. For finding the path of least action, all that needs to be done is to minimize wrt to all the intermediate points $x_1, x_2, ... x_{N-1}$. If one does this, then this can be shown to reduce to a linear system of equations, where the matrix is in tridiagonal form: $$\begin{bmatrix} -2 & 1 & 0 & 0 &... &0\\ 1 & -2 & 1 & 0 & ... &0\\ 0 & 1 & -2 & 1 & ... &0\\ ...\\ 0 & ... & & 1 & -2 &1\\ 0 & ... & & 0 & 1 &-2\\ \end{bmatrix} \begin{bmatrix} x_1\\ x_2\\ x_3\\ ...\\ x_{N-2}\\ x_{N-1}\\ \end{bmatrix}= \begin{bmatrix} -V'(x_1)-x_0\\ -V'(x_2)\\ -V'(x_3)\\ ...\\ -V'(x_{N-2})\\ -V'(x_{N-1})-x_N\\ \end{bmatrix}$$

I have checked my notes in "Computational Physics". This problem is very similar to solving Poisson's equation in 1D. However the function are the coordinates themselves and the Potential depends on the coordinates. I am not sure how to solve this system of equations.

Answer 1

$$\underbrace{\begin{bmatrix} -2 & 1 & 0 & 0 &... &0\\ 1 & -2 & 1 & 0 & ... &0\\ 0 & 1 & -2 & 1 & ... &0\\ ...\\ 0 & ... & & 1 & -2 &1\\ 0 & ... & & 0 & 1 &-2\\ \end{bmatrix}}_{\mathbf A} \underbrace{\begin{bmatrix} x_1\\ x_2\\ x_3\\ ...\\ x_{N-2}\\ x_{N-1}\\ \end{bmatrix}}_{\mathbf x}=- \underbrace{\begin{bmatrix} V'(x_1)+x_0\\ V'(x_2)\\ V'(x_3)\\ ...\\ V'(x_{N-2})\\ V'(x_{N-1})+x_N\\ \end{bmatrix}}_{\mathbf b}$$

you have to solve this vector equation

$$\mathbf f(\mathbf x)=\mathbf A\,\mathbf x+\mathbf b(\mathbf x)=\mathbf 0\tag 1$$

the numerical solution is

$$\mathbf x_{n+1}=\mathbf x_n-\left[\mathbf F(\mathbf x)^{-1}\,\mathbf f(\mathbf x)\right]_{\mathbf x=\mathbf x_n}\quad n=0~,1~,\ldots$$

where $$\mathbf F=\frac{\partial \mathbf f}{\partial \mathbf x}=\mathbf A+\frac{\partial \mathbf b}{\partial \mathbf x}$$

$~\mathbf x_0~$ is the "start solution" and $~n~$ the iteration steps .

you stop the iteration if $ |f_i(\mathbf x_n)|\le \mathbf \epsilon~$ thus the solution is $~\mathbf x_{n}$