Items tagged with ode ode Tagged Items Feed

I have an ODE with L(x,y), which includes partial derivatives of L(x,y) and coefficients y,_y1 and _y1^2. I want to split this equation into seperate equations of these coefficients. So I can solve for the unknown variable L(x,y).

 

ODE:=diff(L(x,y),x,x)+2*_y1(x)*diff(L(x,y),x,y)+_y1(x)^2*diff(L(x,y),y,y)+y(x)/x*diff(L(x,y),x)+(y(x)*_y1(x)/x+_y1(x))*diff(L(x,y),y)-y(x)/x^2*L(x,y)=0;

Could someone help me? 

 

I am currently trying to use the coeffs command but not really getting anywhere.

 

Thank you.

I'm trying to find lypunov exponent for this  system of ODEs. I know I need to take the Jacobian but not sure if it's possible the way it's currently defined. If anyone could provide insight it would be much appreciated.

Eqns:= diff(omega(t),t)=-(G*MSat*beta^(2)*(xH(t)*sin(theta(t))-yH(t)* cos(theta(t)))*(xH(t)*cos(theta(t))+yH(t)*sin(theta(t))))/((xH(t)^(2)+yH(t)^(2))^(2.5)),diff(theta(t),t)=omega(t), diff(xH(t),t)=vxH(t),diff(vxH(t),t)=-(G*M*xH(t))/((xH(t)^(2)+yH(t)^(2))^(1.5)),diff(yH(t),t)=vyH(t),diff(vyH(t),t)=-(G*M*yH(t))/((xH(t)^(2)+yH(t)^(2))^(1.5)): ;

ICs := omega(0) = omega0, theta(0) = theta0, xH(0) = a*(1+e), yH(0) = 0, vxH(0) = 0, vyH(0) = sqrt(G*M*(1-e)/(a*(1+e)));

I want the exponent for omega, I  procedure that takes some initial conditions, changing just w0, and computes the long term value of omega. This plots a sort of bifurcation diagram. I'd like an estimate of the exponent to compare what I see. 

thanks for the help

ft

I don't know what is the problem with this. Please take a look at it. Any help is appreciated. Thanks.

restart

with(RealDomain)

[Im, Re, `^`, arccos, arccosh, arccot, arccoth, arccsc, arccsch, arcsec, arcsech, arcsin, arcsinh, arctan, arctanh, cos, cosh, cot, coth, csc, csch, eval, exp, expand, limit, ln, log, sec, sech, signum, simplify, sin, sinh, solve, sqrt, surd, tan, tanh]

(1)

with(Slode)

[DEdetermine, FPseries, FTseries, candidate_mpoints, candidate_points, dAlembertian_formal_sol, dAlembertian_series_sol, hypergeom_formal_sol, hypergeom_series_sol, mhypergeom_formal_sol, mhypergeom_series_sol, msparse_series_sol, polynomial_series_sol, rational_series_sol, series_by_leastsquare]

(2)

ODE := x*(diff(y(x), x, x))-(diff(y(x), x))-x*y(x)

x*(diff(diff(y(x), x), x))-(diff(y(x), x))-x*y(x)

(3)

SS := convert(sin(x), FormalPowerSeries)

Sum((-1)^k*x^(2*k+1)/factorial(2*k+1), k = 0 .. infinity)

(4)

FPseries(ODE = SS, y(x), v(n))

Error, the right-hand side of the differential equation must be a polynomial or a series in x

 

``

 

Download AA.mw

Dear all;

Good morning everyone.

I solve a simple ode, i want how can I write this program as procedure with output the two coefficient involved in the solution, after solving with ics. i.e coef:=dsolve(ode);  in my example I want as output [4, -1] .

 

restart:
ode:=diff(y(x),x,x)=3*diff(y(x),x)-2*y(x);
coef:=dsolve(ode);
ics:=y(0)=3, D(y)(0)=2;
dsolve({ics,ode});

with best regards

 

Hello,

I have two regimes. Each regime is characterized by system of state differential equations (diff(S(t), t), diff(K(t), t)) and co-state differential equations (diff(psi[S](t), t), diff(psi[Iota](t), t)) as follows: 

(1)

diff(S(t), t) = -eta*K(t)*S(t)/(w*N*(S(t)+K(t))), diff(K(t), t) = eta*K(t)*S(t)/(w*N*(S(t)+K(t)))-upsilon

diff(psi[S](t), t) = eta*K(t)^2*(psi[S](t)-psi[Iota](t))/(w*N*(S(t)+K(t))^2), diff(psi[Iota](t), t) = eta*S(t)^2*(psi[S](t)-psi[Iota](t))/(w*N*(S(t)+K(t))^2)

(2)

diff(S(t), t) = -eta*K(t)*S(t)/(w*N*(S(t)+K(t))), diff(K(t), t) = eta*K(t)*S(t)/(w*N*(S(t)+K(t)))

diff(psi[S](t), t) = eta*K(t)^2*(psi[S](t)-psi[Iota](t))/(w*N*(S(t)+K(t))^2), diff(psi[Iota](t), t) = eta*S(t)^2*(psi[S](t)-psi[Iota](t))/(w*N*(S(t)+K(t))^2)

The first regime is employed from 0 to t1 (where t1 is unknown) and then regime 2, from t1 to T. I know the initial values for the state variables of the first system at t=0, that is, S(0)=S0 and K(0)=K0, as well as boundary conditions for the co-state variables for regime 2 at t=T, that is, psi_S(T)=a and psi_I(T)=b. 

I also know that my unknown t1 should satisfy the algebraic equation: -c - psi_S(t1)=0, where psi_S(t1) is the solution of the co-state diff equation of the regime 2 at t=t1. 

My question is: If I assume that the systems have not analytical solution, how can I found unknown t1 numerically? Moreover, asssume that all other parameters, such as T, eta, upsilon and others are given.

Hello,

My code records the values I need, however, I need to implement a modulo of 2*Pi on my result for theta. But this leads to a graph with no plots and I'm not sure how to fix it. Any help is greatly aprreciated! Thank you in advance!

Kind regards,

Gam

with(plots):

a := 1.501*10^9:

Th := sqrt(4*Pi^2*a^3/(G*(Mh+Msat))):

HyperionOrbit := proc (`θIC`, `ωIC`, n) local a, Mh, Msat, G, e, beta, M, Eqns, ICs, soln; option remember; global `ωH`, Th; a := 1.501*10^9; Mh := 5.5855*10^18; Msat := 5.6832*10^26; G := 6.67259/10^11; e := .232; beta := .89; M := Mh+Msat; Eqns := diff(theta(t), t) = omega(t), diff(omega(t), t) = -G*Msat*beta^2*(xH(t)*sin(theta(t))-yH(t)*cos(theta(t)))*(xH(t)*cos(theta(t))+yH(t)*sin(theta(t)))/(xH(t)^2+yH(t)^2)^2.5, diff(xH(t), t) = vxH(t), diff(vxH(t), t) = -G*M*xH(t)/(xH(t)^2+yH(t)^2)^(3/2), diff(yH(t), t) = vyH(t), diff(vyH(t), t) = -G*M*yH(t)/(xH(t)^2+yH(t)^2)^(3/2); ICs := xH(0) = a*(1+e), yH(0) = 0, vxH(0) = 0, vyH(0) = sqrt(G*M*(1-e)/(a*(1+e))), theta(0) = `θIC`, omega(0) = `ωIC`; soln := dsolve({Eqns, ICs}, numeric, maxfun = 0, output = array([seq(i, i = 0 .. n*Th, Th)])); plots:-odeplot(soln, [modp(theta(t), 2*Pi), omega(t)/`ωH`], 0 .. n*Th, labels = ["θ(t)","ω(t)/ωH"], axes = boxed, style = plottools:-point, size = [.25, .75]) end proc:

plots:-display(HyperionOrbit(.5, 1.8*`ωH`, 10));

Download Poincare_section_Boyd_plot_fixing_theta.mw

bia Man

Hello,

I have a question about poincare sections. I have this piece of code i need to analyse and I want to use a poincare section in order to so. How could I do it? I am interested in theta and omega. Any help is greatly appreciated! Thank you in advance!

Kind regards,

Gambia Man

with(plots):

a := 1.501*10^9:

Th := sqrt(4*Pi^2*a^3/(G*(Mh+Msat)));

1876321.326

 

0.3348672330e-5

(1)

HyperionOrbit := proc (`θIC`, `ωIC`) local a, Mh, Msat, G, e, beta, M, Eqns, ICs; global `ωH`, Th, soln; a := 1.501*10^9; Mh := 5.5855*10^18; Msat := 5.6832*10^26; G := 6.67259/10^11; e := .232; beta := .89; M := Mh+Msat; Eqns := diff(theta(t), t) = omega(t), diff(omega(t), t) = -G*Msat*beta^2*(xH(t)*sin(theta(t))-yH(t)*cos(theta(t)))*(xH(t)*cos(theta(t))+yH(t)*sin(theta(t)))/(xH(t)^2+yH(t)^2)^2.5, diff(xH(t), t) = vxH(t), diff(vxH(t), t) = -G*M*xH(t)/(xH(t)^2+yH(t)^2)^(3/2), diff(yH(t), t) = vyH(t), diff(vyH(t), t) = -G*M*yH(t)/(xH(t)^2+yH(t)^2)^(3/2); ICs := xH(0) = a*(1+e), yH(0) = 0, vxH(0) = 0, vyH(0) = sqrt(G*M*(1-e)/(a*(1+e))), theta(0) = `θIC`, omega(0) = `ωIC`; soln := dsolve({Eqns, ICs}, numeric); odeplot(soln, [theta(t), omega(t)/`ωH`], 0 .. 5*Th, numpoints = 2000, labels = ["θ(t)","ω(t)/ωH"], axes = boxed, size = [.25, .75]) end proc

``

 

Download New_Poincare_section.mw

Hi everyone, I'm working a problem

Given a dog and a man. At t=0, the position of the man and his dog is (0,0) and (0,h), respectively. Then, the man start moving along Ox with constant speed vm. The dog keep running toward its master with constant speed vd. Describe the movement of the dog?

Therefore, let say the position of the dog is (f1(x),f2(x)), I wrote these:

restart;

with(plots);

MD := proc (h, vm, vd) local x, y, ode, ics, sol;

ode := {(diff(f2(t), t))*(vm*t-f1(t))+(diff(f1(t), t))*f2(t), diff(f2(t), t) = -sqrt(vd^2-(diff(f1(t), t))^2)};

ics := {f1(0) = h, f2(0) = 0};

sol := dsolve(ode union ics);

x := unapply(eval(f1(t), sol), t);

y := unapply(eval(f2(t), sol), t);

plot([x(t), y(t), t = 0 .. 10], scaling = constrained) end proc;

MD(10, 1, 5);

However, "sol" is returned NULL, which means the equation has no solution. I supposed I completed the motion part of the problem correctly. Please help me or point out an another method

 

Hello,

I have an error in my code. I don't know where it came from. Earlier today I loaded this and it worked fine and now an error comes up. Any help is greatly appreciated! Thank you in advance!

Kind regards,

Gambia Man

with(plots):

``

HyperionOrbit := proc (`θIC`, `ωIC`) local a, Mh, Msat, G, e, beta, M, Eqns, Th, ICs, soln, `ωH`; a := 1.501*10^9; Mh := 5.5855*10^18; Msat := 5.6832*10^26; G := 6.67259/10^11; e := .232; beta := .89; M := Mh+Msat; Eqns := diff(theta(t), t) = omega(t), diff(omega(t), t) = -G*Msat*beta^2*(xH(t)*sin(theta(t))-yH(t)*cos(theta(t)))*(xH(t)*cos(theta(t))+yH(t)*sin(theta(t)))/(xH(t)^2+yH(t)^2)^2.5, diff(xH(t), t) = vxH(t), diff(vxH(t), t) = -G*M*xH(t)/(xH(t)^2+yH(t)^2)^(3/2), diff(yH(t), t) = vyH(t), diff(vyH(t), t) = -G*M*yH(t)/(xH(t)^2+yH(t)^2)^(3/2); Th := sqrt(4*Pi^2*a^3/(G*(Mh+Msat))); `ωH` := 2*Pi/Th; ICs := xH(0) = a*(1+e), yH(0) = 0, vxH(0) = 0, vyH(0) = sqrt(G*M*(1-e)/(a*(1+e))), theta(0) = `θIC`, omega(0) = `ωIC`; soln := dsolve({Eqns, ICs}, numeric); odeplot(soln, [theta(t), omega(t)/`ωH`], 0 .. 5*Th, numpoints = 2000, labels = ["θ(t)","ω(t)/ωH"], axes = boxed) end proc

HyperionOrbit(.5, 1.8*`ωH`)

Warning, The use of global variables in numerical ODE problems is deprecated, and will be removed in a future release. Use the 'parameters' argument instead (see ?dsolve,numeric,parameters)

 

Error, (in dsolve/numeric/SC/IVPpoints) unable to store '11.3097335529232/Th' when datatype=float[8]

 

``

``

``

 

Download Poincare_section.mw

Hello,

I have an issue with this code. The issue lies in the warning which is outputed. However, I don't know how to resolve it. Any help would be greatly appreciated! Thank you in advance!

Kind regards,

Gambia Man

with(plots):

a := 1.501*10^9:

Eqns := diff(xH(t), t) = vxH(t), diff(vxH(t), t) = -G*(Mh+Msat)*xH(t)/(xH(t)^2+yH(t)^2)^(3/2), diff(yH(t), t) = vyH(t), diff(vyH, t) = -G*(Mh+Msat)*yH(t)/(xH(t)^2+yH(t)^2)^(3/2):

ICs := xH(0) = a*(1+e), yH(0) = 0, vxH(0) = 0, vyH(0) = sqrt(G*(Mh+Msat)*(1-e)/(a*(1+e))):

T := sqrt(4*Pi^2*a^3/(G*(Mh+Msat))):

soln := dsolve({Eqns, ICs}, numeric):

plots:-odeplot(soln, [xH(t)/a, yH(t)/a], 0 .. 10*T, scaling = constrained, labels = ["x/a", "y/a"], numpoints = 2000);

Warning, cannot evaluate the solution further right of 453574.15, probably a singularity

 

 

``

 

Download Hyperion_Orbit.mw

 

Suppose we have a differential equation that it's order is 2. For example

diff(x(t), t, t)+(1000*(x(t)^2-1))*(diff(x(t), t))+x(t) = 0,

and we want to solve it numerically. When we use some initial values and dsolve it, the maple solves it very quickly. We can plot each new equations of x(t) without any problem respect of "t" or else.

But when we want to show the variation of "diff(x(t), t, t)" respect of "t" (when our derivation order is equal or upper than our differential equation order) our answer is 0 however dx(t)/dt is not 0 or constant!

What's happen here?

I tested it with all kind of methods as "rkf45","rosenbrock","lsode" etc.  but non of them showed me a correct answer.

How can I solve it correctly?
 
  
Please help me!

 

This is the code hw2_final.mw

This is the warning 

Warning, The use of global variables in numerical ODE problems is deprecated, and will be removed in a future release. Use the 'parameters' argument instead (see ?dsolve,numeric,parameters )"

 

How to solve it?

Hello,

 I am trying to solve an differential equation in Maple, but I have no response from Maple.Can someone give a look for me? The coefficients a and b cannot be zero.

 Thanks a lot.

 

Download edo_solução_A_U.mw

Hello everybody

I have a system of ODEs that are too large. I am trying to solve this system using Maple. The corrsponding file has been uploaded. I have 2 question. 1- Does anyone know what is the source of error? 2- Is this a high expectation for Maple to solve such a too big system of ODEs. Note: It makes a few minutes to run all the programm.

Thanks in advance

M.mw

Hello

i have an ODE like this:

I sove this ODE with plot order:

with(plots);
odeplot(sol, [x, (3*D1*a+4*D2)*P(x)/((1-q*S(x))*D2)], .5 .. (1/2)*Pi, tickmarks = [[seq((1/10)*i*Pi = (180*i*(1/10))*`°`, i = 1 .. 8)], default]);
my plot work very well. but i need to plot this ODE with five different parameter (q for for instance, q=0.1 & q=0.2 ....) all in one axis. something like this:

2 3 4 5 6 7 8 Last Page 4 of 28