Sistem odvisen od parametra

Needs["Graphics`ImplicitPlot`"]

Ičemo presečiča implicitno podanih krivulj 
  x^2+y^2 = r, 
  x^2-2x y+y = 1,
kjer je r parameter.

Ko se radij kronice spreminja, se spreminja tudi tevilo točk, kjer se sekata krivulja in kronica. Z ukazom NSolve lahko poičemo vse reitve, nas pa zanimajo samo realne.

ImplicitPlot[{x^2 + y^2 == 0.5, 
        x^2 - 2  x y + y == 1}, {x, -2, 2}, {y, -2, 
        2}]

FindRoot[{x^2 + y^2 == 0.5, 
        x^2 - 2  x y + y == 1}, {x, -0.2}, {y, 0.7}]

ImplicitPlot[{x^2 + y^2 == 0.6, 
        x^2 - 2  x y + y == 1}, {x, -2, 2}, {y, -2, 
        2}]

Pri malo spremenjenem radiju kronice so prejnja presečiča dobri začetni pribliki za nova presečiča.

FindRoot[{x^2 + y^2 == 0.6, 
        x^2 - 2  x y + y == 1}, {x, -0.2}, {y, 0.67}]

Med radijem 0.3 in radijem 0.4 pride do preskoka v tevilu presečič. Pri mejni vrednosti r bo Jacobi singularen. Tako dobimo e tretjo enačbo za tri neznanke x,y,r. Za začetni pribliek vzamemo točko x,y,r, ki je presečiče pri taki vrednosti r, da je zelo blizu drugemu presečiču.

ImplicitPlot[{x^2 + y^2 == 0.4, 
        x^2 - 2  x y + y == 1}, {x, -2, 2}, {y, -2, 
        2}]

FindRoot[{x^2 + y^2 == 0.4, 
        x^2 - 2  x y + y == 1}, {x, -0.2}, {y, 0.67}]

jac[x_, y_] = 
      Det[{{2  x, 2  y}, {2  x - 2  y, -2 x + 1}}]

Reimo sistem 3 enačb s tremi neznankami. S parametri WorkingPrecision in AccuracyGoal povečamo natančnost računanja, da dobimo r na dovolj decimalnih mest.

N[FindRoot[{x^2 + y^2 == r, x^2 - 2  x y + y == 1, 
          jac[x, y] == 0}, {x, -0.32}, {y, 0.67}, {r, 0.4}, 
        WorkingPrecision -> 20, AccuracyGoal -> 12], 18]

Tako smo izračunali mejno vrednost r, pri kateri pride do preskoka v tevilu presečič.

Wilkinsonov primer

p[x_] = Expand[Product[(x - i), {i, 1, 20}]]

NSolve[p[x] == 0]

FindRoot[p[x] == 0, {x, 12}, WorkingPrecision -> 60, 
      AccuracyGoal -> 25]

FindRoot[p[x] == 0, {x, 16}, WorkingPrecision -> 60, 
      AccuracyGoal -> 25]

g[x_] = p[x] - 2^(-23) x^19

FindRoot[g[x] == 0, {x, 11.79 + 1.65  I}, 
      WorkingPrecision -> 60, AccuracyGoal -> 25]

FindRoot[g[x] == 0, {x, 16.7 - 2.8  I}, 
      WorkingPrecision -> 60, AccuracyGoal -> 25]

NSolve[g[x] == 0]