Advertisement
MatsGranvik

polynomial roots as a function of z (transformed polynomial)

Jul 24th, 2020
168
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
text 2.03 KB | None | 0 0
  1. (*start*)Clear[A, B, a, b, x];
  2. z = 4;
  3. x = (x1 - z);
  4. Expand[1 + 2 x + x^2/2 + x^3/2 + x^4/6]
  5. a1 = CoefficientList[%, x1]*Range[0, 4]!
  6. a = Flatten[{a1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
  7. 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
  8. 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
  9. 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
  10. 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
  11. 0}];
  12. nn = Length[a];
  13. d = Max[Flatten[Position[Sign[Abs[a]], 1]]]
  14. A = Table[
  15. Table[If[n >= k, Binomial[n - 1, k - 1]*a[[n - k + 1]], 0], {k, 1,
  16. Length[a]}], {n, 1, Length[a]}];
  17. b = Inverse[A][[All, 1]];
  18. x = N[Table[(n - 1)*b[[n - 1]]/b[[n]], {n, nn - 8, nn - 1}], 30]
  19. Sum[a[[k + 1]]/k!*x^k, {k, 0, d}]
  20. -z + x
  21. Clear[x];
  22. Sum[a[[k + 1]]/k!*x^k, {k, 0, d}]
  23. NSolve[1 + 2 x + x^2/2 + x^3/2 + x^4/6, x]
  24. (*end*)
  25.  
  26. This program gives the first root for z =
  27. 4 as - z + x by limiting ratios : -3.0842859641660889018019
  28.  
  29.  
  30. (*start*)Clear[A, B, a, b, x];
  31. z = 0;
  32. x = (x1 - z);
  33. Expand[1 + 2 x + x^2/2 + x^3/2 + x^4/6]
  34. a1 = CoefficientList[%, x1]*Range[0, 4]!
  35. a = Flatten[{a1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
  36. 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
  37. 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
  38. 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
  39. 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
  40. 0}];
  41. nn = Length[a];
  42. d = Max[Flatten[Position[Sign[Abs[a]], 1]]]
  43. A = Table[
  44. Table[If[n >= k, Binomial[n - 1, k - 1]*a[[n - k + 1]], 0], {k, 1,
  45. Length[a]}], {n, 1, Length[a]}];
  46. b = Inverse[A][[All, 1]];
  47. x = N[Table[(n - 1)*b[[n - 1]]/b[[n]], {n, nn - 8, nn - 1}], 30]
  48. Sum[a[[k + 1]]/k!*x^k, {k, 0, d}]
  49. -z + x
  50. Clear[x];
  51. Sum[a[[k + 1]]/k!*x^k, {k, 0, d}]
  52. NSolve[1 + 2 x + x^2/2 + x^3/2 + x^4/6, x]
  53. (*end*)
  54.  
  55. and the second root for z =
  56. 0 as - z +
  57. x also by limiting ratios : -0.5406916617670971775777849727
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement