# HG changeset patch # User John W. Eaton # Date 1440014556 14400 # Node ID 92ac2e05f393048a46e6ad8acff54f2dd4b05c18 # Parent 962724eae09112e261e7eb0099a38ba3475fa2fe * lo-specfun.cc: Reindent comments. diff -r 962724eae091 -r 92ac2e05f393 liboctave/numeric/lo-specfun.cc --- a/liboctave/numeric/lo-specfun.cc Wed Aug 19 16:00:05 2015 -0400 +++ b/liboctave/numeric/lo-specfun.cc Wed Aug 19 16:02:36 2015 -0400 @@ -3004,11 +3004,13 @@ } // This algorithm is due to P. J. Acklam. +// // See http://home.online.no/~pjacklam/notes/invnorm/ -// The rational approximation has relative accuracy 1.15e-9 in the whole region. -// For doubles, it is refined by a single step of Halley's 3rd order method. -// For single precision, the accuracy is already OK, so we skip it to get -// faster evaluation. +// +// The rational approximation has relative accuracy 1.15e-9 in the whole +// region. For doubles, it is refined by a single step of Halley's 3rd +// order method. For single precision, the accuracy is already OK, so +// we skip it to get faster evaluation. static double do_erfinv (double x, bool refine) { @@ -3083,9 +3085,10 @@ return do_erfinv (x, false); } -// The algorthim for erfcinv is an adaptation of the erfinv algorithm above -// from P. J. Acklam. It has been modified to run over the different input -// domain of erfcinv. See the notes for erfinv for an explanation. +// The algorthim for erfcinv is an adaptation of the erfinv algorithm +// above from P. J. Acklam. It has been modified to run over the +// different input domain of erfcinv. See the notes for erfinv for an +// explanation. static double do_erfcinv (double x, bool refine) { @@ -3744,9 +3747,10 @@ static T Lanczos_approximation_psi (const T zc) { - // Coefficients for C.Lanczos expansion of psi function from XLiFE++ gammaFunctions - // psi_coef[k] = - (2k+1) * lg_coef[k] (see melina++ gamma functions) - // -1/12, 3/360,-5/1260, 7/1680,-9/1188, 11*691/360360,-13/156, 15*3617/122400, ? , ? + // Coefficients for C.Lanczos expansion of psi function from XLiFE++ + // gammaFunctions psi_coef[k] = - (2k+1) * lg_coef[k] (see melina++ + // gamma functions -1/12, 3/360,-5/1260, 7/1680,-9/1188, + // 11*691/360360,-13/156, 15*3617/122400, ? , ? static const T dg_coeff[10] = { -0.83333333333333333e-1, 0.83333333333333333e-2, -0.39682539682539683e-2, 0.41666666666666667e-2, @@ -3846,15 +3850,20 @@ unsigned char n = 8 - z_ra; z_m = z + std::complex (n, 0.0); - // Recurrence formula - // for | Re(z) | < 8 , use recursively DiGamma(z) = DiGamma(z+1) - 1/z + // Recurrence formula. For | Re(z) | < 8, use recursively + // + // DiGamma(z) = DiGamma(z+1) - 1/z std::complex z_p = z + P (n - 1); for (unsigned char k = n; k > 0; k--, z_p -= 1.0) dgam -= P (1.0) / z_p; } - // for | Re(z) | > 8, use derivative of C.Lanczos expansion for LogGamma - // psi(z) = log(z) - 1/(2z) - 1/12z^2 + 3/360z^4 - 5/1260z^6 + 7/1680z^8 - 9/1188z^10 + ... + // for | Re(z) | > 8, use derivative of C.Lanczos expansion for + // LogGamma + // + // psi(z) = log(z) - 1/(2z) - 1/12z^2 + 3/360z^4 - 5/1260z^6 + // + 7/1680z^8 - 9/1188z^10 + ... + // // (Abramowitz&Stegun, page 259, formula 6.3.18 dgam += Lanczos_approximation_psi (z_m); }