Advanced Search

MaplePrimes Questions

Deleted Questions/Posts/Answers...

Asked by: Christopher2222 6030
November 02 2024
3  9

I notice a LOT of posts/answers being deleted.  Why?  What is the purpose of this destructive nature?

To mapleprimes administrators - Please restore deleted content not related to spam. 

If it is some disgruntled mapleprimes user, I'm only guessing, but just because someone didn't reply or upvote a particular answer doesn't mean it isn't good.  Lost knowledge on mapleprimes is just the first sign of a collapsing civilization ... well .. maybe that's a bit of a stretch. 

It still begs the question, why?

How collect series solution of PDE?...

Asked by: salim-barzani 1560
pde differential-equations laplace
November 02 2024
1  31

this example is easiest one for getting solution but i can't collect each part and do like elite i can do each part seperatly but it take to much time i want collect solution and do by easier way if possible this is a laplace adomian decomposition methd which contain adomian polynomial too i want upgrade the code, can any one help the  process for get better vision of this topic 
i do upload some picture and my mw. for more undrestanding

and please can any one explan why when i take laplace why is write D[2](u)(x, 0) must be D[1]?

> 

restart
> 

with(inttrans)
> 

with(PDEtools)
> 

with(LinearAlgebra)
> 

NULL
> 

with(SolveTools)
> 

undeclare(prime)

`There is no more prime differentiation variable; all derivatives will be displayed as indexed functions`

 (1)
> 

with(PDEtools)
> 

declare()

`Declared :`

  

u(x, t)*`to be displayed as`*u

  

`The prime differentiation variable has not been declared yet`

  

`Displayed derivatives and declared functions will be copied and pasted "as they were entered"`

 (2)
> 

declare(u(x, t))

u(x, t)*`will now be displayed as`*u

 (3)
> 

eq := diff(u(x, t), `$`(t, 2))+u(x, t)^2-(diff(u(x, t), x))^2 = 0

diff(diff(u(x, t), t), t)+u(x, t)^2-(diff(u(x, t), x))^2 = 0

 (4)
> 

NULL
> 

eqs := laplace(eq, t, s)

s^2*laplace(u(x, t), t, s)-(D[2](u))(x, 0)-s*u(x, 0)+laplace(u(x, t)^2, t, s)-laplace((diff(u(x, t), x))^2, t, s) = 0

 (5)
> 

solve({eqs}, {laplace(u(x, t), t, s)})

{laplace(u(x, t), t, s) = (s*u(x, 0)+(D[2](u))(x, 0)-laplace(u(x, t)^2, t, s)+laplace((diff(u(x, t), x))^2, t, s))/s^2}

 (6)
> 

subs({u(x, 0) = 0, (D[2](u))(x, 0) = exp(x)}, %)

{laplace(u(x, t), t, s) = (exp(x)-laplace(u(x, t)^2, t, s)+laplace((diff(u(x, t), x))^2, t, s))/s^2}

 (7)
> 

eq3 := invlaplace(%, s, t)

{u(x, t) = exp(x)*t-(int(u(x, _U1)^2*(t-_U1), _U1 = 0 .. t))+int((diff(u(x, _U1), x))^2*(t-_U1), _U1 = 0 .. t)}

 (8)
> 

NULL
> 

NULL
> 

NULL
> 

"u[0](x,t):=exp(x)*t"

proc (x, t) options operator, arrow, function_assign; exp(x)*t end proc

 (9)
> 

n := N

N

 (10)
> 

k := K

K

 (11)
> 

f := proc (u) options operator, arrow; u^2 end proc

proc (u) options operator, arrow; u^2 end proc

 (12)
> 

for j from 0 to 3 do A[j] := subs(lambda = 0, (diff(f(seq(sum(lambda^i*u[i](x, t), i = 0 .. 20), m = 1 .. 2)), [`$`(lambda, j)]))/factorial(j)) end do

(exp(x))^2*t^2

  

2*exp(x)*t*u[1](x, t)

  

u[1](x, t)^2+2*exp(x)*t*u[2](x, t)

  

2*u[1](x, t)*u[2](x, t)+2*exp(x)*t*u[3](x, t)

 (13)
> 

NULL
> 

NULL
> 

n := N

N

 (14)
> 

k := K

K

 (15)
> 

f := proc (u) options operator, arrow; (diff(u, x))^2 end proc

proc (u) options operator, arrow; (diff(u, x))^2 end proc

 (16)
> 

for j from 0 to 3 do B[j] := subs(lambda = 0, (diff(f(seq(sum(lambda^i*u[i](x, t), i = 0 .. 20), m = 1 .. 2)), [`$`(lambda, j)]))/factorial(j)) end do

(exp(x))^2*t^2

  

2*exp(x)*t*(diff(u[1](x, t), x))

  

(diff(u[1](x, t), x))^2+2*exp(x)*t*(diff(u[2](x, t), x))

  

2*(diff(u[1](x, t), x))*(diff(u[2](x, t), x))+2*exp(x)*t*(diff(u[3](x, t), x))

 (17)
> 

"#` know we need find all term of u[0]=exp(x)*t` #` u`[1]=-invlaplace(1/(s^(2))(`A__0`+`B__0`)) u[2]=-invlaplace(1/(s^(2))(`A__1`+`B__1`)) ans so on ... at end i want collect all of them and find final result even if is aproximate and want do test of pde too "
> 

NULL

Download explananing_of_get_solution.mw

Typesetting internal error when Physics update is ...

Asked by: nm 11353
libname physics latex
November 02 2024
0  3

[moderator: The Physics update Library fixes this bug with the same error generated and reported by the same Mapleprimes member on another ODESteps example.]

I have removed Physics update from libname path. 

Now I find I get error calling latex command. When Physics update is on libname, no error.

The error is 

         Error, (in Typesetting:-Parse) too many levels of recursion

I am using worksheet with typesetting extended. But also when I change it to typesetting standard, same error. 

Does this mean one must keep Physics update on libname path for Maple to work OK?

Is this error expected if Physics update is not on libname?

Worksheet below that shows this problem

> 

restart;
> 

interface(version);

`Standard Worksheet Interface, Maple 2024.2, Windows 10, October 29 2024 Build ID 1872373`

CASE 1. With Physics update in libname path, no error

 
> 

restart;

> 

libname;

"C:\Users\Owner\maple\toolbox\2024\Physics Updates\lib", "C:\Program Files\Maple 2024\lib"

> 

ode:=[diff(x__1(t),t)=2*x__1(t)+x__2(t),diff(x__2(t),t)=2*x__1(t)+3*x__2(t)];

[diff(x__1(t), t) = 2*x__1(t)+x__2(t), diff(x__2(t), t) = 2*x__1(t)+3*x__2(t)]

> 

the_output:=Student:-ODEs:-ODESteps(ode,output=typeset):

> 

latex(the_output,'output'=string):

> 

 

CASE 2.  Removing Physics from libname path, gives internal error

 
> 

restart;

> 

libname:=libname[2];

"C:\Program Files\Maple 2024\lib"

> 

ode:=[diff(x__1(t),t)=2*x__1(t)+x__2(t),diff(x__2(t),t)=2*x__1(t)+3*x__2(t)];

[diff(x__1(t), t) = 2*x__1(t)+x__2(t), diff(x__2(t), t) = 2*x__1(t)+3*x__2(t)]

> 

the_output:=Student:-ODEs:-ODESteps(ode,output=typeset):

> 

latex(the_output,'output'=string):

Error, (in Typesetting:-Parse) too many levels of recursion

> 

 

 

 

Download internal_error_from_latex_when_libname_changed_nov_2_2024.mw

plotting numerically solved PDE at specific range ...

Asked by: alfarunner 90
pde numeric pdsolve
November 02 2024
0  1

it's easy to plot a numeric solution of a pde, e.g u(x,t) at a specific time but it's more difficult to plot a solution for a specific x-value as a function of time. How do I do that in Maple?

I find an error...

Asked by: Yuri Zamorano 30
regression statistics incomplete-question
November 01 2024
0  1

I find an error when calculating ExponentialFit(X, Y, x, summarize = embed) function. It's have an error in R-squared and Adjusted R-squared.

Issues with Variational Derivative and Commutator...

Asked by: WA573 80
November 01 2024
0  4

I’m working on deriving the equation of motion (16) from equation (9) (see attached image), but I’ve encountered a couple of issues in Maple. i) The variational derivative of Lagrangian, and ii) the commutator is zero, which isn’t expected.

How do we fix these issues? Are there any specific packages regarding the construction of Lagrangian?

VariationDerivative.mw

How to compute with "exports" from dsolve,numeric...

Asked by: C_R 3412
ode procedure dsolve
November 01 2024
1  12

As an example, I have tried to compute the length of a pendulum to verify how good the lagrange multiplier establishes the constraint of a constant length.

Computation of single values from a formula was possible. Using the same formulas as argument of Maple commands (plot in this case) turned out to be complicated.

I tried many variants with unevaluation quotes before turning to a call to unapply within unevaluation quotes.
There must be simpler ways. Please advise.

PS.: I have put exports in quotes because I assume it would have been easier if dsolve,numeric could export functions. Am I correct in this assumption?

From dsolve,numeric,DAE

> 

dsys := {1 - x(t)^2 - y(t)^2, diff(x(t), t, t) = -2*x(t)*lambda(t), diff(y(t), t, t) = -9.8 - 2*y(t)*lambda(t)}

{1-x(t)^2-y(t)^2, diff(diff(x(t), t), t) = -2*x(t)*lambda(t), diff(diff(y(t), t), t) = -9.8-2*y(t)*lambda(t)}

 (1)
> 

dini := {lambda(0) = 5.025000000, x(0) = 0, y(0) = -1, D(x)(0) = 1/2, D(y)(0) = 0}

{lambda(0) = 5.025000000, x(0) = 0, y(0) = -1, (D(x))(0) = 1/2, (D(y))(0) = 0}

 (2)
> 

dsol1 := dsolve(dsys union dini, numeric);

proc (x_rkf45_dae) local _res, _dat, _vars, _solnproc, _xout, _ndsol, _pars, _n, _i; option `Copyright (c) 2000 by Waterloo Maple Inc. All rights reserved.`; if 1 < nargs then error "invalid input: too many arguments" end if; _EnvDSNumericSaveDigits := Digits; Digits := 15; if _EnvInFsolve = true then _xout := evalf[_EnvDSNumericSaveDigits](x_rkf45_dae) else _xout := evalf(x_rkf45_dae) end if; _dat := Array(1..4, {(1) = proc (_xin) local _xout, _dtbl, _dat, _vmap, _x0, _y0, _val, _dig, _n, _ne, _nd, _nv, _pars, _ini, _par, _i, _j, _k, _src; option `Copyright (c) 2002 by Waterloo Maple Inc. All rights reserved.`; table( [( "complex" ) = false ] ) _xout := _xin; _pars := []; _dtbl := array( 1 .. 4, [( 1 ) = (array( 1 .. 28, [( 1 ) = (datatype = float[8], order = C_order, storage = rectangular), ( 2 ) = (datatype = float[8], order = C_order, storage = rectangular), ( 3 ) = ([0, 0, 0, Array(1..0, 1..2, {}, datatype = float[8], order = C_order)]), ( 4 ) = (Array(1..65, {(1) = 5, (2) = 4, (3) = 0, (4) = 0, (5) = 0, (6) = 2, (7) = 1, (8) = 0, (9) = 0, (10) = 0, (11) = 0, (12) = 0, (13) = 0, (14) = 0, (15) = 0, (16) = 0, (17) = 0, (18) = 1, (19) = 30000, (20) = 0, (21) = 0, (22) = 1, (23) = 4, (24) = 0, (25) = 1, (26) = 15, (27) = 1, (28) = 0, (29) = 1, (30) = 3, (31) = 3, (32) = 0, (33) = 1, (34) = 0, (35) = 0, (36) = 0, (37) = 0, (38) = 0, (39) = 0, (40) = 0, (41) = 0, (42) = 0, (43) = 1, (44) = 0, (45) = 0, (46) = 0, (47) = 0, (48) = 0, (49) = 0, (50) = 50, (51) = 1, (52) = 0, (53) = 0, (54) = 0, (55) = 0, (56) = 0, (57) = 0, (58) = 0, (59) = 10000, (60) = 0, (61) = 1000, (62) = 0, (63) = 0, (64) = -1, (65) = 0}, datatype = integer[8])), ( 5 ) = (Array(1..28, {(1) = .0, (2) = 0.10e-5, (3) = .0, (4) = 0.500001e-14, (5) = .0, (6) = 0.5047658755841546e-2, (7) = .0, (8) = 0.10e-5, (9) = .0, (10) = .0, (11) = .0, (12) = .0, (13) = 1.0, (14) = .0, (15) = .49999999999999, (16) = .0, (17) = 1.0, (18) = 1.0, (19) = .0, (20) = .0, (21) = 1.0, (22) = 1.0, (23) = .0, (24) = .0, (25) = 0.10e-14, (26) = .0, (27) = .0, (28) = .0}, datatype = float[8], order = C_order)), ( 6 ) = (Array(1..5, {(1) = .0, (2) = .5, (3) = -1.0, (4) = .0, (5) = 5.025}, datatype = float[8], order = C_order)), ( 7 ) = ([Array(1..4, 1..7, {(1, 1) = .0, (1, 2) = .203125, (1, 3) = .3046875, (1, 4) = .75, (1, 5) = .8125, (1, 6) = .40625, (1, 7) = .8125, (2, 1) = 0.6378173828125e-1, (2, 2) = .0, (2, 3) = .279296875, (2, 4) = .27237892150878906, (2, 5) = -0.9686851501464844e-1, (2, 6) = 0.1956939697265625e-1, (2, 7) = .5381584167480469, (3, 1) = 0.31890869140625e-1, (3, 2) = .0, (3, 3) = -.34375, (3, 4) = -.335235595703125, (3, 5) = .2296142578125, (3, 6) = .41748046875, (3, 7) = 11.480712890625, (4, 1) = 0.9710520505905151e-1, (4, 2) = .0, (4, 3) = .40350341796875, (4, 4) = 0.20297467708587646e-1, (4, 5) = -0.6054282188415527e-2, (4, 6) = -0.4770040512084961e-1, (4, 7) = .77858567237854}, datatype = float[8], order = C_order), Array(1..6, 1..6, {(1, 1) = .0, (1, 2) = .0, (1, 3) = .0, (1, 4) = .0, (1, 5) = .0, (1, 6) = 1.0, (2, 1) = .25, (2, 2) = .0, (2, 3) = .0, (2, 4) = .0, (2, 5) = .0, (2, 6) = 1.0, (3, 1) = .1875, (3, 2) = .5625, (3, 3) = .0, (3, 4) = .0, (3, 5) = .0, (3, 6) = 2.0, (4, 1) = .23583984375, (4, 2) = -.87890625, (4, 3) = .890625, (4, 4) = .0, (4, 5) = .0, (4, 6) = .2681884765625, (5, 1) = .1272735595703125, (5, 2) = -.5009765625, (5, 3) = .44921875, (5, 4) = -0.128936767578125e-1, (5, 5) = .0, (5, 6) = 0.626220703125e-1, (6, 1) = -0.927734375e-1, (6, 2) = .626220703125, (6, 3) = -.4326171875, (6, 4) = .1418304443359375, (6, 5) = -0.861053466796875e-1, (6, 6) = .3131103515625}, datatype = float[8], order = C_order), Array(1..6, {(1) = .0, (2) = .386, (3) = .21, (4) = .63, (5) = 1.0, (6) = 1.0}, datatype = float[8], order = C_order), Array(1..6, {(1) = .25, (2) = -.1043, (3) = .1035, (4) = -0.362e-1, (5) = .0, (6) = .0}, datatype = float[8], order = C_order), Array(1..6, 1..5, {(1, 1) = .0, (1, 2) = .0, (1, 3) = .0, (1, 4) = .0, (1, 5) = .0, (2, 1) = 1.544, (2, 2) = .0, (2, 3) = .0, (2, 4) = .0, (2, 5) = .0, (3, 1) = .9466785280815533, (3, 2) = .25570116989825814, (3, 3) = .0, (3, 4) = .0, (3, 5) = .0, (4, 1) = 3.3148251870684886, (4, 2) = 2.896124015972123, (4, 3) = .9986419139977808, (4, 4) = .0, (4, 5) = .0, (5, 1) = 1.2212245092262748, (5, 2) = 6.019134481287752, (5, 3) = 12.537083329320874, (5, 4) = -.687886036105895, (5, 5) = .0, (6, 1) = 1.2212245092262748, (6, 2) = 6.019134481287752, (6, 3) = 12.537083329320874, (6, 4) = -.687886036105895, (6, 5) = 1.0}, datatype = float[8], order = C_order), Array(1..6, 1..5, {(1, 1) = .0, (1, 2) = .0, (1, 3) = .0, (1, 4) = .0, (1, 5) = .0, (2, 1) = -5.6688, (2, 2) = .0, (2, 3) = .0, (2, 4) = .0, (2, 5) = .0, (3, 1) = -2.4300933568337584, (3, 2) = -.20635991570891224, (3, 3) = .0, (3, 4) = .0, (3, 5) = .0, (4, 1) = -.10735290581452621, (4, 2) = -9.594562251021896, (4, 3) = -20.470286148096154, (4, 4) = .0, (4, 5) = .0, (5, 1) = 7.496443313968615, (5, 2) = -10.246804314641219, (5, 3) = -33.99990352819906, (5, 4) = 11.708908932061595, (5, 5) = .0, (6, 1) = 8.083246795922411, (6, 2) = -7.981132988062785, (6, 3) = -31.52159432874373, (6, 4) = 16.319305431231363, (6, 5) = -6.0588182388340535}, datatype = float[8], order = C_order), Array(1..3, 1..5, {(1, 1) = .0, (1, 2) = .0, (1, 3) = .0, (1, 4) = .0, (1, 5) = .0, (2, 1) = 10.126235083446911, (2, 2) = -7.487995877607633, (2, 3) = -34.800918615557414, (2, 4) = -7.9927717075687275, (2, 5) = 1.0251377232956207, (3, 1) = -.6762803392806898, (3, 2) = 6.087714651678606, (3, 3) = 16.43084320892463, (3, 4) = 24.767225114183653, (3, 5) = -6.5943891257167815}, datatype = float[8], order = C_order)]), ( 9 ) = ([Array(1..5, {(1) = .1, (2) = .1, (3) = .1, (4) = .1, (5) = .1}, datatype = float[8], order = C_order), Array(1..4, {(1) = .0, (2) = .0, (3) = .0, (4) = .0}, datatype = float[8], order = C_order), Array(1..4, {(1) = .0, (2) = .0, (3) = .0, (4) = .0}, datatype = float[8], order = C_order), Array(1..4, {(1) = .0, (2) = .0, (3) = .0, (4) = .0}, datatype = float[8], order = C_order), Array(1..4, {(1) = .0, (2) = .0, (3) = .0, (4) = .0}, datatype = float[8], order = C_order), Array(1..4, 1..4, {(1, 1) = .0, (1, 2) = .0, (1, 3) = .0, (1, 4) = .0, (2, 1) = .0, (2, 2) = .0, (2, 3) = .0, (2, 4) = .0, (3, 1) = .0, (3, 2) = .0, (3, 3) = .0, (3, 4) = .0, (4, 1) = .0, (4, 2) = .0, (4, 3) = .0, (4, 4) = .0}, datatype = float[8], order = C_order), Array(1..4, 1..4, {(1, 1) = .0, (1, 2) = .0, (1, 3) = -2.0, (1, 4) = .0, (2, 1) = .5, (2, 2) = .0, (2, 3) = .0, (2, 4) = -1.0, (3, 1) = .0, (3, 2) = .0, (3, 3) = .0, (3, 4) = .0, (4, 1) = .0, (4, 2) = .0, (4, 3) = .0, (4, 4) = .0}, datatype = float[8], order = C_order), Array(1..4, {(1) = .0, (2) = .0, (3) = .0, (4) = .0}, datatype = float[8], order = C_order), Array(1..4, 1..4, {(1, 1) = .0, (1, 2) = .0, (1, 3) = .0, (1, 4) = .0, (2, 1) = .0, (2, 2) = .0, (2, 3) = .0, (2, 4) = .0, (3, 1) = .0, (3, 2) = .0, (3, 3) = .0, (3, 4) = .0, (4, 1) = .0, (4, 2) = .0, (4, 3) = .0, (4, 4) = .0}, datatype = float[8], order = C_order), Array(1..5, 1..6, {(1, 1) = .0, (1, 2) = .0, (1, 3) = .0, (1, 4) = .0, (1, 5) = .0, (1, 6) = .0, (2, 1) = .0, (2, 2) = .0, (2, 3) = .0, (2, 4) = .0, (2, 5) = .0, (2, 6) = .0, (3, 1) = .0, (3, 2) = .0, (3, 3) = .0, (3, 4) = .0, (3, 5) = .0, (3, 6) = .0, (4, 1) = .0, (4, 2) = .0, (4, 3) = .0, (4, 4) = .0, (4, 5) = .0, (4, 6) = .0, (5, 1) = .0, (5, 2) = .0, (5, 3) = .0, (5, 4) = .0, (5, 5) = .0, (5, 6) = .0}, datatype = float[8], order = C_order), Array(1..4, {(1) = 0, (2) = 0, (3) = 0, (4) = 0}, datatype = integer[8]), Array(1..5, {(1) = .0, (2) = .0, (3) = .0, (4) = .0, (5) = .0}, datatype = float[8], order = C_order), Array(1..5, {(1) = .0, (2) = .0, (3) = .0, (4) = .0, (5) = .0}, datatype = float[8], order = C_order), Array(1..5, {(1) = .0, (2) = .0, (3) = .0, (4) = .0, (5) = .0}, datatype = float[8], order = C_order), Array(1..5, {(1) = .0, (2) = .0, (3) = .0, (4) = .0, (5) = .0}, datatype = float[8], order = C_order), Array(1..4, {(1) = .0, (2) = .0, (3) = .0, (4) = .0}, datatype = float[8], order = C_order), Array(1..10, {(1) = .0, (2) = .0, (3) = .0, (4) = .0, (5) = .0, (6) = .0, (7) = .0, (8) = .0, (9) = .0, (10) = .0}, datatype = float[8], order = C_order), Array(1..4, {(1) = 0, (2) = 0, (3) = 0, (4) = 0}, datatype = integer[8])]), ( 8 ) = ([Array(1..5, {(1) = .0, (2) = .5, (3) = -1.0, (4) = .0, (5) = undefined}, datatype = float[8], order = C_order), Array(1..5, {(1) = .0, (2) = .0, (3) = .0, (4) = .0, (5) = .0}, datatype = float[8], order = C_order), Array(1..4, {(1) = .5, (2) = -.0, (3) = .0, (4) = .25}, datatype = float[8], order = C_order), 0, 0]), ( 11 ) = (Array(1..6, 0..5, {(1, 1) = .0, (1, 2) = .0, (1, 3) = .0, (1, 4) = .0, (1, 5) = .0, (2, 0) = .0, (2, 1) = .0, (2, 2) = .0, (2, 3) = .0, (2, 4) = .0, (2, 5) = .0, (3, 0) = .0, (3, 1) = .0, (3, 2) = .0, (3, 3) = .0, (3, 4) = .0, (3, 5) = .0, (4, 0) = .0, (4, 1) = .0, (4, 2) = .0, (4, 3) = .0, (4, 4) = .0, (4, 5) = .0, (5, 0) = .0, (5, 1) = .0, (5, 2) = .0, (5, 3) = .0, (5, 4) = .0, (5, 5) = .0, (6, 0) = .0, (6, 1) = .0, (6, 2) = .0, (6, 3) = .0, (6, 4) = .0, (6, 5) = .0}, datatype = float[8], order = C_order)), ( 10 ) = ([proc (N, X, Y, YP) local Z1, Z2, Z3, Z4, Z5; Z1 := Y[2]^2; Z2 := Y[4]^2; Z3 := Z1+Z2; Z4 := Y[1]^2; Z5 := Y[3]^2+Z4; Z5 := 1/Z5; Y[5] := -(1/10)*(-5*Z3+49*Y[3])*Z5; if N < 1 then return 0 end if; YP[2] := -(1/5)*(5*Z3-49*Y[3])*Y[1]*Z5; YP[4] := -(1/5)*(5*Y[3]*(Z1+Z2)+49*Z4)*Z5; YP[1] := Y[2]; YP[3] := Y[4]; 0 end proc, -1, proc (X, Y, R) R[1] := Y[1]^2+Y[3]^2-1; R[2] := Y[1]*Y[2]+Y[3]*Y[4]; 0 end proc, proc (X, Y, J) J[1 .. 2, 1 .. 4] := 0; J[1, 1] := 2*Y[1]; J[1, 3] := 2*Y[3]; J[2, 1] := Y[2]; J[2, 2] := Y[1]; J[2, 3] := Y[4]; J[2, 4] := Y[3]; 0 end proc, 0, 0, 0, 0, 0, 0, 0, 0]), ( 13 ) = (), ( 12 ) = (), ( 15 ) = ("rkf45"), ( 14 ) = ([0, 0]), ( 18 ) = ([]), ( 19 ) = (0), ( 16 ) = ([proc (X, Y, R) R[1] := Y[1]^2+Y[3]^2-1; R[2] := Y[1]*Y[2]+Y[3]*Y[4]; 0 end proc, proc (X, Y, J) J[1 .. 2, 1 .. 4] := 0; J[1, 1] := 2*Y[1]; J[1, 3] := 2*Y[3]; J[2, 1] := Y[2]; J[2, 2] := Y[1]; J[2, 3] := Y[4]; J[2, 4] := Y[3]; 0 end proc, 0, 0, Array(1..4, {(1) = .0, (2) = .0, (3) = .0, (4) = .0}, datatype = float[8], order = C_order), Array(1..2, {(1) = .0, (2) = .0}, datatype = float[8], order = C_order), [-1+x(t)^2+y(t)^2, x(t)*(diff(x(t), t))+y(t)*(diff(y(t), t))]]), ( 17 ) = ([proc (N, X, Y, YP) local Z1, Z2, Z3, Z4, Z5; Z1 := Y[2]^2; Z2 := Y[4]^2; Z3 := Z1+Z2; Z4 := Y[1]^2; Z5 := Y[3]^2+Z4; Z5 := 1/Z5; Y[5] := -(1/10)*(-5*Z3+49*Y[3])*Z5; if N < 1 then return 0 end if; YP[2] := -(1/5)*(5*Z3-49*Y[3])*Y[1]*Z5; YP[4] := -(1/5)*(5*Y[3]*(Z1+Z2)+49*Z4)*Z5; YP[1] := Y[2]; YP[3] := Y[4]; 0 end proc, -1, proc (X, Y, R) R[1] := Y[1]^2+Y[3]^2-1; R[2] := Y[1]*Y[2]+Y[3]*Y[4]; 0 end proc, proc (X, Y, J) J[1 .. 2, 1 .. 4] := 0; J[1, 1] := 2*Y[1]; J[1, 3] := 2*Y[3]; J[2, 1] := Y[2]; J[2, 2] := Y[1]; J[2, 3] := Y[4]; J[2, 4] := Y[3]; 0 end proc, 0, 0, 0, 0, 0, 0, 0, 0]), ( 22 ) = (0), ( 23 ) = (0), ( 20 ) = ([]), ( 21 ) = (0), ( 27 ) = (""), ( 26 ) = (Array(1..0, {})), ( 25 ) = (Array(1..0, {})), ( 24 ) = (0), ( 28 ) = (0) ] )) ] ); _y0 := Array(0..5, {(1) = 0., (2) = 0., (3) = .5000000000000000000000, (4) = -1., (5) = 0.}); _vmap := array( 1 .. 5, [( 1 ) = (5), ( 2 ) = (1), ( 3 ) = (2), ( 4 ) = (3), ( 5 ) = (4) ] ); _x0 := _dtbl[1][5][5]; _n := _dtbl[1][4][1]; _ne := _dtbl[1][4][3]; _nd := _dtbl[1][4][4]; _nv := _dtbl[1][4][16]; if not type(_xout, 'numeric') then if member(_xout, ["start", "left", "right"]) then if _Env_smart_dsolve_numeric = true or _dtbl[1][4][10] = 1 then if _xout = "left" then if type(_dtbl[2], 'table') then return _dtbl[2][5][1] end if elif _xout = "right" then if type(_dtbl[3], 'table') then return _dtbl[3][5][1] end if end if end if; return _dtbl[1][5][5] elif _xout = "method" then return _dtbl[1][15] elif _xout = "storage" then return evalb(_dtbl[1][4][10] = 1) elif _xout = "leftdata" then if not type(_dtbl[2], 'array') then return NULL else return eval(_dtbl[2]) end if elif _xout = "rightdata" then if not type(_dtbl[3], 'array') then return NULL else return eval(_dtbl[3]) end if elif _xout = "enginedata" then return eval(_dtbl[1]) elif _xout = "enginereset" then _dtbl[2] := evaln(_dtbl[2]); _dtbl[3] := evaln(_dtbl[3]); return NULL elif _xout = "initial" then return procname(_y0[0]) elif _xout = "laxtol" then return _dtbl[`if`(member(_dtbl[4], {2, 3}), _dtbl[4], 1)][5][18] elif _xout = "numfun" then return `if`(member(_dtbl[4], {2, 3}), _dtbl[_dtbl[4]][4][18], 0) elif _xout = "parameters" then return [seq(_y0[_n+_i], _i = 1 .. nops(_pars))] elif _xout = "initial_and_parameters" then return procname(_y0[0]), [seq(_y0[_n+_i], _i = 1 .. nops(_pars))] elif _xout = "last" then if _dtbl[4] <> 2 and _dtbl[4] <> 3 or _x0-_dtbl[_dtbl[4]][5][1] = 0. then error "no information is available on last computed point" else _xout := _dtbl[_dtbl[4]][5][1] end if elif _xout = "function" then if _dtbl[1][4][33]-2. = 0 then return eval(_dtbl[1][10], 1) else return eval(_dtbl[1][10][1], 1) end if elif _xout = "map" then return copy(_vmap) elif type(_xin, `=`) and type(rhs(_xin), 'list') and member(lhs(_xin), {"initial", "parameters", "initial_and_parameters"}) then _ini, _par := [], []; if lhs(_xin) = "initial" then _ini := rhs(_xin) elif lhs(_xin) = "parameters" then _par := rhs(_xin) elif select(type, rhs(_xin), `=`) <> [] then _par, _ini := selectremove(type, rhs(_xin), `=`) elif nops(rhs(_xin)) < nops(_pars)+1 then error "insufficient data for specification of initial and parameters" else _par := rhs(_xin)[-nops(_pars) .. -1]; _ini := rhs(_xin)[1 .. -nops(_pars)-1] end if; _xout := lhs(_xout); _i := false; if _par <> [] then _i := `dsolve/numeric/process_parameters`(_n, _pars, _par, _y0) end if; if _ini <> [] then _i := `dsolve/numeric/process_initial`(_n-_ne, _ini, _y0, _pars, _vmap) or _i end if; if _i then `dsolve/numeric/SC/reinitialize`(_dtbl, _y0, _n, procname, _pars); if _Env_smart_dsolve_numeric = true and type(_y0[0], 'numeric') and _dtbl[1][4][10] <> 1 then procname("right") := _y0[0]; procname("left") := _y0[0] end if end if; if _xout = "initial" then return [_y0[0], seq(_y0[_vmap[_i]], _i = 1 .. _n-_ne)] elif _xout = "parameters" then return [seq(_y0[_n+_i], _i = 1 .. nops(_pars))] else return [_y0[0], seq(_y0[_vmap[_i]], _i = 1 .. _n-_ne)], [seq(_y0[_n+_i], _i = 1 .. nops(_pars))] end if elif _xin = "eventstop" then if _nv = 0 then error "this solution has no events" end if; _i := _dtbl[4]; if _i <> 2 and _i <> 3 then return 0 end if; if _dtbl[_i][4][10] = 1 and assigned(_dtbl[5-_i]) and _dtbl[_i][4][9] < 100 and 100 <= _dtbl[5-_i][4][9] then _i := 5-_i; _dtbl[4] := _i; _j := round(_dtbl[_i][4][17]); return round(_dtbl[_i][3][1][_j, 1]) elif 100 <= _dtbl[_i][4][9] then _j := round(_dtbl[_i][4][17]); return round(_dtbl[_i][3][1][_j, 1]) else return 0 end if elif _xin = "eventstatus" then if _nv = 0 then error "this solution has no events" end if; _i := [selectremove(proc (a) options operator, arrow; _dtbl[1][3][1][a, 7] = 1 end proc, {seq(_j, _j = 1 .. round(_dtbl[1][3][1][_nv+1, 1]))})]; return ':-enabled' = _i[1], ':-disabled' = _i[2] elif _xin = "eventclear" then if _nv = 0 then error "this solution has no events" end if; _i := _dtbl[4]; if _i <> 2 and _i <> 3 then error "no events to clear" end if; if _dtbl[_i][4][10] = 1 and assigned(_dtbl[5-_i]) and _dtbl[_i][4][9] < 100 and 100 < _dtbl[5-_i][4][9] then _dtbl[4] := 5-_i; _i := 5-_i end if; if _dtbl[_i][4][9] < 100 then error "no events to clear" elif _nv < _dtbl[_i][4][9]-100 then error "event error condition cannot be cleared" else _j := _dtbl[_i][4][9]-100; if irem(round(_dtbl[_i][3][1][_j, 4]), 2) = 1 then error "retriggerable events cannot be cleared" end if; _j := round(_dtbl[_i][3][1][_j, 1]); for _k to _nv do if _dtbl[_i][3][1][_k, 1] = _j then if _dtbl[_i][3][1][_k, 2] = 3 then error "range events cannot be cleared" end if; _dtbl[_i][3][1][_k, 8] := _dtbl[_i][3][1][_nv+1, 8] end if end do; _dtbl[_i][4][17] := 0; _dtbl[_i][4][9] := 0; if _dtbl[1][4][10] = 1 then if _i = 2 then try procname(procname("left")) catch: end try else try procname(procname("right")) catch: end try end if end if end if; return elif type(_xin, `=`) and member(lhs(_xin), {"eventdisable", "eventenable"}) then if _nv = 0 then error "this solution has no events" end if; if type(rhs(_xin), {('list')('posint'), ('set')('posint')}) then _i := {op(rhs(_xin))} elif type(rhs(_xin), 'posint') then _i := {rhs(_xin)} else error "event identifiers must be integers in the range 1..%1", round(_dtbl[1][3][1][_nv+1, 1]) end if; if select(proc (a) options operator, arrow; _nv < a end proc, _i) <> {} then error "event identifiers must be integers in the range 1..%1", round(_dtbl[1][3][1][_nv+1, 1]) end if; _k := {}; for _j to _nv do if member(round(_dtbl[1][3][1][_j, 1]), _i) then _k := `union`(_k, {_j}) end if end do; _i := _k; if lhs(_xin) = "eventdisable" then _dtbl[4] := 0; _j := [evalb(assigned(_dtbl[2]) and member(_dtbl[2][4][17], _i)), evalb(assigned(_dtbl[3]) and member(_dtbl[3][4][17], _i))]; for _k in _i do _dtbl[1][3][1][_k, 7] := 0; if assigned(_dtbl[2]) then _dtbl[2][3][1][_k, 7] := 0 end if; if assigned(_dtbl[3]) then _dtbl[3][3][1][_k, 7] := 0 end if end do; if _j[1] then for _k to _nv+1 do if _k <= _nv and not type(_dtbl[2][3][4][_k, 1], 'undefined') then userinfo(3, {'events', 'eventreset'}, `reinit #2, event code `, _k, ` to defined init `, _dtbl[2][3][4][_k, 1]); _dtbl[2][3][1][_k, 8] := _dtbl[2][3][4][_k, 1] elif _dtbl[2][3][1][_k, 2] = 0 and irem(iquo(round(_dtbl[2][3][1][_k, 4]), 32), 2) = 1 then userinfo(3, {'events', 'eventreset'}, `reinit #2, event code `, _k, ` to rate hysteresis init `, _dtbl[2][5][24]); _dtbl[2][3][1][_k, 8] := _dtbl[2][5][24] elif _dtbl[2][3][1][_k, 2] = 0 and irem(iquo(round(_dtbl[2][3][1][_k, 4]), 2), 2) = 0 then userinfo(3, {'events', 'eventreset'}, `reinit #2, event code `, _k, ` to initial init `, _x0); _dtbl[2][3][1][_k, 8] := _x0 else userinfo(3, {'events', 'eventreset'}, `reinit #2, event code `, _k, ` to fireinitial init `, _x0-1); _dtbl[2][3][1][_k, 8] := _x0-1 end if end do; _dtbl[2][4][17] := 0; _dtbl[2][4][9] := 0; if _dtbl[1][4][10] = 1 then procname(procname("left")) end if end if; if _j[2] then for _k to _nv+1 do if _k <= _nv and not type(_dtbl[3][3][4][_k, 2], 'undefined') then userinfo(3, {'events', 'eventreset'}, `reinit #3, event code `, _k, ` to defined init `, _dtbl[3][3][4][_k, 2]); _dtbl[3][3][1][_k, 8] := _dtbl[3][3][4][_k, 2] elif _dtbl[3][3][1][_k, 2] = 0 and irem(iquo(round(_dtbl[3][3][1][_k, 4]), 32), 2) = 1 then userinfo(3, {'events', 'eventreset'}, `reinit #3, event code `, _k, ` to rate hysteresis init `, _dtbl[3][5][24]); _dtbl[3][3][1][_k, 8] := _dtbl[3][5][24] elif _dtbl[3][3][1][_k, 2] = 0 and irem(iquo(round(_dtbl[3][3][1][_k, 4]), 2), 2) = 0 then userinfo(3, {'events', 'eventreset'}, `reinit #3, event code `, _k, ` to initial init `, _x0); _dtbl[3][3][1][_k, 8] := _x0 else userinfo(3, {'events', 'eventreset'}, `reinit #3, event code `, _k, ` to fireinitial init `, _x0+1); _dtbl[3][3][1][_k, 8] := _x0+1 end if end do; _dtbl[3][4][17] := 0; _dtbl[3][4][9] := 0; if _dtbl[1][4][10] = 1 then procname(procname("right")) end if end if else for _k in _i do _dtbl[1][3][1][_k, 7] := 1 end do; _dtbl[2] := evaln(_dtbl[2]); _dtbl[3] := evaln(_dtbl[3]); _dtbl[4] := 0; if _dtbl[1][4][10] = 1 then if _x0 <= procname("right") then try procname(procname("right")) catch: end try end if; if procname("left") <= _x0 then try procname(procname("left")) catch: end try end if end if end if; return elif type(_xin, `=`) and lhs(_xin) = "eventfired" then if not type(rhs(_xin), 'list') then error "'eventfired' must be specified as a list" end if; if _nv = 0 then error "this solution has no events" end if; if _dtbl[4] <> 2 and _dtbl[4] <> 3 then error "'direction' must be set prior to calling/setting 'eventfired'" end if; _i := _dtbl[4]; _val := NULL; if not assigned(_EnvEventRetriggerWarned) then _EnvEventRetriggerWarned := false end if; for _k in rhs(_xin) do if type(_k, 'integer') then _src := _k elif type(_k, 'integer' = 'anything') and type(evalf(rhs(_k)), 'numeric') then _k := lhs(_k) = evalf[max(Digits, 18)](rhs(_k)); _src := lhs(_k) else error "'eventfired' entry is not valid: %1", _k end if; if _src < 1 or round(_dtbl[1][3][1][_nv+1, 1]) < _src then error "event identifiers must be integers in the range 1..%1", round(_dtbl[1][3][1][_nv+1, 1]) end if; _src := {seq(`if`(_dtbl[1][3][1][_j, 1]-_src = 0., _j, NULL), _j = 1 .. _nv)}; if nops(_src) <> 1 then error "'eventfired' can only be set/queried for root-finding events and time/interval events" end if; _src := _src[1]; if _dtbl[1][3][1][_src, 2] <> 0. and _dtbl[1][3][1][_src, 2]-2. <> 0. then error "'eventfired' can only be set/queried for root-finding events and time/interval events" elif irem(round(_dtbl[1][3][1][_src, 4]), 2) = 1 then if _EnvEventRetriggerWarned = false then WARNING(`'eventfired' has no effect on events that retrigger`) end if; _EnvEventRetriggerWarned := true end if; if _dtbl[_i][3][1][_src, 2] = 0 and irem(iquo(round(_dtbl[_i][3][1][_src, 4]), 32), 2) = 1 then _val := _val, undefined elif type(_dtbl[_i][3][4][_src, _i-1], 'undefined') or _i = 2 and _dtbl[2][3][1][_src, 8] < _dtbl[2][3][4][_src, 1] or _i = 3 and _dtbl[3][3][4][_src, 2] < _dtbl[3][3][1][_src, 8] then _val := _val, _dtbl[_i][3][1][_src, 8] else _val := _val, _dtbl[_i][3][4][_src, _i-1] end if; if type(_k, `=`) then if _dtbl[_i][3][1][_src, 2] = 0 and irem(iquo(round(_dtbl[_i][3][1][_src, 4]), 32), 2) = 1 then error "cannot set event code for a rate hysteresis event" end if; userinfo(3, {'events', 'eventreset'}, `manual set event code `, _src, ` to value `, rhs(_k)); _dtbl[_i][3][1][_src, 8] := rhs(_k); _dtbl[_i][3][4][_src, _i-1] := rhs(_k) end if end do; return [_val] elif type(_xin, `=`) and lhs(_xin) = "direction" then if not member(rhs(_xin), {-1, 1, ':-left', ':-right'}) then error "'direction' must be specified as either '1' or 'right' (positive) or '-1' or 'left' (negative)" end if; _src := `if`(_dtbl[4] = 2, -1, `if`(_dtbl[4] = 3, 1, undefined)); _i := `if`(member(rhs(_xin), {1, ':-right'}), 3, 2); _dtbl[4] := _i; _dtbl[_i] := `dsolve/numeric/SC/IVPdcopy`(_dtbl[1], `if`(assigned(_dtbl[_i]), _dtbl[_i], NULL)); if 0 < _nv then for _j to _nv+1 do if _j <= _nv and not type(_dtbl[_i][3][4][_j, _i-1], 'undefined') then userinfo(3, {'events', 'eventreset'}, `reinit #4, event code `, _j, ` to defined init `, _dtbl[_i][3][4][_j, _i-1]); _dtbl[_i][3][1][_j, 8] := _dtbl[_i][3][4][_j, _i-1] elif _dtbl[_i][3][1][_j, 2] = 0 and irem(iquo(round(_dtbl[_i][3][1][_j, 4]), 32), 2) = 1 then userinfo(3, {'events', 'eventreset'}, `reinit #4, event code `, _j, ` to rate hysteresis init `, _dtbl[_i][5][24]); _dtbl[_i][3][1][_j, 8] := _dtbl[_i][5][24] elif _dtbl[_i][3][1][_j, 2] = 0 and irem(iquo(round(_dtbl[_i][3][1][_j, 4]), 2), 2) = 0 then userinfo(3, {'events', 'eventreset'}, `reinit #4, event code `, _j, ` to initial init `, _x0); _dtbl[_i][3][1][_j, 8] := _x0 else userinfo(3, {'events', 'eventreset'}, `reinit #4, event code `, _j, ` to fireinitial init `, _x0-2*_i+5.0); _dtbl[_i][3][1][_j, 8] := _x0-2*_i+5.0 end if end do end if; return _src elif _xin = "eventcount" then if _dtbl[1][3][1] = 0 or _dtbl[4] <> 2 and _dtbl[4] <> 3 then return 0 else return round(_dtbl[_dtbl[4]][3][1][_nv+1, 12]) end if elif type(_xin, `=`) and lhs(_xin) = "setdatacallback" then if not type(rhs(_xin), 'nonegint') then error "data callback must be a nonnegative integer (address)" end if; _dtbl[1][28] := rhs(_xin) else return "procname" end if end if; if _xout = _x0 then return [_x0, seq(evalf(_dtbl[1][6][_vmap[_i]]), _i = 1 .. _n-_ne)] end if; _i := `if`(_x0 <= _xout, 3, 2); if _xin = "last" and 0 < _dtbl[_i][4][9] and _dtbl[_i][4][9] < 100 then _dat := eval(_dtbl[_i], 2); _j := _dat[4][20]; return [_dat[11][_j, 0], seq(_dat[11][_j, _vmap[_i]], _i = 1 .. _n-_ne-_nd), seq(_dat[8][1][_vmap[_i]], _i = _n-_ne-_nd+1 .. _n-_ne)] end if; if not type(_dtbl[_i], 'array') then _dtbl[_i] := `dsolve/numeric/SC/IVPdcopy`(_dtbl[1], `if`(assigned(_dtbl[_i]), _dtbl[_i], NULL)); if 0 < _nv then for _j to _nv+1 do if _j <= _nv and not type(_dtbl[_i][3][4][_j, _i-1], 'undefined') then userinfo(3, {'events', 'eventreset'}, `reinit #5, event code `, _j, ` to defined init `, _dtbl[_i][3][4][_j, _i-1]); _dtbl[_i][3][1][_j, 8] := _dtbl[_i][3][4][_j, _i-1] elif _dtbl[_i][3][1][_j, 2] = 0 and irem(iquo(round(_dtbl[_i][3][1][_j, 4]), 32), 2) = 1 then userinfo(3, {'events', 'eventreset'}, `reinit #5, event code `, _j, ` to rate hysteresis init `, _dtbl[_i][5][24]); _dtbl[_i][3][1][_j, 8] := _dtbl[_i][5][24] elif _dtbl[_i][3][1][_j, 2] = 0 and irem(iquo(round(_dtbl[_i][3][1][_j, 4]), 2), 2) = 0 then userinfo(3, {'events', 'eventreset'}, `reinit #5, event code `, _j, ` to initial init `, _x0); _dtbl[_i][3][1][_j, 8] := _x0 else userinfo(3, {'events', 'eventreset'}, `reinit #5, event code `, _j, ` to fireinitial init `, _x0-2*_i+5.0); _dtbl[_i][3][1][_j, 8] := _x0-2*_i+5.0 end if end do end if end if; if _xin <> "last" then if 0 < 0 then if `dsolve/numeric/checkglobals`(op(_dtbl[1][14]), _pars, _n, _y0) then `dsolve/numeric/SC/reinitialize`(_dtbl, _y0, _n, procname, _pars, _i) end if end if; if _dtbl[1][4][7] = 0 then error "parameters must be initialized before solution can be computed" end if end if; _dat := eval(_dtbl[_i], 2); _dtbl[4] := _i; try _src := `dsolve/numeric/SC/IVPrun`(_dat, _xout) catch: userinfo(2, `dsolve/debug`, print(`Exception in solnproc:`, [lastexception][2 .. -1])); error end try; if _dat[17] <> _dtbl[1][17] then _dtbl[1][17] := _dat[17]; _dtbl[1][10] := _dat[10] end if; if _src = 0 and 100 < _dat[4][9] then _val := _dat[3][1][_nv+1, 8] else _val := _dat[11][_dat[4][20], 0] end if; if _src <> 0 or _dat[4][9] <= 0 then _dtbl[1][5][1] := _xout else _dtbl[1][5][1] := _val end if; if _i = 3 and _val < _xout then Rounding := -infinity; if _dat[4][9] = 1 then error "cannot evaluate the solution further right of %1, probably a singularity", evalf[8](_val) elif _dat[4][9] = 2 then error "cannot evaluate the solution further right of %1, maxfun limit exceeded (see <a href='http://www.maplesoft.com/support/help/search.aspx?term=dsolve,maxfun' target='_new'>?dsolve,maxfun</a> for details)", evalf[8](_val) elif _dat[4][9] = 3 then if _dat[4][25] = 3 then error "cannot evaluate the solution past the initial point, problem may be initially singular or improperly set up" else error "cannot evaluate the solution past the initial point, problem may be complex, initially singular or improperly set up" end if elif _dat[4][9] = 4 then error "cannot evaluate the solution further right of %1, accuracy goal cannot be achieved with specified 'minstep'", evalf[8](_val) elif _dat[4][9] = 5 then error "cannot evaluate the solution further right of %1, too many step failures, tolerances may be too loose for problem", evalf[8](_val) elif _dat[4][9] = 6 then error "cannot evaluate the solution further right of %1, cannot downgrade delay storage for problems with delay derivative order > 1, try increasing delaypts", evalf[8](_val) elif _dat[4][9] = 10 then error "cannot evaluate the solution further right of %1, interrupt requested", evalf[8](_val) elif 100 < _dat[4][9] then if _dat[4][9]-100 = _nv+1 then error "constraint projection failure on event at t=%1", evalf[8](_val) elif _dat[4][9]-100 = _nv+2 then error "index-1 and derivative evaluation failure on event at t=%1", evalf[8](_val) elif _dat[4][9]-100 = _nv+3 then error "maximum number of event iterations reached (%1) at t=%2", round(_dat[3][1][_nv+1, 3]), evalf[8](_val) else if _Env_dsolve_nowarnstop <> true then `dsolve/numeric/warning`(StringTools:-FormatMessage("cannot evaluate the solution further right of %1, event #%2 triggered a halt", evalf[8](_val), round(_dat[3][1][_dat[4][9]-100, 1]))) end if; Rounding := 'nearest'; _xout := _val end if else error "cannot evaluate the solution further right of %1", evalf[8](_val) end if elif _i = 2 and _xout < _val then Rounding := infinity; if _dat[4][9] = 1 then error "cannot evaluate the solution further left of %1, probably a singularity", evalf[8](_val) elif _dat[4][9] = 2 then error "cannot evaluate the solution further left of %1, maxfun limit exceeded (see <a href='http://www.maplesoft.com/support/help/search.aspx?term=dsolve,maxfun' target='_new'>?dsolve,maxfun</a> for details)", evalf[8](_val) elif _dat[4][9] = 3 then if _dat[4][25] = 3 then error "cannot evaluate the solution past the initial point, problem may be initially singular or improperly set up" else error "cannot evaluate the solution past the initial point, problem may be complex, initially singular or improperly set up" end if elif _dat[4][9] = 4 then error "cannot evaluate the solution further left of %1, accuracy goal cannot be achieved with specified 'minstep'", evalf[8](_val) elif _dat[4][9] = 5 then error "cannot evaluate the solution further left of %1, too many step failures, tolerances may be too loose for problem", evalf[8](_val) elif _dat[4][9] = 6 then error "cannot evaluate the solution further left of %1, cannot downgrade delay storage for problems with delay derivative order > 1, try increasing delaypts", evalf[8](_val) elif _dat[4][9] = 10 then error "cannot evaluate the solution further right of %1, interrupt requested", evalf[8](_val) elif 100 < _dat[4][9] then if _dat[4][9]-100 = _nv+1 then error "constraint projection failure on event at t=%1", evalf[8](_val) elif _dat[4][9]-100 = _nv+2 then error "index-1 and derivative evaluation failure on event at t=%1", evalf[8](_val) elif _dat[4][9]-100 = _nv+3 then error "maximum number of event iterations reached (%1) at t=%2", round(_dat[3][1][_nv+1, 3]), evalf[8](_val) else if _Env_dsolve_nowarnstop <> true then `dsolve/numeric/warning`(StringTools:-FormatMessage("cannot evaluate the solution further left of %1, event #%2 triggered a halt", evalf[8](_val), round(_dat[3][1][_dat[4][9]-100, 1]))) end if; Rounding := 'nearest'; _xout := _val end if else error "cannot evaluate the solution further left of %1", evalf[8](_val) end if end if; if _EnvInFsolve = true then _dig := _dat[4][26]; if type(_EnvDSNumericSaveDigits, 'posint') then _dat[4][26] := _EnvDSNumericSaveDigits else _dat[4][26] := Digits end if; _Env_dsolve_SC_native := true; if _dat[4][25] = 1 then _i := 1; _dat[4][25] := 2 else _i := _dat[4][25] end if; _val := `dsolve/numeric/SC/IVPval`(_dat, _xout, _src); _dat[4][25] := _i; _dat[4][26] := _dig; [_xout, seq(_val[_vmap[_i]], _i = 1 .. _n-_ne)] else Digits := _dat[4][26]; _val := `dsolve/numeric/SC/IVPval`(eval(_dat, 2), _xout, _src); [_xout, seq(_val[_vmap[_i]], _i = 1 .. _n-_ne)] end if end proc, (2) = Array(0..0, {}), (3) = [t, lambda(t), x(t), diff(x(t), t), y(t), diff(y(t), t)], (4) = []}); _vars := _dat[3]; _pars := map(lhs, _dat[4]); _n := nops(_vars)-1; _solnproc := _dat[1]; if not type(_xout, 'numeric') then if member(x_rkf45_dae, ["start", 'start', "method", 'method', "left", 'left', "right", 'right', "leftdata", "rightdata", "enginedata", "eventstop", 'eventstop', "eventclear", 'eventclear', "eventstatus", 'eventstatus', "eventcount", 'eventcount', "laxtol", 'laxtol', "numfun", 'numfun', NULL]) then _res := _solnproc(convert(x_rkf45_dae, 'string')); if 1 < nops([_res]) then return _res elif type(_res, 'array') then return eval(_res, 1) elif _res <> "procname" then return _res end if elif member(x_rkf45_dae, ["last", 'last', "initial", 'initial', "parameters", 'parameters', "initial_and_parameters", 'initial_and_parameters', NULL]) then _xout := convert(x_rkf45_dae, 'string'); _res := _solnproc(_xout); if _xout = "parameters" then return [seq(_pars[_i] = _res[_i], _i = 1 .. nops(_pars))] elif _xout = "initial_and_parameters" then return [seq(_vars[_i+1] = [_res][1][_i+1], _i = 0 .. _n), seq(_pars[_i] = [_res][2][_i], _i = 1 .. nops(_pars))] else return [seq(_vars[_i+1] = _res[_i+1], _i = 0 .. _n)] end if elif type(_xout, `=`) and member(lhs(_xout), ["initial", 'initial', "parameters", 'parameters', "initial_and_parameters", 'initial_and_parameters', NULL]) then _xout := convert(lhs(x_rkf45_dae), 'string') = rhs(x_rkf45_dae); if type(rhs(_xout), 'list') then _res := _solnproc(_xout) else error "initial and/or parameter values must be specified in a list" end if; if lhs(_xout) = "initial" then return [seq(_vars[_i+1] = _res[_i+1], _i = 0 .. _n)] elif lhs(_xout) = "parameters" then return [seq(_pars[_i] = _res[_i], _i = 1 .. nops(_pars))] else return [seq(_vars[_i+1] = [_res][1][_i+1], _i = 0 .. _n), seq(_pars[_i] = [_res][2][_i], _i = 1 .. nops(_pars))] end if elif type(_xout, `=`) and member(lhs(_xout), ["eventdisable", 'eventdisable', "eventenable", 'eventenable', "eventfired", 'eventfired', "direction", 'direction', NULL]) then return _solnproc(convert(lhs(x_rkf45_dae), 'string') = rhs(x_rkf45_dae)) elif _xout = "solnprocedure" then return eval(_solnproc) elif _xout = "sysvars" then return _vars end if; if procname <> unknown then return ('procname')(x_rkf45_dae) else _ndsol := 1; _ndsol := _ndsol; _ndsol := pointto(_dat[2][0]); return ('_ndsol')(x_rkf45_dae) end if end if; try _res := _solnproc(_xout); [seq(_vars[_i+1] = _res[_i+1], _i = 0 .. _n)] catch: error end try end proc

 (3)
> 

dsol1(1/2);

[t = .500000000000000, lambda(t) = HFloat(4.837512079685444), x(t) = HFloat(0.15920394558985643), diff(x(t), t) = HFloat(0.003966955001051789), y(t) = HFloat(-0.9872457162788973), diff(y(t), t) = HFloat(6.396972656506275e-4)]

 (4)
> 

r:=sqrt(x(t)^2+y(t)^2);#expression to compute

(x(t)^2+y(t)^2)^(1/2)

 (5)
> 

plots:-odeplot(dsol1,[t,r],t=0..10)

 
> 

eval(r,dsol1(0.5));

HFloat(1.0000000003012055)

 (6)
> 

plot(eval(r,dsol1(t)),t=0..10)

Error, invalid input: eval received dsol1(t), which is not valid for its 2nd argument, eqns

  
> 

plot('eval(r,dsol1(t))',0..10)

Error, (in plot) procedure expected, as range contains no plotting variable

  
> 

plot(unapply('eval(r,dsol1(t))',t),0..10);#complicated way

 
> 

r:=t->sqrt(x(t)^2+y(t)^2);#function to compute

proc (t) options operator, arrow; sqrt(x(t)^2+y(t)^2) end proc

 (7)
> 

f:=x-> dsol1(x),'r(t)'

proc (x) options operator, arrow; dsol1(x) end proc, r(t)

 (8)
> 

subs[eval](f(0.5))

HFloat(1.0000000003012055)

 (9)
> 

plot(unapply('subs[eval](f(t))',t),0..10);# equally complicated

 

NULL

Download computation_with_dsolve_exports.mw

Edit: I found these after sending but I still think that there is more to improve. Is eval aor subs needed at all?

plot('eval(r,dsol1(t))',t=0..10)
plot(''eval(r,dsol1(t))'',t=0..10)

Equation of motion of a vibrational elastic plate...

Asked by: alamir 10
differential-equations incomplete-question
November 01 2024
1  10

Hi, I am trying to build a code to get the equation of motion of a vibrational elastic plate in Maple 18. Could anybody help me?

why order of odetest makes difference in how it wo...

Asked by: nm 11353
ode remember differential-equations
October 31 2024
3  8

I think I just found one of the most serious problems in Maple I've seen (other than timelimit hanging).

This is using Maple 2024.2 on windows 10.

I'll explain in words the problem, then give worksheet below to reproduce this. I can reproduce this all the time.

I have implicit solution in y(x) to an ode.

If I first solve for y(x) from the solution, so solution is now explicit, then call odetest to check if this explicit solution is correct, and if I use assumptions on the odetest call, and then after that call odetest on the original implicit solution, then odetest fails to verify the implicit solution.

But, if I change the order, and first call odetest to verify the implicit solution first, it verifies it OK !  

So the problem ONLY happens if I change the order of calling odetest and if I use assumptions on the odetest call that was used before on the explicit solution.

This tells me that Maple remembers something from earlier call. Does it remember the assumptions used? If so, this is very risky. As some part of program might call odetest with some assumptions, and another part of the program can use no assumptions.  I thought assuming is only applied to the call it is used at only and will not affect future calls.

Is there a way then to clear all assumptions used on earlier call to Maple command before using the command again? Or to tell Maple not to remember assumptions used on a call?

This is a big problem. It took me 14 hrs of debuging to find it. Order of calls to odetest should not make it behave different.

I hope someone could find solution to this, since now I have no idea what to do as I need to use odetest on explicit and implicit solutions and I do not want the order of calling Maple command to make difference in results.

This worksheet shows the problem. 3 cases are given. Notice that when using assumptions on earlier call to odetest, how it fails to verify the implicit solution in later call.

> 

restart;
> 

interface(version);

`Standard Worksheet Interface, Maple 2024.2, Windows 10, October 29 2024 Build ID 1872373`

> 

Physics:-Version();

`The "Physics Updates" version in the MapleCloud is 1824 and is the same as the version installed in this computer, created 2024, October 31, 14:22 hours Pacific Time.`

> 

kernelopts('assertlevel'=2):

CASE 1. Calling odetest on explicit first with assumptions, make odetest hang when calling on implicit after

 
> 

restart;

> 

IMPLICIT_SOL:=ln(10*x^2 + (10*y(x) + 12)*x + 5*y(x)^2 + 8*y(x) + 4)/5 - (4*arctan((5*y(x) + 5*x + 4)/(5*x + 2)))/5 = (3*ln(2))/5 + (4*arctan(3))/5;
ode:=2*x+3*y(x)+2+(y(x)-x)*diff(y(x),x) = 0;
IC:=y(0) = -2;

(1/5)*ln(10*x^2+(10*y(x)+12)*x+5*y(x)^2+8*y(x)+4)-(4/5)*arctan((5*y(x)+5*x+4)/(5*x+2)) = (3/5)*ln(2)+(4/5)*arctan(3)

2*x+3*y(x)+2+(y(x)-x)*(diff(y(x), x)) = 0

y(0) = -2

> 

#find explicit solution first
EXPLICIT_SOL:=solve(IMPLICIT_SOL,y(x)):

> 

#NOw check the explicit solution. Using assumptions to see the problem
timelimit(30,  ( odetest(y(x)=EXPLICIT_SOL,[ode,IC]) assuming positive, y(x)::positive) );

Error, (in evalr) time expired

> 

#Now odetest hangs on the implicit solution

> 

timelimit(30,  odetest(IMPLICIT_SOL,[ode,IC]) );

Error, (in is/internal/rename) time expired

> 

 

> 

 

(1/5)*ln(10*x^2+(10*y(x)+12)*x+5*y(x)^2+8*y(x)+4)-(4/5)*arctan((5*y(x)+5*x+4)/(5*x+2)) = (3/5)*ln(2)+(4/5)*arctan(3)

y(x)

2*x+3*y(x)+2+(y(x)-x)*(diff(y(x), x)) = 0

y(0) = -2

CASE 2. Calling odetest on implicit solution first, then it DOES NOT hang !!

 
> 

restart;

> 

IMPLICIT_SOL:=ln(10*x^2 + (10*y(x) + 12)*x + 5*y(x)^2 + 8*y(x) + 4)/5 - (4*arctan((5*y(x) + 5*x + 4)/(5*x + 2)))/5 = (3*ln(2))/5 + (4*arctan(3))/5;
ode:=2*x+3*y(x)+2+(y(x)-x)*diff(y(x),x) = 0;
IC:=y(0) = -2;

(1/5)*ln(10*x^2+(10*y(x)+12)*x+5*y(x)^2+8*y(x)+4)-(4/5)*arctan((5*y(x)+5*x+4)/(5*x+2)) = (3/5)*ln(2)+(4/5)*arctan(3)

2*x+3*y(x)+2+(y(x)-x)*(diff(y(x), x)) = 0

y(0) = -2

> 

#notice, no hang now, since called before
timelimit(30,  odetest(IMPLICIT_SOL,[ode,IC]) );

[0, 0]

> 

#NOw check the explicit solution. This will timeout which is OK
EXPLICIT_SOL:=solve(IMPLICIT_SOL,y(x)):
timelimit(30,  ( odetest(y(x)=EXPLICIT_SOL,[ode,IC]) assuming positive, y(x)::positive) );

Error, (in type/complex) time expired

> 

#check again that odetest still verifies the implicit solution OK

> 

timelimit(30,  odetest(IMPLICIT_SOL,[ode,IC]) );

[0, 0]

> 

 

CASE 3. Calling odetest with no assumptions on explicit solution, then it also does not hang

 
> 

restart;

> 

IMPLICIT_SOL:=ln(10*x^2 + (10*y(x) + 12)*x + 5*y(x)^2 + 8*y(x) + 4)/5 - (4*arctan((5*y(x) + 5*x + 4)/(5*x + 2)))/5 = (3*ln(2))/5 + (4*arctan(3))/5;
ode:=2*x+3*y(x)+2+(y(x)-x)*diff(y(x),x) = 0;
IC:=y(0) = -2;

(1/5)*ln(10*x^2+(10*y(x)+12)*x+5*y(x)^2+8*y(x)+4)-(4/5)*arctan((5*y(x)+5*x+4)/(5*x+2)) = (3/5)*ln(2)+(4/5)*arctan(3)

2*x+3*y(x)+2+(y(x)-x)*(diff(y(x), x)) = 0

y(0) = -2

> 

#Now check the explicit solution. but DO NOT use assumptions
EXPLICIT_SOL:=solve(IMPLICIT_SOL,y(x)):
timelimit(30, odetest(y(x)=EXPLICIT_SOL,[ode,IC]));

Error, (in factor/remember) time expired

> 

#check again that odetest still verifies the implicit solution OK

> 

timelimit(30,  odetest(IMPLICIT_SOL,[ode,IC]) );

[0, 0]

> 

 


 

Download order_of_ode_test_makes_difference_oct_31_2024.mw

 

update NOV 2, 2024 6 AM

I found the cause.

Removing PHYSICS from libname, then the problem goes away !

So this is caused by PHYSICS package. For some reason, having Physics package in the libname causes odetest to hang/fail compared to when the physics package is not in the libname path. Worksheet below.

> 

restart;
> 

interface(version);

`Standard Worksheet Interface, Maple 2024.2, Windows 10, October 29 2024 Build ID 1872373`

> 

Physics:-Version();

`The "Physics Updates" version in the MapleCloud is 1824 and is the same as the version installed in this computer, created 2024, October 31, 14:22 hours Pacific Time.`

> 

libname;

"C:\Users\Owner\maple\toolbox\2024\Physics Updates\lib", "C:\Program Files\Maple 2024\lib"

CASE 1. without PHYSICS on libname, it works !!

 
> 

restart;

> 

kernelopts('assertlevel'=2):

> 

libname;

"C:\Users\Owner\maple\toolbox\2024\Physics Updates\lib", "C:\Program Files\Maple 2024\lib"

> 

libname:=libname[2]; #remove PHYSICS

"C:\Program Files\Maple 2024\lib"

> 

libname;

"C:\Program Files\Maple 2024\lib"

> 

IMPLICIT_SOL:=ln(10*x^2 + (10*y(x) + 12)*x + 5*y(x)^2 + 8*y(x) + 4)/5 - (4*arctan((5*y(x) + 5*x + 4)/(5*x + 2)))/5 = (3*ln(2))/5 + (4*arctan(3))/5;
ode:=2*x+3*y(x)+2+(y(x)-x)*diff(y(x),x) = 0;
IC:=y(0) = -2;

(1/5)*ln(10*x^2+(10*y(x)+12)*x+5*y(x)^2+8*y(x)+4)-(4/5)*arctan((5*y(x)+5*x+4)/(5*x+2)) = (3/5)*ln(2)+(4/5)*arctan(3)

2*x+3*y(x)+2+(y(x)-x)*(diff(y(x), x)) = 0

y(0) = -2

> 

#find explicit solution first
EXPLICIT_SOL:=solve(IMPLICIT_SOL,y(x)):

> 

#NOw check the explicit solution. Using assumptions to see the problem
timelimit(30,  ( odetest(y(x)=EXPLICIT_SOL,[ode,IC]) assuming positive, y(x)::positive) );

Error, (in evalr/ln) time expired

> 

#Now try odetest  on the implicit solution

> 

infolevel[odetest]:=5;

5

> 

timelimit(30,  odetest(IMPLICIT_SOL,[ode,IC]) );

odetest: Performing an implicit solution test

odetest: Performing an explicit (try soft) solution test

odetest: Performing an implicit solution (II) test

[0, 0]

> 

 

> 

 

CASE 2. With Physics on libname, it fails !

 
> 

restart;

> 

kernelopts('assertlevel'=2):

> 

libname;

"C:\Users\Owner\maple\toolbox\2024\Physics Updates\lib", "C:\Program Files\Maple 2024\lib"

> 

IMPLICIT_SOL:=ln(10*x^2 + (10*y(x) + 12)*x + 5*y(x)^2 + 8*y(x) + 4)/5 - (4*arctan((5*y(x) + 5*x + 4)/(5*x + 2)))/5 = (3*ln(2))/5 + (4*arctan(3))/5;
ode:=2*x+3*y(x)+2+(y(x)-x)*diff(y(x),x) = 0;
IC:=y(0) = -2;

(1/5)*ln(10*x^2+(10*y(x)+12)*x+5*y(x)^2+8*y(x)+4)-(4/5)*arctan((5*y(x)+5*x+4)/(5*x+2)) = (3/5)*ln(2)+(4/5)*arctan(3)

2*x+3*y(x)+2+(y(x)-x)*(diff(y(x), x)) = 0

y(0) = -2

> 

#find explicit solution first
EXPLICIT_SOL:=solve(IMPLICIT_SOL,y(x)):

> 

timelimit(30,  ( odetest(y(x)=EXPLICIT_SOL,[ode,IC]) assuming positive, y(x)::positive) );

Error, (in evalr/shake) time expired

> 

#Now try odetest  on the implicit solution

> 

infolevel[odetest]:=5;

5

> 

timelimit(30,  odetest(IMPLICIT_SOL,[ode,IC]) );

odetest: Performing an implicit solution test

odetest: Performing an explicit (try soft) solution test

odetest: Performing an implicit solution (II) test

Error, (in is/internal/rename) time expired

 

 

Download order_of_ode_test_makes_difference_NOV_2_2024.mw

my question is: Is it safe to permanently remove Physics package from libname? Why is having physics package in libname (which is by default) causes this problem?

I do not use Physics package explicitly in my code. i.e. I do not do Physics:- calls. 

Any one knows what effect not having Physics on libname cause? Will Maple still work OK for everything if one is not calling Physics package explicitly?

Content of exported SVG files is always outside th...

Asked by: mafitch 10
export svg
October 31 2024
2  0

When I export images as SVG the resulting file always has the image much larger than the viewbox resulting in only part of it showing. My normal workflow is to then load in Inkscape and correct the error, however, I would love to not have an extra step. Is there any other fix for this?

How do you convert this system of equations into m...

Asked by: GFY 65
matrix ode differential-equations
October 31 2024
0  3

How do you convert this system of equations into matrix form? The decoupling process is performed.Convert the equations into photographs.fuxian1030.mw

Assign variables to curve fit parameters...

Asked by: priemnotas 5
parameters curve-fitting incomplete-question
October 31 2024
1  1

Low power user (occasional) here! I would like to extract the parameters from a curve fit and assign to another variable?

Change of coordinates: Are there other ways to loo...

Asked by: C_R 3412
coords
October 31 2024
1  2

What I do at the moment

u-> changecoords(u,[x,y,z],spherical_physics,[r,theta,phi]):
map(%,[x,y,z])
  [r cos(phi) sin(theta), r sin(phi) sin(theta), r cos(theta)]

Any other (preferably shorter) ways to look up transformations defined in ?coords

What is the probability that the total value of th...

Asked by: yangtheary 280
probability
October 30 2024
1  3

How do i define a tensor with an equation of three...

Asked by: Steve255 15
tensor differentialgeometry physics
October 30 2024
3  4

I am studying linear operators that have vanishing Nijenhuis Torsion, that is a (1-1) tensor L whose corresponding (1-2) tensor N given by a tensor equation of L (see attached) is identically zero.

I am new to maple, i have used it to plot vector fields and solve systems of equations in the past but i am unfamiliar with the DifferentialGeometry and Physics packages.

Attached is my best effort at solving this problem directly for a simple 2d case without the use of any packages and i am wondering if it is possible to do it all in one line without having to define tensor components one by one.

nijenhuis_torsion.mw

First 57 58 59 60 61 62 63 Last Page 59 of 2425