Artifact 2ee438c57bed8d860fbd2a1ced5ff33b730f088e1e532796e3fb9505526ecfb8:
- Executable file
r38/packages/odesolve/zimopbas.rlg
— part of check-in
[f2fda60abd]
at
2011-09-02 18:13:33
on branch master
— Some historical releases purely for archival purposes
git-svn-id: https://svn.code.sf.net/p/reduce-algebra/code/trunk/historical@1375 2bfe0521-f11c-4a00-b80e-6202646ff360 (user: arthurcnorman@users.sourceforge.net, size: 22944) [annotate] [blame] [check-ins using] [more...]
REDUCE Development Version, Wed Sep 13 20:40:41 2000 ... ODESolve 1.065 % -*- REDUCE -*- % The Postel/Zimmermann (11/4/96) ODE test examples. % Equation names from Postel/Zimmermann. % This version uses Maple-style functional notation wherever possible. % It outputs general solutions of linear ODEs in basis format. % It also checks all solutions. on odesolve_basis, odesolve_check; on div, intstr; off allfac; % to look prettier % 1 Single equations without initial conditions % ============================================== % 1.1 Linear equations % ==================== operator y; % (1) Linear Bernoulli 1 odesolve((x^4-x^3)*df(y(x),x) + 2*x^4*y(x) = x^3/3 + C, y(x), x); - 2*x e {{--------------}, 2 x - 2*x + 1 1 -2 1 1 ---*c*x + ---*x - --- 2 6 4 -------------------------} 2 x - 2*x + 1 % (2) Linear Bernoulli 2 odesolve(-1/2*df(y(x),x) + y(x) = sin x, y(x), x); 2*x 2 4 {{e },---*cos(x) + ---*sin(x)} 5 5 % (3) Linear change of variables (FJW: shifted Euler equation) odesolve(df(y(x),x,2)*(a*x+b)^2 + 4df(y(x),x)*(a*x+b)*a + 2y(x)*a^2 = 0, y(x), x); 2 a a {{----------------------,---------}} 2 2 2 a*x + b a *x + 2*a*b*x + b % (4) Adjoint odesolve((x^2-x)*df(y(x),x,2) + (2x^2+4x-3)*df(y(x),x) + 8x*y(x) = 1, y(x), x); 2 2 {df(y(x),x,2)*x - df(y(x),x,2)*x + 2*df(y(x),x)*x + 4*df(y(x),x)*x - 3*df(y(x),x) + 8*y(x)*x - 1=0} % (5) Polynomial solutions % (FJW: Currently very slow, and fails anyway!) % odesolve((x^2-x)*df(y(x),x,2) + (1-2x^2)*df(y(x),x) + (4x-2)*y(x) = 0, % y(x), x); % (6) Dependent variable missing odesolve(df(y(x),x,2) + 2x*df(y(x),x) = 2x, y(x), x); 1 {{---*sqrt(pi)*erf(x),1},x} 2 % (7) Liouvillian solutions % (FJW: INTEGRATION IMPOSSIBLY SLOW WITHOUT EITHER ALGINT OR NOINT OPTION) begin scalar !*allfac; !*allfac := t; return odesolve((x^3/2-x^2)*df(y(x),x,2) + (2x^2-3x+1)*df(y(x),x) + (x-1)*y(x) = 0, y(x), x, algint); end; -1 1/x - 1/2 - x - 1/2 e *sqrt(x - 2) {{x *e *(x - 2) *int(--------------------------,x), 2 sqrt(x)*x - 2*sqrt(x)*x -1 - 1/2 - x - 1/2 x *e *(x - 2) }} % NB: DO NOT RE-EVALUATE RESULT WITHOUT TURNING ON ALGINT OR NOINT SWITCH % (8) Reduction of order % (FJW: Attempting to make explicit currently too slow.) odesolve(df(y(x),x,2) - 2x*df(y(x),x) + 2y(x) = 3, y(x), x); ***** Cannot convert nonlinear combination solution to basis! {arbconst(2) + sqrt(pi)*arbconst(1) erf(i*x) *int(-----------------------------------------------------------,x)*i 2 x sqrt(pi)*arbconst(1)*erf(i*x)*i*x + e *arbconst(1) - 2*x 1 - 2*int(-----------------------------------------------------------,x) 2 x sqrt(pi)*arbconst(1)*erf(i*x)*i*x + e *arbconst(1) - 2*x 3 - log(y(x) - ---)=0} 2 % (9) Integrating factors % (FJW: Currently very slow, and fails anyway!) % odesolve(sqrt(x)*df(y(x),x,2) + 2x*df(y(x),x) + 3y(x) = 0, y(x), x); % (10) Radical solution (FJW: omitted for now) % (11) Undetermined coefficients odesolve(df(y(x),x,2) - 2/x^2*y(x) = 7x^4 + 3*x^3, y(x), x); -1 2 1 6 1 5 {{x ,x },---*x + ---*x } 4 6 % (12) Variation of parameters odesolve(df(y(x),x,2) + y(x) = csc(x), y(x), x); {{cos(x),sin(x)}, - cos(x)*x + log(sin(x))*sin(x)} % (13) Linear constant coefficients << factor exp(x); write odesolve(df(y(x),x,7) - 14df(y(x),x,6) + 80df(y(x),x,5) - 242df(y(x),x,4) + 419df(y(x),x,3) - 416df(y(x),x,2) + 220df(y(x),x) - 48y(x) = 0, y(x), x); remfac exp(x) >>; 3*x {{e , 4*x e , 2*x e *x, 2*x e , x 2 e *x , x e *x, x e }} % (14) Euler odesolve(df(y(x),x,4) - 4/x^2*df(y(x),x,2) + 8/x^3*df(y(x),x) - 8/x^4*y(x) = 0, y(x), x); -1 2 4 {{x ,x,x ,x }} % (15) Exact n-th order odesolve((1+x+x^2)*df(y(x),x,3) + (3+6x)*df(y(x),x,2) + 6df(y(x),x) = 6x, y(x), x); 1 2 ---*x 2 {{------------, 2 x + x + 1 x ------------, 2 x + x + 1 1 ------------}, 2 x + x + 1 1 4 ---*x 4 ------------} 2 x + x + 1 % 1.2 Nonlinear equations % ======================= % (16) Integrating factors 1 odesolve(df(y(x),x) = y(x)/(y(x)*log y(x) + x), y(x), x); 1 2 {x=arbconst(4)*y(x) + ---*log(y(x)) *y(x)} 2 % (17) Integrating factors 2 odesolve(2y(x)*df(y(x),x)^2 - 2x*df(y(x),x) - y(x) = 0, y(x), x); 4 2 - 1/3 - 2/3 1/3 {{y(x)=2*(4*arbparam(1) - 12*arbparam(1) + 9) *arbparam(1) *2 *arbconst(5)*arbparam(1), 4 2 - 1/3 - 2/3 1/3 x=2*(4*arbparam(1) - 12*arbparam(1) + 9) *arbparam(1) *2 2 4 2 - 1/3 *arbconst(5)*arbparam(1) - (4*arbparam(1) - 12*arbparam(1) + 9) - 2/3 1/3 *arbparam(1) *2 *arbconst(5), arbparam(1)}} % This parametric solution is correct, cf. Zwillinger (1989) p.168 (41.10) % (except that first edition is missing the constant C)! % (18) Bernoulli 1 odesolve(df(y(x),x) + y(x) = y(x)^3*sin x, y(x), x, explicit); {y(x) 2*x - 1/2 =(5*e *arbconst(6) + 2*cos(x) + 4*sin(x)) *sqrt(5)*plus_or_minus(tag_1)} expand_plus_or_minus ws; 2*x - 1/2 {y(x)=(5*e *arbconst(6) + 2*cos(x) + 4*sin(x)) *sqrt(5), 2*x - 1/2 y(x)= - (5*e *arbconst(6) + 2*cos(x) + 4*sin(x)) *sqrt(5)} % (19) Bernoulli 2 operator P, Q; begin scalar soln, !*exp, !*allfac; % for a neat solution on allfac; soln := odesolve(df(y(x),x) + P(x)*y(x) = Q(x)*y(x)^n, y(x), x); off allfac; return soln end; - n (n - 1)*int(p(x),x) {y(x) *y(x)= - e - int(p(x),x)*n + int(p(x),x) *((n - 1)*int(e *q(x),x) - arbconst(7))} odesolve(df(y(x),x) + P(x)*y(x) = Q(x)*y(x)^(2/3), y(x), x); 1/3 - 1/3*int(p(x),x) {y(x) =e *arbconst(8) 1 - 1/3*int(p(x),x) 1/3*int(p(x),x) + ---*e *int(e *q(x),x)} 3 % (20) Clairaut 1 odesolve((x^2-1)*df(y(x),x)^2 - 2x*y(x)*df(y(x),x) + y(x)^2 - 1 = 0, y(x), x, explicit); 2 {y(x)=arbconst(9)*x + sqrt(arbconst(9) + 1), 2 y(x)=arbconst(9)*x - sqrt(arbconst(9) + 1), 2 y(x)=sqrt( - x + 1), 2 y(x)= - sqrt( - x + 1)} % (21) Clairaut 2 operator f, g; odesolve(f(x*df(y(x),x)-y(x)) = g(df(y(x),x)), y(x), x); {f(arbconst(10)*x - y(x)) - g(arbconst(10))=0} % (22) Equations of the form y' = f(x,y) odesolve(df(y(x),x) = (3x^2-y(x)^2-7)/(exp(y(x))+2x*y(x)+1), y(x), x); y(x) 2 3 {arbconst(11) + e + y(x) *x + y(x) - x + 7*x=0} % (23) Homogeneous odesolve(df(y(x),x) = (2x^3*y(x)-y(x)^4)/(x^4-2x*y(x)^3), y(x), x); 3 3 {arbconst(12)*y(x)*x + y(x) + x =0} % (24) Factoring the equation odesolve(df(y(x),x)*(df(y(x),x)+y(x)) = x*(x+y(x)), y(x), x); - x {y(x)=e *arbconst(13) - x + 1, 1 2 y(x)=arbconst(14) + ---*x } 2 % (25) Interchange variables % (NB: Soln in Zwillinger (1989) wrong, as is last eqn in Table 68!) odesolve(df(y(x),x) = x/(x^2*y(x)^2+y(x)^5), y(x), x); 3 2 2/3*y(x) 3 3 {x =e *arbconst(15) - y(x) - ---} 2 % (26) Lagrange 1 odesolve(y(x) = 2x*df(y(x),x) - a*df(y(x),x)^3, y(x), x); -1 1 3 {{y(x)=2*arbconst(16)*arbparam(2) + ---*arbparam(2) *a, 2 -2 3 2 x=arbconst(16)*arbparam(2) + ---*arbparam(2) *a, 4 arbparam(2)}} odesolve(y(x) = 2x*df(y(x),x) - a*df(y(x),x)^3, y(x), x, implicit); 3 2 2 2 2 {64*arbconst(17) *a + 128*arbconst(17) *a*x - 144*arbconst(17)*y(x) *a*x 4 4 2 3 + 64*arbconst(17)*x + 27*y(x) *a - 16*y(x) *x =0} % root_of quartic is VERY slow if explicit option used! % (27) Lagrange 2 odesolve(y(x) = 2x*df(y(x),x) - df(y(x),x)^2, y(x), x); -1 1 2 {{y(x)=2*arbconst(18)*arbparam(3) + ---*arbparam(3) , 3 -2 2 x=arbconst(18)*arbparam(3) + ---*arbparam(3), 3 arbparam(3)}} odesolve(y(x) = 2x*df(y(x),x) - df(y(x),x)^2, y(x), x, implicit); 2 3 3 { - 9*arbconst(19) + 18*arbconst(19)*y(x)*x - 12*arbconst(19)*x - 4*y(x) 2 2 + 3*y(x) *x =0} % (28) Riccati 1 odesolve(df(y(x),x) = exp(x)*y(x)^2 - y(x) + exp(-x), y(x), x); - x - x e *arbconst(20)*sin(x) - e *cos(x) {y(x)=------------------------------------------} arbconst(20)*cos(x) + sin(x) % (29) Riccati 2 << factor x; write odesolve(df(y(x),x) = y(x)^2 - x*y(x) + 1, y(x), x); remfac x >>; 2 1/2*x 2*e *arbconst(21) {y(x)=x + ------------------------------------------------------} - 1/2 sqrt(pi)*sqrt(2)*arbconst(21)*erf(2 *x*i)*i - 2 % (30) Separable odesolve(df(y(x),x) = (9x^8+1)/(y(x)^2+1), y(x), x); 3 9 {3*arbconst(22) + y(x) + 3*y(x) - 3*x - 3*x=0} % (31) Solvable for x odesolve(y(x) = 2x*df(y(x),x) + y(x)*df(y(x),x)^2, y(x), x); -1 {{y(x)= - 2*arbconst(23)*arbparam(4) , -2 x= - arbconst(23)*arbparam(4) + arbconst(23), arbparam(4)}} odesolve(y(x) = 2x*df(y(x),x) + y(x)*df(y(x),x)^2, y(x), x, implicit); 2 2 { - 4*arbconst(24) + 4*arbconst(24)*x + y(x) =0} % (32) Solvable for y begin scalar !*allfac; !*allfac := t; return odesolve(x = y(x)*df(y(x),x) - x*df(y(x),x)^2, y(x), x) end; 2 - 1/2*arbparam(5) 2 {{y(x)=e *arbconst(25)*(arbparam(5) + 1), 2 - 1/2*arbparam(5) x=e *arbconst(25)*arbparam(5), arbparam(5)}} % (33) Autonomous 1 odesolve(df(y(x),x,2)-df(y(x),x) = 2y(x)*df(y(x),x), y(x), x, explicit); {y(x) 1 arbconst(27)*arbconst(26) - arbconst(26)*x 1 = - ---*arbconst(26)*tan(--------------------------------------------) - ---, 2 2 2 y(x)=arbconst(28)} % (34) Autonomous 2 (FJW: Slow without either algint or noint option.) odesolve(df(y(x),x,2)/y(x) - df(y(x),x)^2/y(x)^2 - 1 + 1/y(x)^3 = 0, y(x), x, algint); {arbconst(31)*plus_or_minus(tag_4) + sqrt(3) 3 3 - 1/2 *int(sqrt(y(x))*(3*arbconst(30)*y(x) + 6*log(y(x))*y(x) + 2) ,y(x)) - plus_or_minus(tag_4)*x=0} % (35) Differentiation method odesolve(2y(x)*df(y(x),x,2) - df(y(x),x)^2 = 1/3(df(y(x),x) - x*df(y(x),x,2))^2, y(x), x, explicit); 2 2 2 {y(x)=arbconst(33) *x + 2*sqrt(3)*arbconst(33)*arbconst(32)*x + 4*arbconst(32) , 2 2 2 y(x)=arbconst(34) *x - 2*sqrt(3)*arbconst(34)*arbconst(32)*x + 4*arbconst(32) , y(x)=arbconst(35)} % (36) Equidimensional in x odesolve(x*df(y(x),x,2) = 2y(x)*df(y(x),x), y(x), x, explicit); 1 arbconst(37)*arbconst(36) - arbconst(36)*log(x) {y(x)= - ---*arbconst(36)*tan(-------------------------------------------------) 2 2 1 - ---, 2 y(x)=arbconst(38)} % (37) Equidimensional in y odesolve((1-x)*(y(x)*df(y(x),x,2)-df(y(x),x)^2) + x^2*y(x)^2 = 0, y(x), x); 3 2 arbconst(40) + arbconst(39)*x + 1/6*x + 1/2*x - x x e *(x - 1) {y(x)=---------------------------------------------------------------} x - 1 % (38) Exact second order odesolve(x*y(x)*df(y(x),x,2) + x*df(y(x),x)^2 + y(x)*df(y(x),x) = 0, y(x), x, explicit); {y(x)=sqrt( - arbconst(42) + log(x))*sqrt(arbconst(41))*sqrt(2), y(x)= - sqrt( - arbconst(42) + log(x))*sqrt(arbconst(41))*sqrt(2), y(x)=arbconst(43)} % (39) Factoring differential operator odesolve(df(y(x),x,2)^2 - 2df(y(x),x)*df(y(x),x,2) + 2y(x)*df(y(x),x) - y(x)^2 = 0, y(x), x); x x {y(x)=e *arbconst(45) + e *arbconst(44)*x, x - x y(x)=e *arbconst(47) + e *arbconst(46)} % (40) Scale invariant (fails with algint option) odesolve(x^2*df(y(x),x,2) + 3x*df(y(x),x) = 1/(y(x)^3*x^4), y(x), x); {2*arbconst(49)*plus_or_minus(tag_7) + log( 2 - 1/2 2 - 1/2 - 2*(4*arbconst(48) + 1) *arbconst(48) + (4*arbconst(48) + 1) 2 2 4 4 *sqrt( - 4*arbconst(48)*y(x) *x + y(x) *x - 1) 2 - 1/2 2 2 + (4*arbconst(48) + 1) *y(x) *x ) - 2*log(x)*plus_or_minus(tag_7)=0} % Revised scale-invariant example (hangs with algint option): ode := x^2*df(y(x),x,2) + 3x*df(y(x),x) + 2*y(x) = 1/(y(x)^3*x^4); 2 -3 -4 ode := df(y(x),x,2)*x + 3*df(y(x),x)*x + 2*y(x)=y(x) *x % Choose full (explicit and expanded) solution: odesolve(ode, y(x), x, full); 1 {y(x)= - ---*sqrt(15*arbconst(50) 2 2 - 1/2 -1 - sqrt(225*arbconst(50) - 64)*sin(2*(arbconst(51) - log(x))))*2 *x , 1 y(x)= - ---*sqrt(15*arbconst(50) 2 2 - 1/2 -1 + sqrt(225*arbconst(50) - 64)*sin(2*(arbconst(51) - log(x))))*2 *x , 1 y(x)=---*sqrt(15*arbconst(50) 2 2 - sqrt(225*arbconst(50) - 64)*sin(2*(arbconst(51) - log(x)))) - 1/2 -1 *2 *x , 1 y(x)=---*sqrt(15*arbconst(50) 2 2 + sqrt(225*arbconst(50) - 64)*sin(2*(arbconst(51) - log(x)))) - 1/2 -1 *2 *x } % or "explicit, expand" % Check it -- each solution should simplify to 0: foreach soln in ws collect trigsimp sub(soln, num(lhs ode - rhs ode)); {0,0,0,0} % (41) Autonomous, 3rd order odesolve((df(y(x),x)^2+1)*df(y(x),x,3) - 3df(y(x),x)*df(y(x),x,2)^2 = 0, y(x), x); 2 2 {y(x)=arbconst(55) + sqrt(arbconst(53) *arbconst(52) 2 2 2 - 2*arbconst(53)*arbconst(52) *x + 2*arbconst(53) + arbconst(52) *x - 2*x) -1 *arbconst(52) *i, y(x)=arbconst(56) + i*x, y(x)=arbconst(57) - i*x, y(x)=arbconst(58) + arbconst(54)*x} % odesolve((df(y(x),x)^2+1)*df(y(x),x,3) - 3df(y(x),x)*df(y(x),x,2)^2 = 0, % y(x), x, implicit); % Implicit form is currently too messy! % (42) Autonomous, 4th order odesolve(3*df(y(x),x,2)*df(y(x),x,4) - 5df(y(x),x,3)^2 = 0, y(x), x); {y(x)=arbconst(63)*x + arbconst(62) -3 - 3*sqrt(arbconst(60) - x)*sqrt(6)*arbconst(59) , y(x)=arbconst(65)*x + arbconst(64) -3 + 3*sqrt(arbconst(60) - x)*sqrt(6)*arbconst(59) , 1 2 y(x)=arbconst(67)*x + arbconst(66) + ---*arbconst(61)*x } 2 % 1.3 Special equations % ===================== % (43) Delay odesolve(df(y(x),x) + a*y(x-1) = 0, y(x), x); ***** Arguments of y differ -- solving delay equations is not implemented. % (44) Functions with several parameters odesolve(df(y(x,a),x) = a*y(x,a), y(x,a), x); a*x {{e }} % 2 Single equations with initial conditions % =========================================== % (45) Exact 4th order odesolve(df(y(x),x,4) = sin x, y(x), x, {x=0, y(x)=0, df(y(x),x)=0, df(y(x),x,2)=0, df(y(x),x,3)=0}); 1 3 {y(x)=sin(x) + ---*x - x} 6 % (46) Linear polynomial coefficients -- Bessel J0 odesolve(x*df(y(x),x,2) + df(y(x),x) + 2x*y(x) = 0, y(x), x, {x=0, y(x)=1, df(y(x),x)=0}); {y(x)=besselj(0,sqrt(2)*x)} % (47) Second-degree separable soln := odesolve(x*df(y(x),x)^2 - y(x)^2 + 1 = 0, y(x)=1, x=0, explicit); 1 2*sqrt(x)*plus_or_minus(tag_9) soln := {y(x)=---*e 2 1 - 2*sqrt(x)*plus_or_minus(tag_9) + ---*e } 2 % Alternatively ... soln where e^~x => cosh x + sinh x; {y(x)=cosh(2*sqrt(x)*plus_or_minus(tag_9))} % but this works ONLY with `on div, intstr; off allfac;' % A better alternative is ... trigsimp(soln, hyp, combine); {y(x)=cosh(2*sqrt(x)*plus_or_minus(tag_9))} expand_plus_or_minus ws; {y(x)=cosh(2*sqrt(x))} % (48) Autonomous odesolve(df(y(x),x,2) + y(x)*df(y(x),x)^3 = 0, y(x), x, {x=0, y(x)=0, df(y(x),x)=2}); 3 {y(x) + 3*y(x) - 6*x=0} %% Only one explicit solution satisfies the conditions: begin scalar !*trode, !*fullroots; !*fullroots := t; return odesolve(df(y(x),x,2) + y(x)*df(y(x),x)^3 = 0, y(x), x, {x=0, y(x)=0, df(y(x),x)=2}, explicit); end; 2 1/3 2 - 1/3 {y(x)=(sqrt(9*x + 1) + 3*x) - (sqrt(9*x + 1) + 3*x) } % 3 Systems of equations % ======================= % (49) Integrable combinations operator x, z; odesolve({df(x(t),t) = -3y(t)*z(t), df(y(t),t) = 3x(t)*z(t), df(z(t),t) = -x(t)*y(t)}, {x(t),y(t),z(t)}, t); odesolve-system({df(x(t),t) + 3*y(t)*z(t), df(y(t),t) - 3*x(t)*z(t), df(z(t),t) + x(t)*y(t)},{x(t),y(t),z(t)},t) % (50) Matrix Riccati operator a, b; odesolve({df(x(t),t) = a(t)*(y(t)^2-x(t)^2) + 2b(t)*x(t)*y(t) + 2c*x(t), df(y(t),t) = b(t)*(y(t)^2-x(t)^2) - 2a(t)*x(t)*y(t) + 2c*y(t)}, {x(t),y(t)}, t); 2 2 odesolve-system({a(t)*x(t) - a(t)*y(t) - 2*b(t)*x(t)*y(t) + df(x(t),t) - 2*c*x(t), 2 2 2*a(t)*x(t)*y(t) + b(t)*x(t) - b(t)*y(t) + df(y(t),t) - 2*c*y(t)},{x(t),y(t)},t) % (51) Triangular odesolve({df(x(t),t) = x(t)*(1 + cos(t)/(2+sin(t))), df(y(t),t) = x(t) - y(t)}, {x(t),y(t)}, t); odesolve-system({( - cos(t)*x(t) + df(x(t),t)*sin(t) + 2*df(x(t),t) - sin(t)*x(t) - 2*x(t))/(sin(t) + 2), df(y(t),t) - x(t) + y(t)},{x(t),y(t)},t) % (52) Vector odesolve({df(x(t),t) = 9x(t) + 2y(t), df(y(t),t) = x(t) + 8y(t)}, {x(t),y(t)}, t); odesolve-system({df(x(t),t) - 9*x(t) - 2*y(t), df(y(t),t) - x(t) - 8*y(t)},{x(t),y(t)},t) % (53) Higher order odesolve({df(x(t),t) - x(t) + 2y(t) = 0, df(x(t),t,2) - 2df(y(t),t) = 2t - cos(2t)}, {x(t),y(t)}, t); odesolve-system({df(x(t),t) - x(t) + 2*y(t), cos(2*t) + df(x(t),t,2) - 2*df(y(t),t) - 2*t},{x(t),y(t)},t) % (54) Inhomogeneous system equ := {df(x(t),t) = -1/(t*(t^2+1))*x(t) + 1/(t^2*(t^2+1))*y(t) + 1/t, df(y(t),t) = -t^2/(t^2+1)*x(t) + (2t^2+1)/(t*(t^2+1))*y(t) + 1}; -1 -2 -1 - x(t)*t + y(t)*t + t + t equ := {df(x(t),t)=----------------------------------, 2 t + 1 2 -1 2 - x(t)*t + 2*y(t)*t + y(t)*t + t + 1 df(y(t),t)=-------------------------------------------} 2 t + 1 odesolve(equ, {x(t),y(t)}, t); 2 -1 -1 -2 df(x(t),t)*t + df(x(t),t) - t + t *x(t) - t - y(t)*t odesolve-system({------------------------------------------------------------, 2 t + 1 2 2 2 (df(y(t),t)*t + df(y(t),t) + t *x(t) - t - 2*t*y(t) -1 2 - y(t)*t - 1)/(t + 1)},{x(t),y(t)},t) end; Time for test: 23953 ms, plus GC time: 1807 ms