6
$\begingroup$

Description

Recently, I have been learning a couse called "Numerical Analysis". The fixed point iteration theory was introducted to solve the the approximate root of the equations in this course.

Fixed point iteration theory: $x=Ψ(x) \Rightarrow x_{k+1}=Ψ(x_k)$

Example

Solving the approximate root of $f(x)=x^3+4x^2-10 = 0$

So I give three styles of $x=Ψ(x)$, shown as below:

  • $Ψ_1(x)=\frac{1}{2}\sqrt{10-x^3}$

  • $Ψ_2(x)=\sqrt{\frac{10}{x+4}}$

  • $Ψ_3(x)=x-\frac{f(x)}{f'(x)}=x-\frac{x^3+4x^2-10}{3x^2+8x}$

Using Mathematica to implement it and verify the results of textbook

My trial: I set the precision to 10

 result1 = FixedPointList[N[1/2 Sqrt[10 - #^3] &, 10], 1.5`10]

{1.500000000, 1.286953768, 1.402540804, 1.345458374, 1.375170253, 
   1.360094193, 1.367846968, 1.363887004, 1.365916733, 1.364878217, 
   1.365410061, 1.365137821, 1.365277209, 1.365205850, 1.365242384, 
   1.365223680, 1.365233256, 1.365228353, 1.365230863, 1.365229578, 
   1.365230236, 1.365229899, 1.365230072, 1.365229984, 1.365230029, 
   1.365230006, 1.365230017, 1.365230011, 1.365230014, 1.365230013, 
   1.365230014}

  result2 = FixedPointList[N[Sqrt[10/(# + 4)] &, 10], 1.5`10]

{1.500000000, 1.348399725, 1.367376372, 1.364957015, 1.365264748, 
   1.365225594, 1.365230576, 1.365229942, 1.365230023, 1.365230012, 
   1.365230014, 1.365230013}

 result3 = FixedPointList[N[# - (#^3 + 4 #^2 - 10)/(3 #^2 + 8 #) &, 10], 1.5`10]

 {1.500000000, 1.373333333, 1.365262015, 1.36523001, 1.36523001}

However, the last two results of the result3 is 1.36523001, 1.36523001, whose precision is 9 rather than 10

TableForm[
  Flatten[{result1, result2, result3}, {{2}, {1}}], 
  TableHeadings -> {None, Style[#, 15, Red] & /@ {"\[Phi]1 steps", "\[Phi]2 steps", 
   "\[Phi]3 steps"}}]

enter image description here

As the picture shown, the result2 ia same as the result of textbook, however, the result1 and result3 are a little different from the result of textbook. In addition, the precison of result3 is 9, not 10.

Question

Can someone give me a explalation about this trial? For me, I cannot understand this result.

|improve this question
$\endgroup$
10
$\begingroup$

A nice trick to force Mathematica to use a given precision is to use Block and make $MinPrecision equal to $MaxPrecision. So you can write your result1 as:

Block[{$MinPrecision = 10, $MaxPrecision = 10}, 
         FixedPointList[N[1/2 Sqrt[10 - #^3] &, 10], 1.5`10]]

{1.500000000, 1.286953768, 1.402540804, 1.345458374, 1.375170253, 1.360094193,
 1.367846968, 1.363887004, 1.365916733, 1.364878217, 1.365410061, 1.365137821,
 1.365277209, 1.365205850, 1.365242384, 1.365223680, 1.365233256, 1.365228353,
 1.365230863, 1.365229578,  1.365230236, 1.365229899, 1.365230072, 1.365229984,
 1.365230029, 1.365230006, 1.365230017, 1.365230011, 1.365230014, 1.365230013,      
 1.365230014, 1.365230013, 1.365230013, 1.365230013}

And similarly for result3

Block[{$MinPrecision = 10, $MaxPrecision = 10}, 
     FixedPointList[N[# - (#^3 + 4 #^2 - 10)/(3 #^2 + 8 #) &, 10], 1.5`10]]

{1.500000000, 1.373333333, 1.365262015, 1.365230014, 1.365230013, 1.365230013}
|improve this answer
$\endgroup$
4
$\begingroup$

Mathematical operations on numbers of a given precision cannot guarantee the output to maintain the precision of the input numbers. In general precision will decrease. The amount of decrease depends on the operations and some operations will cause precision to decrease significantly. If you want a specific precision on the final result then the working precision should be higher than the desired output precision. Frequently, the working precision is done at twice the desired output precision to allow for significant lose of precision during complex calculations. In your example just using machine precision numbers will perform well.

result1 =
  FixedPointList[N[1/2 Sqrt[10 - #^3] &], 1.5];
{Length[result1], result1[[-1]] // InputForm}

{53, InputForm[1.3652300134140969`]}

result2 =
  FixedPointList[N[Sqrt[10/(# + 4)] &], 1.5];
{Length[result2], result2[[-1]] // InputForm}

{19, InputForm[1.3652300134140969`]}

result3 =
  FixedPointList[N[# - (#^3 + 4 #^2 - 10)/(3 #^2 + 8 #) &], 1.5];
{Length[result3], result3[[-1]] // InputForm}

{6, InputForm[1.3652300134140969`]}

You can shorten the displayed numbers using NumberForm

|improve this answer
$\endgroup$
  • 3
    $\begingroup$ From the explanation above you might conclude that the result of a computation always has less precision than the precision than the numbers used in the computation. However, you increase precision in the next two examples. ------- In[1]:= x=SetPrecision[6.2,4.0]; Precision[Erf[x]] Out[2]= 19.9801 In[3]:= Precision[2+(x-62/10-10^-30)^6] Out[3]=19.5467 $\endgroup$ – Ted Ersek Sep 21 '14 at 19:44
  • $\begingroup$ @Ted Ersek - I believe that I qualified my statements: "cannot guarantee the output to maintain the precision of the input numbers" and " In general precision will decrease" so I believe none of my statements imply "always". On the other hand, I have somewhat frequently encountered situations where increased WorkingPrecision was required to get acceptable results. $\endgroup$ – Bob Hanlon Sep 21 '14 at 23:45
  • $\begingroup$ I didn't intend to say you are wrong. However, one who does not read carefully might think precision always decreases as computation continues. I wanted to make it clear that isn't always the case. $\endgroup$ – Ted Ersek Sep 23 '14 at 0:11

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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