3

EDIT: solved. ans var needed to be set to 0 in every iteration.That was a mistake I overlooked.A mpz_set_ui(ans,0); at the beginning of each loop solves the problem.Thanks to Oli Charlesworth for his effort!

I'm coding a fibonacci function in C using the GMP library,I use the matrix multiplication algorithm,but for some unexpected reason,I get wrong results(I strongly believe that the algorithm is right).It must have something to do with Gmp structs and inits.Any ideas?

/*********************************************************
Find n th Fibonacci number.Matrix algorithm.
+--          --+     +--     --+
|f(n+1) f(n)   |     | 1     1 |  ^n
|              |  == |         |
|f(n)   f(n-1) |     | 1     0 |
+--          --+     +--     --+
*********************************************************/

int
fastfib(mpz_t result,int n)
{


    int cons[2][2]={{1,1},{1,0}};      //"constant" matrix
    mpz_t mat[2][2];                   //actual matrix holding fibonacci numbers
    mpz_t ans;                         //holds f(n+1)   
    mpz_t return_val;                  //used for calculations
    mpz_t temp;                        //used for calculations as well

    //initialize (automatically set to zero)
    mpz_init(ans);
    mpz_init(return_val);
    mpz_init(temp);
    mpz_init(mat[0][0]);
    mpz_init(mat[0][1]);
    mpz_init(mat[1][0]);
    mpz_init(mat[1][1]);

    //start with n=1
    mpz_set_ui(mat[0][0],1);
    mpz_set_ui(mat[1][0],1);
    mpz_set_ui(mat[0][1],1);
    mpz_set_ui(mat[1][1],0);

    //some trivial cases
    if(n==0){
        mpz_set_ui(result,0);
        return 0;
    }
    if(n==1){
        mpz_set_ui(result,1);
        return 0;
    }
    if(n==2){
        mpz_set_ui(result,1);
        return 0;
    }

    n--;


    while(n>1){

        //fib[n+1]
        //ans=mat[0][0]*cons[0][0]+mat[0][1]*cons[1][0];

        mpz_set_ui(ans,0);

        mpz_mul_ui(temp,mat[0][0],cons[0][0]);
        mpz_add(ans,ans,temp);
        mpz_mul_ui(temp,mat[0][1],cons[1][0]);
        mpz_add(ans,ans,temp);

        //update matrix


        mpz_set(mat[1][1],mat[1][0]);   //mat[1][1]=mat[1][0];
        mpz_set(mat[1][0],mat[0][0]);   //mat[0][1]=mat[1][0]=mat[0][0];
        mpz_set(mat[0][1],mat[0][0]);
        mpz_set(mat[0][0],ans);         //mat[0][0]=ans;

        n--;
    }


    //clear vars
    mpz_clear(ans);
    mpz_clear(return_val);
    mpz_clear(temp);

    mpz_set(result,mat[0][0]);
    return 0;

}

Some debug info

results for n=2       // expected
2 1                   //  2 1
1 1                   //  1 1

results for n=3
5 2                   //  3 2
2 1                   //  2 1

results for n=4
12 5                  //  5 3
5 2                   //  3 2

results for n=5
29 12                 //  8 5
12 5                  //  5 3

results for n=6
70 29                 //  13 8
29 12                 //  8  5

70
10
  • Well, at what point does it go wrong? After the first iteration? What are the results you get? Sep 5, 2011 at 10:21
  • I can't break and debug with the gmp structs...I don't know how to get their values.I can only get the final result with a gmp_print.For instance fib 30 gives 44560482149 which certainly is not the case. Sep 5, 2011 at 10:25
  • Why can't you just print out the matrix on every iteration of your while loop? Sep 5, 2011 at 10:31
  • results for n=2 2 1 1 1 results for n=3 5 2 2 1 results for n=4 12 5 5 2 results for n=5 29 12 12 5 results for n=6 70 29 29 12 where I print the first row and then the second Sep 5, 2011 at 10:43
  • That's kind of hard to read. Can you update your answer to include this information, so it can be formatted. Also, you should include a comparison against what the expected output should be, so people don't have to figure it out themselves. Sep 5, 2011 at 10:45

1 Answer 1

0

solved. ans var needed to be set to 0 in every iteration.That was a mistake I overlooked.A mpz_set_ui(ans,0); at the beginning of each loop solves the problem.Thanks to Oli Charlesworth for his effort!

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.