In my project, I am computing some Quadrature Rule for Truncated Normal Distribution using the python code here (http://people.sc.fsu.edu/~jburkardt/py_src/truncated_normal_rule/truncated_normal_rule.html). This code involves a Cholesky decomposition of some moment matrix. In some cases, it works and gives me the desired nodes and weights, but in many cases mainly when I increase the number of quadrature points, it says that moment matrix is not longer PSD. I checked the eigenvalues of that matrices and the most of the time I get only one negative value of small order (~ 10^(-6)). Is there any trick to fix this and make the Cholesky decomposition works for my case. Thank you.
Is it possible to slightly modify your moment matrix so that it is PSD, which will make the computation slightly less accurate but still preserve an approximate solution?
The book "Numerical Optimization" by Jorge Nocedal and Stephen J. Wright, which can be found here, provide some techniques for modifying symmetric indefinite matrices. These methods include: eigenvalue modification (you could flip the sign of your negative eigenvalue since it is of extremely small order), adding a multiple of the identity (adding $\tau I$ to your moment matrix, where $\tau$ is a tiny positive real-valued number and $I$ is the identity matrix of the appropriate size), and the Modified Cholesky Factorization method, which has many variants.