Re: Differentiating Functions and Root objects [ was Re: ArcCos[]]
- To: mathgroup at smc.vnet.net
- Subject: [mg24857] Re: [mg24833] Differentiating Functions and Root objects [ was Re: ArcCos[]]
- From: Daniel Lichtblau <danl at wolfram.com>
- Date: Wed, 16 Aug 2000 03:24:17 -0400 (EDT)
- References: <8mqvd7$18j@smc.vnet.net> <8mtc0m$7qb@smc.vnet.net> <8n7q9r$1j7@smc.vnet.net> <200008150704.DAA15216@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
Allan Hayes wrote:
>
> John's posting stirred me into taking another look a this problem.
> Here are four observations
>
> 1. We can differentiating a pure function with respect to a parameter
>
> D[Function[x, Sin[y] x^2], y]
>
> Function[x, x^2*Cos[y]]
>
> D[Sin[y] #^2 &, y]
>
> Cos[y]*#1^2 &
>
> 2. The original problem, differentiating Root objects
>
> a = Root[-t + 2*#1 + 2*t^2*#1 + #1^3 &, 1]
>
> Root[-t + 2*#1 + 2*t^2*#1 + #1^3 & , 1]
>
> (b = D[a, t]) // InputForm
>
> (1 - 4*t*Root[-t + 2*#1 + 2*t^2*#1 + #1^3 & , 1])/
> (2 + 2*t^2 + 3*Root[-t + 2*#1 + 2*t^2*#1 + #1^3 & , 1]^
> 2)
>
> Make this output explicit
>
> (rb = ToRadicals[b]) // InputForm
>
> (1 - 4*t*(-(2^(1/3)*(6 + 6*t^2))/
> (3*(27*t + Sqrt[729*t^2 + 4*(6 + 6*t^2)^3])^
> (1/3)) + (27*t + Sqrt[729*t^2 + 4*(6 + 6*t^2)^3])^
> (1/3)/(3*2^(1/3))))/(2 + 2*t^2 +
> 3*(-(2^(1/3)*(6 + 6*t^2))/
> (3*(27*t + Sqrt[729*t^2 + 4*(6 + 6*t^2)^3])^
> (1/3)) +
> (27*t + Sqrt[729*t^2 + 4*(6 + 6*t^2)^3])^(1/3)/
> (3*2^(1/3)))^2)
>
> Maybe Solve and NSolve will work on rb == 0, -- I gave up after a a brief
> wait -- but we can quickly get
>
> FindRoot[rb == 0, {t, 0}]
>
> {t -> -818.986}
>
> Is b correct for the rate of change of a respect to t?
> We can check
>
> rb - D[ToRadicals[a], t] // Simplify
>
> 0
>
> We might expect with sol the corresponding polynomial would have a multiple
> root.
>
> Check:
>
> (poly = -t + 2*#1 + 2*t^2*#1 + #1^3 &[x] /. sol[[1]]) // InputForm
>
> 818.985893312185 + 1.3414777868887153*^6*x + x^3
>
> Solve[% == 0, x] // InputForm
>
> {{x -> -0.0006105102159100438},
> {x -> 0.0003052551079550219 - 1158.2218211072502*I},
> {x -> 0.0003052551079550219 + 1158.2218211072502*I}}
>
> 3. Can we still get the solution without using ToRadicals?
>
> sol2 = FindRoot[b == 0, {t, 0}]
>
> FindRoot::"jsing": "Encountered a singular Jacobian at the point \!\(t\)
> = \!\
> \(0.`\). Try perturbing the initial point(s)."
>
> FindRoot::"jsing": "Encountered a singular Jacobian at the point \!\(t\)
> = \!\
> \(0.`\). Try perturbing the initial point(s)."
>
> FindRoot[b == 0, {t, 0}]
>
> However
>
> sol2 = FindRoot[b == 0, {t, 0.1}]
>
> {t -> -1.0150511440294088}
>
> A different solution.
>
> And we get
>
> Solve[(-t + 2*#1 + 2*t^2*#1 + #1^3 &[x] /. sol2[[1]]) == 0, x]
>
> {{x -> -0.2462928577523092},
> {x -> 0.1231464288761546 - 2.0263644239932934*I},
> {x -> 0.1231464288761546 + 2.0263644239932934*I}}
>
> 4. Is Differentiation of Root objects always valid? Is this method a general
> one for finding which parameter values give multiple roots?
>
> Allan
> ---------------------
> Allan Hayes
> Mathematica Training and Consulting
> Leicester UK
> www.haystack.demon.co.uk
> hay at haystack.demon.co.uk
> Voice: +44 (0)116 271 4198
> Fax: +44 (0)870 164 0565
> [...]
Not sure if this answers all the questions raised above, but at least it
will illustrate another approach to handling these derivatives.
We can differentiate as you did, or by explicitly forming an algebraic
equation and using implicit differentiation. This latter approach gives
an equivalent result but to me it looks a bit cleaner. We then moreover
use it directly and avoid radicals altogether in solving various
equations. Among other things this may help to avoid parasite results
that can arise in the presence of radicals.
In[29]:= b = Root[-t + 2*#1 + 2*t^2*#1 + #1^3 &, 1];
In[30]:= InputForm[ee = b[[1]] [a[t]] ]
Out[30]//InputForm= -t + 2*a[t] + 2*t^2*a[t] + a[t]^3
In[31]:= InputForm[deriv = (a'[t] /. First[Solve[D[ee,t]==0, a'[t]]]) ]
Out[31]//InputForm= (1 - 4*t*a[t])/(2 + 2*t^2 + 3*a[t]^2)
In[32]:= InputForm[Solve[{ee,deriv}==0, {a[t],t}]]
Out[32]//InputForm=
{{t -> (-3*Sqrt[(-4 + 3*Sqrt[2])/2])/2 - Sqrt[-4 + 3*Sqrt[2]],
a[t] -> -Sqrt[-4 + 3*Sqrt[2]]/2},
{t -> (3*Sqrt[(-4 + 3*Sqrt[2])/2])/2 + Sqrt[-4 + 3*Sqrt[2]],
a[t] -> Sqrt[-4 + 3*Sqrt[2]]/2},
{t -> ((3*I)/2)*Sqrt[(4 + 3*Sqrt[2])/2] - I*Sqrt[4 + 3*Sqrt[2]],
a[t] -> (-I/2)*Sqrt[4 + 3*Sqrt[2]]},
{t -> ((-3*I)/2)*Sqrt[(4 + 3*Sqrt[2])/2] + I*Sqrt[4 + 3*Sqrt[2]],
a[t] -> (I/2)*Sqrt[4 + 3*Sqrt[2]]}}
In[33]:= N[%]
Out[33]= {{t -> -1.01505, a[t] -> -0.246293},
> {t -> 1.01505, a[t] -> 0.246293},
> {t -> 0. + 0.174155 I, a[t] -> 0. - 1.4355 I},
> {t -> 0. - 0.174155 I, a[t] -> 0. + 1.4355 I}}
This does not give values where the roots (as functions of t) cross, but
rather values where the functions become horizontal (so the answer to
the second question in (4) is "No"). Here is a brute force approach to
etting crossing points. We define a pair of these root functions, equate
them, and insist that their derivatives be unequal at the points where
they cross (so tangential crossings would need separate attention).
In[67]:= ee1 = b[[1]] [a1[t]];
In[68]:= ee2 = b[[1]] [a2[t]];
In[69]:= deriv1 = (a1'[t] /. First[Solve[D[ee1,t]==0, a1'[t]]]);
In[70]:= deriv2 = (a2'[t] /. First[Solve[D[ee2,t]==0, a2'[t]]]);
In[72]:= InputForm[nsol = NSolve[{ee1,ee2,a1[t]-a2[t],
Numerator[Together[(deriv1-deriv2)*recip-1]]} /.
{a1[t]->a1,a2[t]->a2},
{a1,t}, {recip,a2}]]
Out[72]//InputForm=
{{a1 -> 0.6977990670838363 - 0.5167396037367114*I,
t -> 0.43840749850541516 + 1.233716626570067*I},
{a1 -> 0.6977990670838377 + 0.51673960373671*I,
t -> 0.43840749850541516 - 1.233716626570067*I},
{a1 -> -0.6977990670838368 - 0.5167396037367105*I,
t -> -0.43840749850541516 + 1.2337166265700668*I},
{a1 -> -0.6977990670838368 + 0.5167396037367111*I,
t -> -0.43840749850541516 - 1.2337166265700668*I},
{a1 -> 0.6631797945551418*I, t -> 0.5833428152816322*I},
{a1 -> -0.6631797945551418*I, t -> -0.5833428152816322*I}}
Let's check one of these.
In[77]:= eqn = (ee1 /. a1[t]->a) /. nsol[[1]];
In[78]:= Solve[eqn==0, a]
Out[78]= {{a -> -1.3956 + 1.03348 I}, {a -> 0.697799 - 0.51674 I},
> {a -> 0.697799 - 0.51674 I}}
Works as advertised.
Daniel Lichtblau
Wolfram Research
- References:
- Differentiating Functions and Root objects [ was Re: ArcCos[]]
- From: "Allan Hayes" <hay@haystack.demon.co.uk>
- Differentiating Functions and Root objects [ was Re: ArcCos[]]