[Top][All Lists]

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
## Re: [Help-gsl] Recovering from failed Cholesky decomp.

**From**: |
Brian Gough |

**Subject**: |
Re: [Help-gsl] Recovering from failed Cholesky decomp. |

**Date**: |
Tue, 06 Apr 2010 10:38:51 +0100 |

**User-agent**: |
Wanderlust/2.14.0 (Africa) Emacs/22.2 Mule/5.0 (SAKAKI) |

At Wed, 24 Mar 2010 12:38:53 -0600,
Daniel Neilson wrote:
>* I have an application in which I'm using gsl_linalg_cholesky_decomp() *
>* to decompose a symmetric positive-definite matrix, A. However, the *
>* decomp fails in a very small percentage of the matrices that I feed to *
>* this function due to round-off error and all of that fun stuff.*
>* *
>* Does anyone happen to know of a good method that will still give me *
>* the L L^T decomposition for A when one of these errors is detected? For *
>* example, is there a method whereby I can add a small-magnitude diagonal *
>* matrix to A, Cholesky decomp that, and modify the result to obtain the *
>* decomp as though I'd done it with A instead of the modified A?*
In linalg/cholesky.c there is a function which handles all the square roots
double
quiet_sqrt (double x)
/* avoids runtime error, for checking matrix for positive definiteness
*/
{
return (x >= 0) ? sqrt(x) : GSL_NAN;
}
You make a copy of cholesky.c in your application and modify this to
return 0 if x is negative but "close enough" to zero.
--
Brian Gough
GNU Scientific Library -
http://www.gnu.org/software/gsl/

[Prev in Thread] |
**Current Thread** |
[Next in Thread] |

**Re: [Help-gsl] Recovering from failed Cholesky decomp.**,
*Brian Gough* **<=**