Search code examples
calgorithmfunctionformulalogarithm

How to implement natural logarithm with continued fraction in C?


Here I have a little problem. Create something from this formula:

enter image description here

This is what I have, but it doesn't work. Franky, I really don't understand how it should work.. I tried to code it with some bad instructions. N is number of iteration and parts of fraction. I think it leads somehow to recursion but don't know how.

Thanks for any help.

double contFragLog(double z, int n)
{
    double cf = 2 * z;
    double a, b;
    for(int i = n; i >= 1; i--)
    {
        a = sq(i - 2) * sq(z);
        b = i + i - 2;
        cf = a / (b - cf);

    }
    return (1 + cf) / (1 - cf);
}

Solution

  • The central loop is messed. Reworked. Recursion not needed either. Just compute the deepest term first and work your way out.

    double contFragLog(double z, int n) {
      double zz = z*z;
      double cf = 1.0;  // Important this is not 0
      for (int i = n; i >= 1; i--) {
        cf = (2*i -1) - i*i*zz/cf;
      }
      return 2*z/cf;
    }
    
    void testln(double z) {
      double y = log((1+z)/(1-z));
      double y2 = contFragLog(z, 8);
      printf("%e %e %e\n", z, y, y2);
    }
    
    int main() {
      testln(0.2);
      testln(0.5);
      testln(0.8);
      return 0;
    }
    

    Output

    2.000000e-01 4.054651e-01 4.054651e-01
    5.000000e-01 1.098612e+00 1.098612e+00
    8.000000e-01 2.197225e+00 2.196987e+00
    

    [Edit]

    As prompted by @MicroVirus, I found double cf = 1.88*n - 0.95; to work better than double cf = 1.0;. As more terms are used, the value used makes less difference, yet a good initial cf requires fewer terms for a good answer, especially for |z| near 0.5. More work could be done here as I studied 0 < z <= 0.5. @MicroVirus suggestion of 2*n+1 may be close to my suggestion due to an off-by-one of what n is.

    This is based on reverse computing and noting the value of CF[n] as n increased. I was surprised the "seed" value did not appear to be some nice integer equation.