I wonder where is my mistake. This system is linear has only two equations. I want the analytical solution but Sagemath returns an error
var('X Y X0 Y0 t a0 b0 A0 a1 b1 A1')
X = function('X')(t)
Y = function('Y')(t)
de1=(diff(X,t) - a0*X - b0*Y + A0 == 0)
de2=(diff(Y,t) - a1*X - b1*Y + A1 == 0)
show(de1)
show(de2)
sol=desolve_system([de1,de2],[X,Y],[0,X0,Y0],t)
show(sol)
Add some assumptions, for example
assume(b1^2-2*a0*b1+4*a1*b0>0)
assume(b1^2-2*a0*b1+4*a1*b0+a0^2>0)
Alternative to @achresz's answer : change tools. Running :
X0, Y0, t, a0, b0, A0, a1, b1, A1, = SR.var('X0, Y0, t, a0, b0, A0, a1, b1, A1')
X, Y = function("X, Y")
de1 = X(t).diff(t) == a0*X(t)+b0*Y(t)-A0
de2 = Y(t).diff(t) == a1*X(t)+b1*Y(t)-A1
import sympy
%time foo = sympy.dsolve(*map(sympy.sympify, ([de1, de2], [X(t), Y(t)])))
%time tfoo=[u.rewrite("sin") for u in foo]
%time sfoo=[u._sage_() for u in foo]
%time stfoo=[u._sage_() for u in tfoo]
gives :
CPU times: user 1.78 s, sys: 8.09 ms, total: 1.79 s
Wall time: 1.79 s
CPU times: user 24.9 ms, sys: 56 µs, total: 24.9 ms
Wall time: 24.9 ms
CPU times: user 10 ms, sys: 0 ns, total: 10 ms
Wall time: 10 ms
CPU times: user 7.47 ms, sys: 3.89 ms, total: 11.4 ms
Wall time: 11.4 ms
The expressions are somewhat bulky, and neither Sympy nor Sage can lighten them significatively.
In this case, "the competition" (the gratis-but-not-free Wolfram Engine) does better. Running :
M1=mathematica("Sys={X'[t]==a0*X[t]+b0*Y[t]-A0, Y'[t]==a0*X[t]+b1*Y[t]-A1}")
M2=mathematica("TSol=Timing[FullSimplify[DSolve[Sys, {X[t], Y[t]}, t]]]")
gives
CPU times: user 215 µs, sys: 13 µs, total: 228 µs
Wall time: 12.2 s
(half a magnitude slower), but :
sage: M2[2]
{{X[t] -> (A1*b0 - A0*b1)/(a0*b0 - a0*b1) + E^(((a0 + b1)*t)/2)*
(C[1]*Cosh[(Sqrt[a0^2 + 4*a0*b0 - 2*a0*b1 + b1^2]*t)/2] +
(((a0 - b1)*C[1] + 2*b0*C[2])*
Sinh[(Sqrt[a0^2 + 4*a0*b0 - 2*a0*b1 + b1^2]*t)/2])/
Sqrt[a0^2 + 4*a0*b0 - 2*a0*b1 + b1^2]),
Y[t] -> (A0 - A1)/(b0 - b1) + E^(((a0 + b1)*t)/2)*
(C[2]*Cosh[(Sqrt[a0^2 + 4*a0*b0 - 2*a0*b1 + b1^2]*t)/2] +
((2*a0*C[1] - a0*C[2] + b1*C[2])*
Sinh[(Sqrt[a0^2 + 4*a0*b0 - 2*a0*b1 + b1^2]*t)/2])/
Sqrt[a0^2 + 4*a0*b0 - 2*a0*b1 + b1^2])}}
waaaay lighter...
Checking the identity of the solutions symbolically is problematic : (M2[2][1][1][2]-stfoo[0].rhs()._mathematica_()).FullSimplify() doesn't return after several minutes.
Checking them numerically is (lazily) left as an exercise to the reader.
EDIT : Fixed a typo (mis-pasting) in Mathematica's code.
HTH,
Asked: 2024-06-04 18:37:41 +0100
Seen: 265 times
Last updated: Jun 08
Copyright Sage, 2010. Some rights reserved under creative commons license. Content on this site is licensed under a Creative Commons Attribution Share Alike 3.0 license.