Re: Re: Solving simple equations
- To: mathgroup at smc.vnet.net
- Subject: [mg83248] Re: [mg83206] Re: Solving simple equations
- From: DrMajorBob <drmajorbob at bigfoot.com>
- Date: Fri, 16 Nov 2007 05:24:22 -0500 (EST)
- References: <fheh0u$mm7$1@smc.vnet.net> <18907167.1195135705102.JavaMail.root@m35>
- Reply-to: drmajorbob at bigfoot.com
Things are not hopeless, even without numerical values for everything.
meqn = {x'[t] == beta (x[t] + (Pf - p[t])),
p'[t] == (Tanh[a1 (x[t] + (Pf - p[t])) + a2 x[t]] - x[t])};
eqns = meqn /. {p'[t] -> 0, x'[t] -> 0}
{0 == beta (Pf - p[t] + x[t]),
0 == Tanh[a2 x[t] + a1 (Pf - p[t] + x[t])] - x[t]}
Notice the first equation has a trivial solution (beta == 0) and a second
solution that simplifies everything:
Quiet@Solve[eqns[[1]]]
{{beta -> 0}, {Pf -> p[t] - x[t]}}
{first, simple} =
Flatten@DeleteCases[eqns /. %, True, {2}] /. x[t] -> x // Simplify
{x == Tanh[a2 x + a1 (Pf + x - p[t])], x == Tanh[a2 x]}
Either of these can be solved numerically, but I'll just take the second.
Solving for x is no help:
Solve[simple, {x}]
Solve::tdep: The equations appear to involve the variables to be \
solved for in an essentially non-algebraic way. >>
Solve[x == Tanh[a2 x], {x}]
but solving for a2 gives us some idea where to start:
Quiet@Solve[simple, {a2}]
{{a2 -> ArcTanh[x]/x}}
NMinimize[ArcTanh[x]/x, {x}]
{1., {x -> -1.45306*10^-8}}
Plot[ArcTanh[x]/x, {x, -5, 5}]
So we see that a2 >= 1.
The ratio is flat at x == 0, so inverting it numerically can be unstable
there. I deal with this using Series:
seriesSoln = x /. Solve[x == Normal@Series[Tanh[a x], {x, 0, 9}], x];
seriesSoln /. a -> 1.1
{0, -0.981186 + 1.01603 \[ImaginaryI],
0.981186- 1.01603 \[ImaginaryI], -0.981186 - 1.01603 \[ImaginaryI],
0.981186+
1.01603 \[ImaginaryI], -0.503006, 0.503006, -1.38767, 1.38767}
0.503006 is the correct root; the others are spurious:
Plot[{x, Tanh[1.1 x]}, {x, 0, 2}]
Finally, x in terms of a can be calculated efficiently by
Clear[xnum]
xnum[a_] /; 1 <= a < 1.1 = seriesSoln[[-3]];
xnum[a_] /; a >= 1.1 := x /. FindRoot[x == Tanh[a x], {x, 1}]
xnum[a_] /; a < 1 := xnum[2 - a]
xint = Quiet@FunctionInterpolation[xnum[a], {a, -5, 5}];
Plot[xint[a], {a, -5, 5}, PlotRange -> All]
If the beta = 0 case is important, I'd start by simplifying
first // ExpandAll
x == Tanh[a1 Pf + a1 x + a2 x - a1 p[t]]
to just
x == Tanh[a x + b]
and applying similar techniques to those above, first finding allowable
ranges for a and b, then using FindRoot and Series to get x values.
I hope that helps!
Bobby
On Thu, 15 Nov 2007 04:32:48 -0600, Jean-Marc Gulliet
<jeanmarc.gulliet at gmail.com> wrote:
> Holger wrote:
>
>> I'm trying to solve these two simple equations but it doesn't work.
>>
>> meqn = {
>> x'[t]==beta(x[t]+(Subscript[P,f]-p[t])),
>> p'[t]==(Tanh[Subscript[a,1](x[t]+(Subscript[P,f]-p[t]))+
>> Subscript[a,2]x[t]]-x[t])
>> };
>>
>> eqp = NSolve[
>> meqn/.{
>> p'[t]->0,x'[t]->0},{p[t], x[t]
>> }
>> ]
>>
>> I guess I'm doing a mistake somewhere. Does anyone have idea where the
>> mistake is?
>
> So you want to solve *numerically* the system meqn (note that I have
> modified the subscripted stuff).
>
> In[1]:= meqn = {x'[t] == beta (x[t] + (Pf - p[t])),
> p'[t] == (Tanh[a1 (x[t] + (Pf - p[t])) + a2 x[t]] - x[t])};
> meqn /. {p'[t] -> 0, x'[t] -> 0}
>
> Out[2]= {0 == beta (Pf - p[t] + x[t]),
> 0 == Tanh[a2 x[t] + a1 (Pf - p[t] + x[t])] - x[t]}
>
> To do so, you must provide some numeric values for *all* the parameter=
s,
> that is beta, Pf, a1, and a2, (though in the above example the value o=
f
> beta does not matter and beta can be discarded before solving).
>
> Of course some (many?) other things might be wrong, but since you did
> tell us neither what you mean by "not working" nor reported any error
> messages returned by by Mathematica, it is hard to tell more (perhaps
> DSolve or NDSolve might be more appropriate).
>
> Regards,
-- =
DrMajorBob at bigfoot.com