Mathematica Numerical Questions

  1. How can I coax N[] to give the precision I asked for?

  1. Sometimes N[,] doesn't appear to work. Like here:

       x = {2.0, 3.0, 5.0};
       A = { {6.0, 2.0, 1.0}, {2.0, 3.0, 1.0}, {1.0, 1.0, 1.0} };
       For[ k=0, k < 15, ++k,
          lambda = x.A.x/(x.x);
          y = LinearSolve[A,x];
          x = y / Norm[y,Infinity];
       ]
       N[lambda, 30]
    
    The output is:
       Out[5]= 0.578933
    

    I was expecting 30 digits. Why did N[] ignore my request for 30 digits?


    N can not create additional precission. You calculated lambda with machine precission!

    I rather like "Mathematica for Scientists & Engineers" by Thomas B. Bahder, Addison-Wesley, 1995, 846pp, ISBN 0-201-54090-8. Chapter 8 of this book has a section on vectors and tensors, including MathTensor from Mathsolutions Inc

    FullForm[lambda] will reveal more digits, but in general, your example works with machine precision numbers, and this depends on your machine. One possibility is to set the precision of the numbers, that is 1.0`30 (30 digits precision) or to work with infinite precision, that is, integers and rationals, with no decimal point. For the second case:

    x = {2, 3, 5};
    A = {{6, 2, 1}, {2, 3, 1}, {1, 1, 1}};
    For[k = 0, k < 15, ++k, lambda = x.A.x/(x.x);
      y = LinearSolve[A, x];
      x = y/Norm[y, Infinity];]
    N[lambda, 30]
    
    returns
    0.578933385691052787623495851172
    

    For machine precision input, the output will be a machine precision number.

    x={2.0,3.0,5.0};
    A={{6.0,2.0,1.0},{2.0,3.0,1.0},{1.0,1.0,1.0}};
    For[k=0, k < 15, ++k, lambda=x.A.x/(x.x);
      y=LinearSolve[A,x];
      x=y/Norm[y,Infinity];]
    lambda
    
    0.578933
    
    %//InputForm
    
    0.5789333856910527
    

    Providing higher precision input numbers will give a more precise output. However, some precision is lost in the calculation.

    x={2.0`30,3.0`30,5.0`30};
    A={{6.0`30,2.0`30,1.0`30},{2.0`30,3.0`30,1.0`30},{1.0`30,1.0`30,1.0`30}};
    For[k=0,k<15,++k,lambda=x.A.x/(x.x);
      y=LinearSolve[A,x];
      x=y/Norm[y,Infinity];]
    lambda
    
    0.5789333856910527876234959
    

    For very high precision use exact numbers (Rationalize)

    x = Rationalize[{2.0,3.0,5.0},0];
    A = Rationalize[{{6.0,2.0,1.0},{2.0,3.0,1.0},{1.0,1.0,1.0}},0];
    For[k=0, k < 15, ++k, lambda=x.A.x/(x.x);
      y=LinearSolve[A,x];
      x=y/Norm[y,Infinity];]
    lambda
    
    1102158619423970036337/1903774504398915184457
    
    N[lambda,30]
    
    0.578933385691052787623495851172
    

    Mathematica can't create precision out of nowhere, you only specified 2 digits of precision for the elements of x and A, fortunately Mathematica actually used machine precision numbers and gave you 16, but it can't read your mind and give you 30. Since lambda has 16 digits of precision it's treated as a machine precision number and prints out as such. Here's how to put 30 digits of precision on x and A:

    In[1]:=
    x = SetPrecision[{2.0, 3.0, 5.0},30];
    A = SetPrecision[{ {6.0, 2.0, 1.0}, {2.0, 3.0, 1.0}, {1.0, 1.0, 1.0} },30];
    For[ k=0, k<15, ++k,
        lambda = x.A.x/(x.x);
        y = LinearSolve[A,x];
        x = y / Norm[y,Infinity];
    ]
    N[lambda, 30]
    
    Out[4]=
    0.5789333856910527876234959
    

    N[] is N[] not Clairvoiance[] ;-)

    x = N[{2, 3, 5}, 30];
    A = N[{{6, 2, 1}, {2, 3, 1}, {1, 1, 1}}, 30];
    For[k = 0, k < 15, ++k, lambda = x.A.x/(x.x);
       y = LinearSolve[A, x];
       x = y/Norm[y, Infinity];]
    
    lambda
    
    0.57893338569105278762349585117179881397`24.96959227682693