2012-01-24 40 views
2

稍後我問了一下question,這有助於達成解決方案。我已經達成了一個有點可以接受的方法,但仍然沒有完全放在我想要的地方。假設有兩個函數f1[x]g1[y],我想確定公用切線的值爲xy。我至少可以確定切線的一個xy例如下列要求:Mathematica:FindRoot公用切線

f1[x_]:=(5513.12-39931.8x+23307.5x^2+(-32426.6+75662.x-43235.4x^2)Log[(1.-1.33333x)/(1.-1.x)]+x(-10808.9+10808.9x)Log[x/(1.-1.x)])/(-1.+x) 
g1[y_]:=(3632.71+3806.87y-51143.6y^2+y(-10808.9+10808.9y)Log[y/(1.-1.y)]+(-10808.9+32426.6y-21617.7y^2)Log[1.-(1.y)/(1.-1.y)])/(-1.+y) 

Show[ 
Plot[f1[x],{x,0,.75},PlotRange->All], 
Plot[g1[y],{y,0,.75},PlotRange->All] 
] 

Chop[FindRoot[ 
{ 
(f1[x]-g1[y])/(x-y)==D[f1[x],x]==D[g1[y],y] 
}, 
{x,0.0000001},{y,.00000001} 
] 
[[All,2]] 
] 

但是,你會從存在於稍大的xy值另一種常見的切線(積通知例如x〜4和y〜5)。現在,有趣的是,如果我略有f1[x]g1[y]上述表述更改爲類似如下:

f2[x_]:=(7968.08-59377.8x+40298.7x^2+(-39909.6+93122.4x-53212.8x^2)Log[(1.-1.33333x)/(1.-1.x)]+x(-13303.2+13303.2x)Log[x/(1.-1.x)])/(-1.+x) 
    g2[y_]:=(5805.16-27866.2y-21643.y^2+y(-13303.2+13303.2y)Log[y/(1.-1.y)]+(-13303.2+39909.6y-26606.4y^2)Log[1.-(1.y)/(1.-1.y)])/(-1.+y) 

    Show[ 
    Plot[f2[x],{x,0,.75},PlotRange->All], 
    Plot[g2[y],{y,0,.75},PlotRange->All] 
    ] 

    Chop[FindRoot[ 
    { 
    (f2[x]-g2[y])/(x-y)==D[f2[x],x]==D[g2[y],y] 
    }, 
    {x,0.0000001},{y,.00000001} 
    ] 
    [[All,2]] 
    ] 

,並使用相同的方法來確定公切線,數學選擇找到xy較大值正傾斜切線。

最後,我的問題:是否有可能讓Mathematica找到公共切線的高低xy值,並以類似的方式存儲這些值,以便我可以製作列表圖?上述函數fg是另一個變量z的所有複雜函數,而且我正在使用類似下面的內容來繪製切點(應該是兩個x和兩個y),作爲z的函數。

ex[z_]:=Chop[FindRoot[ 
{ 
(f[x,z]-g[y,z])/(x-y)==D[f[x],x]==D[g[y],y] 
}, 
{x,0.0000001},{y,.00000001} 
] 
[[All,2]] 
] 

ListLinePlot[ 
Table[{ex[z][[i]],z},{i,1,2},{z,1300,1800,10}] 
] 

回答

1

好了,讓我們趕緊重寫到目前爲止,你都幹了些什麼:

使用您f1g1,我們有情節

plot = Plot[{f1[x], g1[x]}, {x, 0, .75}] 

da plot

和第一共享正切在

sol1 = Chop[FindRoot[{(f1[x] - g1[y])/(x - y) == D[f1[x], x] == D[g1[y], y]}, 
    {x, 0.0000001}, {y, .00000001}]] 

(* {x -> 0.00840489, y -> 0.105801} *) 

定義功能

l1[t_] = (1 - t) {x, f1[x]} + t {y, g1[y]} /. sol1 

然後,您可以使用

Show[plot, Graphics[Point[{l1[0], l1[1]}]], 
    ParametricPlot[l1[t], {t, -1, 2}], 
    PlotRange -> {{-.2, .4}, {-10000, 10000}}] 

combined plot


我簡要說明(爲了我自己)的方程式您繪製的切線使用 (例如,以上產生sol1) 來自於要求f1的切線x 在某點y切向地命中g1,即,

LogicalExpand[{x, f[x]} + t {1, f'[x]} == {y, g[y]} && f'[x] == g'[y]] 

探討其中共享的切線說謊,你可以使用一個Manipulate

Manipulate[Show[plot, ParametricPlot[{x, f1[x]} + t {1, f1'[x]}, {t, -1, 1}]], 
     {x, 0, .75, Appearance -> "Labeled"}] 

產生類似

animation!

使用直勾勾值xy,您可以使用

sol = Chop[Table[ 
    FindRoot[{(f1[x] - g1[y])/(x - y) == D[f1[x], x] == D[g1[y], y]}, 
    {x, xy[[1]]}, {y, xy[[2]]}], {xy, {{0.001, 0.01}, {0.577, 0.4}}}]] 

使用

l[t_] = (1 - t) {x, f1[x]} + t {y, g1[y]} /. sol 

然後

Show[plot, Graphics[Point[Flatten[{l[0], l[1]}, 1]]], 
ParametricPlot[l[t], {t, -1, 2}, PlotStyle -> Dotted]] 

anudda plot


這個過程湊定義兩條切線得到實際的解決方案ld是自動的,但我不知道如何有效地做到這一點。

3

要找到解決方程的{x, y}的估計值,可以在ContourPlot中繪製它們並查找相交點。例如

f1[x_]:=(5513.12-39931.8 x+23307.5 x^2+(-32426.6+75662. x- 
    43235.4 x^2)Log[(1.-1.33333 x)/(1.-1.x)]+ 
    x(-10808.9+10808.9 x) Log[x/(1.-1.x)])/(-1.+x) 
g1[y_]:=(3632.71+3806.87 y-51143.6 y^2+y (-10808.9+10808.9y) Log[y/(1.-1.y)]+ 
    (-10808.9+32426.6 y-21617.7 y^2) Log[1.-(1.y)/(1.-1.y)])/(-1.+y) 

plot = ContourPlot[{f1'[x] == g1'[y], f1[x] + f1'[x] (y - x) == g1[y]}, 
    {x, 0, 1}, {y, 0, 1}, PlotPoints -> 40] 

Mathematica graphics

正如你可以看到有在間隔(0,1) 2個交點。然後,您可以從圖中讀出點,並利用這些作爲起始值爲FindRoot

seeds = {{.6,.4}, {.05, .1}}; 
sol = FindRoot[{f1'[x] == g1'[y], f1[x] + f1'[x] (y - x) == g1[y]}, 
    {x, #1}, {y, #2}] & @@@ seeds 

拿分的由溶膠的對,你可以使用ReplaceAll

points = {{x, f1[x]}, {y, g1[y]}} /. sol 

(* 
==> {{{0.572412, 19969.9}, {0.432651, 4206.74}}, 
     {{0.00840489, -5747.15}, {0.105801, -7386.68}}} 
*) 

要說明的是,這些是正確的點:

Show[Plot[{f1[x], g1[x]}, {x, 0, 1}], 
{ParametricPlot[#1 t + (1 - t) #2, {t, -5, 5}, PlotStyle -> {Gray, Dashed}], 
    Graphics[{PointSize[Medium], Point[{##}]}]} & @@@ points] 

Mathematica graphics

+0

使用「ContourPlot」的+1 – Simon