List Archives > 
Maple User Group List Archive > 
Archive by date > 
This Month By Date > 
This Month By Topic
[MUG] Re: Numerical solution of ODE
| [MUG] Re: Numerical solution of ODE |
|
Author: Maple User Group
Posted: Wed, 22 May 2002 16:39:37 -0400
|
>> From: Maple User Group "maple_gr"
-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-
Date: Sat, 18 May 2002 22:10:29 -0700
From: Allan Wittkopf "wittkopf"
Subject: Numerical solution of ODE
To: "maple-list"
Hello PierLuigi, MUG,
The problem here is that dsolve(...,numeric) (in Maple 7 at least, perhaps 6)
needs to know that all functions in the system are either known (like sin(t),
exp(t)) or unknown (any other functions).
It does so by querying `type/PDEtools/known_function`.
The problem here is that you have defined refx1 and refx2, but the type check
does not know this.
You can work around this, knowing one fact, namely that the known function check
looks for a differentiation rule.
You can tell dsolve that 'refx1' and 'refx2' are known then by defining
differentiation rules.
I will assume you are using Maple 7, in which case some of the above can be
simplified - the rules for unevaluated dsolve/numeric procedures were refined
for Maple 7 so the code for refx1 and refx2 is unnecessary:
> ref1:= diff(x1(t),t) = x2(t):
> ref2:= diff(x2(t),t) = (1-x1(t)^2)*x2(t)-x1(t)+
PIECEWISE([1,t<=.5],[-1,t<=1],[1,otherwise]):
# A given solution corresponding to fixed initial condition
> refsol:=dsolve([ref1,ref2,x1(0)=-4/10,x2(0)=6/10],
[x1(t),x2(t)],type=numeric,method=rkf45,
relerr=Float(1,2-Digits),output=listprocedure):
# the two components of the solution as functions of t (NOTE THE CHANGE)
> refx1 := rhs(refsol[2]):
> refx2 := rhs(refsol[3]):
> # The linearization along the given solution
> lin1:=diff(xx1(t),t)=xx2(t):
> lin2:=diff(xx2(t),t)=(-2*refx1(t)*refx2(t)-1)*xx1(t)+(1-refx1(t)^2)*xx2(t):
# Now dsolve/numeric thinks that refx1 is unknown
> `type/PDEtools/known_function`(refx1(t));
false
# We can change this by adding a differentiation rule - note I set it to
undefined intentionally:
> forget(type);
> forget(`type/PDEtools/known_function`);
> `diff/refx1` := proc() return undefined; end:
> `type/PDEtools/known_function`(refx1(t));
true
# Now the same rule for refx2:
> `diff/refx2` := proc() return undefined; end:
# And obtain the solution:
> X_star:=dsolve([lin1,lin2,xx1(0)=1,xx2(0)=0],[xx1(t),xx2(t)],type=numeric,
> method=rkf45,relerr=Float(1,2-Digits));
X_star := proc(rkf45_x) ... end proc
> X_star(1);
bytes used=3012464, alloc=1507052, time=0.10
...
bytes used=16132732, alloc=2489912, time=0.66
[t = 1., xx1(t) = .381455642254406, xx2(t) = -1.68940332483715]
> X_star(2);
bytes used=17132892, alloc=2489912, time=0.70
...
bytes used=48312276, alloc=2752008, time=2.13
[t = 2., xx1(t) = -.606466227250658, xx2(t) = -.0148086678906940]
Additional note: it is not necessary to define the differentiation rules for
Maple 8, the solution is computed directly, but since the type of functions that
refx1,refx2 are is unknown, the code assumes they might be complex-valued, and
we obtain the result:
> X_star(1);
bytes used=4000616, alloc=3538296, time=0.24
[t(1) = 1, xx1(t)(1) = 0.381457714734347 + 0. I,
xx2(t)(1) = -1.68939911251356 + 0. I]
> X_star(2);
[t(2) = 2, xx1(t)(2) = -0.606465529082917 + 0. I,
xx2(t)(2) = -0.0148105934650757 + 0. I]
Where the complex assumption is clear from the presence of the zero imaginary
parts.
Sincerely,
Allan Wittkopf
-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-
From: "Willard, Daniel Dr DUSA-OR" "daniel.willard"
To: "''"
Subject: Numerical solution of ODE
Date: Mon, 20 May 2002 11:43:23 -0400
First of all, eliminate ref1 and x2 to get:
diff(x1(t),t,t)=(1-x1(t)^2)*diff(x1(t),t)-x1(t)+PIECEWISE([1, t <=
> .5],[-1, t <= 1],[1, otherwise]);
>
>From there you are on your own with a DE text, solving it first with
the piecewiae piece = 0.
|
[View Complete Thread]
Previous by date: [MUG] Re: Coupled, nonlinear ODEs., Edgardo S Cheb-Terrab
Next by date: [MUG] Re: Question about texturing using Maple, Maple User Group
Previous thread: [MUG] evaluating _F1 /_F2 in pdsolve o/p, Nagaraj Mahavir
Next thread: [MUG] Question about texturing using Maple, Bguerrieri
|