Items tagged with numeric numeric Tagged Items Feed

Hi everyone!

I have a problem solving the nonlinear ode (as attached below). I got this error ---> Error, (in fproc) unable to store '-1.32352941215398+(-0.441176470717993e-1, -0.)' when datatype=float[8]

1) Could someone please explain to me what does the unable store .... error means? 

and i will be grateful if you could help me finding the solution out. Thanks in advance


problem4.mw

I am working in Document mode, 2D input.  I have doing some complex number evalutions

Lets say I have a complex number Z= 2+3i

when I ask maple for the "abs(Z)" I get the display of "abs(2+3i)" rather than the numberic answer. This occurs with all the complex operation arg, abs, polar conversion, 

polar(abs(2+3*i), argument(2+3*i))   this is from Maple. 

How do I display the numeric value for these functions in my document

 

Thanks Bill

HI.please help me for solve differenrtial equation with finite difference method not dsolve numeric solver in maple

thanks

FDM2.mw

L := 1/50000000; -1; eta := 1; -1; PDE[111] := 7.65692307692309*10^(-8)*(diff(f1(x), x, x, x, x))-1.56784615384616*10^12*(diff(f1(x), x, x))+220.592307692308*(diff(f2(x), x, x, x))-3.52947692307693*10^21*(diff(f2(x), x))+43.7538461538462*(diff(f3(x), x, x, x))+4.81292307692309*10^20*(diff(f3(x), x))+6.50473846153848*10^30*f1(x)-7.90000000000000*10^(-8)*eta*f1(x)

0.7656923077e-7*(diff(diff(diff(diff(f1(x), x), x), x), x))-0.1567846154e13*(diff(diff(f1(x), x), x))+220.592307692308*(diff(diff(diff(f2(x), x), x), x))-0.3529476923e22*(diff(f2(x), x))+43.7538461538462*(diff(diff(diff(f3(x), x), x), x))+0.4812923077e21*(diff(f3(x), x))+0.6504738462e31*f1(x)

(1)

PDE[222] := 2.14211538461539*10^(-8)*(diff(f2(x), x, x, x, x))-1.64988461538462*10^12*(diff(f2(x), x, x))+7.90486153846156*10^30*f2(x)-220.592307692308*(diff(f1(x), x, x, x))+3.52947692307693*10^21*(diff(f1(x), x))-5.94323076923080*10^11*(diff(f3(x), x, x))+5.13378461538463*10^30*f3(x)-7.90000000000000*10^(-8)*eta*f2(x)

0.2142115385e-7*(diff(diff(diff(diff(f2(x), x), x), x), x))-0.1649884615e13*(diff(diff(f2(x), x), x))+0.7904861538e31*f2(x)-220.592307692308*(diff(diff(diff(f1(x), x), x), x))+0.3529476923e22*(diff(f1(x), x))-0.5943230769e12*(diff(diff(f3(x), x), x))+0.5133784615e31*f3(x)

(2)

PDE[333] := -6.38076923076924*10^(-31)*(diff(f3(x), x, x, x, x, x, x))+9.66537046153848*10^(-8)*(diff(f3(x), x, x, x, x))-3.10154753538461*10^12*(diff(f3(x), x, x))-43.7538461538462*(diff(f1(x), x, x, x))-4.81292307692309*10^20*(diff(f1(x), x))-5.94323076923080*10^11*(diff(f2(x), x, x))+5.13378461538463*10^30*f2(x)+2.29989058707693*10^31*f3(x)-7.90105333333333*10^(-8)*omega^2*f3(x)+6.58333333333333*10^(-31)*eta*(diff(f3(x), x, x))

-0.6380769231e-30*(diff(diff(diff(diff(diff(diff(f3(x), x), x), x), x), x), x))+0.9665370462e-7*(diff(diff(diff(diff(f3(x), x), x), x), x))-0.3101547535e13*(diff(diff(f3(x), x), x))-43.7538461538462*(diff(diff(diff(f1(x), x), x), x))-0.4812923077e21*(diff(f1(x), x))-0.5943230769e12*(diff(diff(f2(x), x), x))+0.5133784615e31*f2(x)+0.2299890587e32*f3(x)-0.7901053333e-7*omega^2*f3(x)

(3)

bcs := {f1(0) = 0, f1(L) = 0, f2(0) = 0, f2(L) = 0, f3(0) = 0, f3(L) = 0, ((D@@1)(f1))(0) = 0, ((D@@1)(f1))(L) = 0, ((D@@1)(f2))(0) = 0, ((D@@1)(f2))(L) = 0, ((D@@1)(f3))(0) = 0, ((D@@1)(f3))(L) = 0, ((D@@2)(f3))(0) = 0, ((D@@2)(f3))(L) = 0}

{f1(0) = 0, f1(1/50000000) = 0, f2(0) = 0, f2(1/50000000) = 0, f3(0) = 0, f3(1/50000000) = 0, (D(f1))(0) = 0, (D(f1))(1/50000000) = 0, (D(f2))(0) = 0, (D(f2))(1/50000000) = 0, (D(f3))(0) = 0, (D(f3))(1/50000000) = 0, ((D@@2)(f3))(0) = 0, ((D@@2)(f3))(1/50000000) = 0}

(4)

``


Download FDM2.mw

Also always i try to solve that strange ecuation maple give me like an aswer in terms of Z, with an integral on it, i don't get nothing. I evaluate a point with right click, evaluate a point. Choose the variables value, and the time at wich i want to evaluate the ecuation. Later it give me the same strange solution but with the variables replaced for the values i gave to the program. Then i choose right click, aproximate, and don't matter if i select 5, 10 or more digits, it launchs to me the same strange ecuation, not the number i'm waiting for.

So that i want to know now is how can i manipulate a numerical solution to a diff ecuation.

For example if i got this system of diff ecuation

dx/dt = (-k/m) (x srqt(x^2 + y^2))

dy/dt = (-k/m) (y srqt(x^2 + y^2)) - g

Which is a projectile trajectory taking care of air drag force, and i want to know at which initial speed i have to fire the body to reach 300m if it is fired from the initial point at (0,0), and it's initial velocity on x is 3 times it's initial velocity on y (initial launch angle condition)? 

How can i input that conditions to make maple give me a solution for something like that?

 

Hi

Please find the attachment

It seems a singularity has been occurred at left end point for n less than unit (say n=0.5, n>0)

Is there a way to fix it?

 

n.mw

DEAR SIR,

PLEASE HELP ME WITH THAT QUESTION

 

DEAR SIR

ANYONE CAN HELP TO COMPUTE TIME IN DSOLVE COMMAND?

restart

with(plots)

Nb := 0.1e-4; Nt := 0.1e-4; Sc := 3.0; Sb := 15; Pe := 1; Bi := .5; Pr := 6.8; c[4] := 0; c[6] := .3; c[8] := .4; k[1] := 0; k[2] := 1; k[3] := 0; Un := .1; M := .5

Eq1 := (101-100*lambda)*(1+c[2]*phi(eta))*(diff(f(eta), `$`(eta, 3)))+(diff(f(eta), `$`(eta, 2)))*(f(eta)+g(eta)+c[2]*(diff(phi(eta), eta)))-(diff(f(eta), eta))^2-M*(diff(f(eta), eta)-k[2])+k[2]+Un*(k[2]-(1/2)*eta*(diff(f(eta), `$`(eta, 2)))-(diff(f(eta), eta))) = 0

(101-100*lambda)*(1+c[2]*phi(eta))*(diff(diff(diff(f(eta), eta), eta), eta))+(diff(diff(f(eta), eta), eta))*(f(eta)+g(eta)+c[2]*(diff(phi(eta), eta)))-(diff(f(eta), eta))^2-.6*(diff(f(eta), eta))+1.6-0.5000000000e-1*eta*(diff(diff(f(eta), eta), eta)) = 0

(1)

Eq2 := (101-100*lambda)*(1+c[2]*phi(eta))*(diff(g(eta), `$`(eta, 3)))+(diff(g(eta), `$`(eta, 2)))*(f(eta)+g(eta)+c[2]*(diff(phi(eta), eta)))-(diff(g(eta), eta))^2-M*(diff(g(eta), eta)-k[2])+k[2]+Un*(k[2]-(1/2)*eta*(diff(g(eta), `$`(eta, 2)))-(diff(g(eta), eta))) = 0

(101-100*lambda)*(1+c[2]*phi(eta))*(diff(diff(diff(g(eta), eta), eta), eta))+(diff(diff(g(eta), eta), eta))*(f(eta)+g(eta)+c[2]*(diff(phi(eta), eta)))-(diff(g(eta), eta))^2-.6*(diff(g(eta), eta))+1.6-0.5000000000e-1*eta*(diff(diff(g(eta), eta), eta)) = 0

(2)

Eq3 := (1+c[4]*phi(eta))*(diff(theta(eta), `$`(eta, 2)))+Pr*(diff(theta(eta), eta))*(f(eta)+g(eta))+Nb*Pr*(diff(theta(eta), eta))*(diff(phi(eta), eta))*(1+c[6]*(2*phi(eta)+1))+Nt*Pr*(diff(theta(eta), eta))^2+c[4]*(diff(theta(eta), eta))*(diff(phi(eta), eta))-(1/2)*Pr*eta*Un*(diff(theta(eta), eta)) = 0

diff(diff(theta(eta), eta), eta)+6.8*(diff(theta(eta), eta))*(f(eta)+g(eta))+0.68e-4*(diff(theta(eta), eta))*(diff(phi(eta), eta))*(1.3+.6*phi(eta))+0.68e-4*(diff(theta(eta), eta))^2-.3400000000*eta*(diff(theta(eta), eta)) = 0

(3)

Eq4 := (1+c[6]*phi(eta))*(diff(phi(eta), `$`(eta, 2)))+Sc*(f(eta)+g(eta))*(diff(phi(eta), eta))+c[6]*(diff(phi(eta), eta))^2+Nt*(diff(theta(eta), `$`(eta, 2)))/Nb-(1/2)*Sc*eta*Un*(diff(phi(eta), eta)) = 0

(1+.3*phi(eta))*(diff(diff(phi(eta), eta), eta))+3.0*(f(eta)+g(eta))*(diff(phi(eta), eta))+.3*(diff(phi(eta), eta))^2+1.000000000*(diff(diff(theta(eta), eta), eta))-.1500000000*eta*(diff(phi(eta), eta)) = 0

(4)

Eq5 := (1+c[8]*phi(eta))*(diff(chi(eta), `$`(eta, 2)))+c[8]*(diff(phi(eta), eta))*(diff(chi(eta), eta))-Pe*(chi(eta)*(diff(phi(eta), `$`(eta, 2)))+(diff(phi(eta), eta))*(diff(chi(eta), eta)))+Sb*(diff(chi(eta), eta))*(f(eta)+g(eta))-(1/2)*Sb*eta*Un*(diff(chi(eta), eta)) = 0

(1+.4*phi(eta))*(diff(diff(chi(eta), eta), eta))-.6*(diff(phi(eta), eta))*(diff(chi(eta), eta))-chi(eta)*(diff(diff(phi(eta), eta), eta))+15*(diff(chi(eta), eta))*(f(eta)+g(eta))-.7500000000*eta*(diff(chi(eta), eta)) = 0

(5)

Vc[2] := [.2, .3, .4]

etainf := 1.85

bcs := (D(f))(0) = k[1], (D(g))(0) = k[3], f(0) = 0, g(0) = 0, (D(theta))(0) = -Bi*(1-theta(0))/(1+c[4]*phi(0)), Nb*(D(phi))(0)*(1+c[6]*(2*phi(0)+1))+Nt*(D(theta))(0) = 0, chi(0) = 1, (D(f))(etainf) = k[2], (D(g))(etainf) = k[2], theta(etainf) = 0, phi(etainf) = 0, chi(etainf) = 0

(D(f))(0) = 0, (D(g))(0) = 0, f(0) = 0, g(0) = 0, (D(theta))(0) = -.5+.5*theta(0), 0.1e-4*(D(phi))(0)*(1.3+.6*phi(0))+0.1e-4*(D(theta))(0) = 0, chi(0) = 1, (D(f))(1.85) = 1, (D(g))(1.85) = 1, theta(1.85) = 0, phi(1.85) = 0, chi(1.85) = 0

(6)

dsys := {Eq1, Eq2, Eq3, Eq4, Eq5, bcs}

for i to 3 do c[2] := Vc[2][i]; dsol[i] := dsolve(dsys, numeric, continuation = lambda); print(c[2]); print(dsol[i](0)) end do

.2

 

[eta = 0., chi(eta) = HFloat(1.0), diff(chi(eta), eta) = HFloat(3.888290578689045), f(eta) = HFloat(0.0), diff(f(eta), eta) = HFloat(0.0), diff(diff(f(eta), eta), eta) = HFloat(2.244199282192492), g(eta) = HFloat(0.0), diff(g(eta), eta) = HFloat(0.0), diff(diff(g(eta), eta), eta) = HFloat(2.244199282192492), phi(eta) = HFloat(-2.044191234673432), diff(phi(eta), eta) = HFloat(5.227515304629519), theta(eta) = HFloat(0.2317093657771352), diff(theta(eta), eta) = HFloat(-0.38414531711143246)]

 

.3

 

[eta = 0., chi(eta) = HFloat(1.0), diff(chi(eta), eta) = HFloat(4.148187853914835), f(eta) = HFloat(0.0), diff(f(eta), eta) = HFloat(0.0), diff(diff(f(eta), eta), eta) = HFloat(3.1086884419918364), g(eta) = HFloat(0.0), diff(g(eta), eta) = HFloat(0.0), diff(diff(g(eta), eta), eta) = HFloat(3.1086884419918364), phi(eta) = HFloat(-2.049332060722701), diff(phi(eta), eta) = HFloat(5.527786294980874), theta(eta) = HFloat(0.2216792480031605), diff(theta(eta), eta) = HFloat(-0.38916037599841974)]

 

.4

 

[eta = 0., chi(eta) = HFloat(0.9999999999999998), diff(chi(eta), eta) = HFloat(4.580796631072469), f(eta) = HFloat(0.0), diff(f(eta), eta) = HFloat(0.0), diff(diff(f(eta), eta), eta) = HFloat(5.687607599246298), g(eta) = HFloat(0.0), diff(g(eta), eta) = HFloat(0.0), diff(diff(g(eta), eta), eta) = HFloat(5.687607599246298), phi(eta) = HFloat(-2.0568809171520708), diff(phi(eta), eta) = HFloat(6.0203396482123575), theta(eta) = HFloat(0.20686299926628238), diff(theta(eta), eta) = HFloat(-0.39656850036685887)]

(7)

NULL

 

 

Download 3DAKc2w_-_Copy.mw

hi .please help me for solve this equations.

bbb2.mw

restart; d[11] := 1; mu[11] := 1; q[311] := 1; d[33] := 1; mu[33] := 1; a[11] := 1; e[311] := 1; a[33] := 1; A := 1; g[111111] := 1; c[1111] := 1; g[113113] := 1; f[3113] := 1; beta[11] := 1; `ΔT` := 1; II := 1; L := 1

J := d[11]*(diff(Phi(x, z), x, x))+mu[11]*(diff(psi(x, z), x, x))+q[311]*(diff(w(x), x, x))+d[33]*(diff(Phi(x, z), z, z))+mu[33]*(diff(psi(x, z), z, z));

diff(diff(Phi(x, z), x), x)+diff(diff(psi(x, z), x), x)+diff(diff(w(x), x), x)+diff(diff(Phi(x, z), z), z)+diff(diff(psi(x, z), z), z)

(1)

B := a[11]*(diff(Phi(x, z), x, x))+d[11]*(diff(psi(x, z), x, x))+e[311]*(diff(w(x), x, x))+a[33]*(diff(Phi(x, z), z, z))+d[33]*(diff(psi(x, z), z, z));

diff(diff(Phi(x, z), x), x)+diff(diff(psi(x, z), x), x)+diff(diff(w(x), x), x)+diff(diff(Phi(x, z), z), z)+diff(diff(psi(x, z), z), z)

(2)

R := A*(g[111111]*(diff(u[0](x), x, x, x, x))-c[1111]*(diff(u[0](x), x, x)+(1/2)*(diff((diff(w(x), x))^2, x)))+e[311]*(diff(diff(Phi(x, z), z), x))+q[311]*(diff(diff(psi(x, z), z), x)));

diff(diff(diff(diff(u[0](x), x), x), x), x)-(diff(diff(u[0](x), x), x))-(diff(w(x), x))*(diff(diff(w(x), x), x))+diff(diff(Phi(x, z), x), z)+diff(diff(psi(x, z), x), z)

(3)

S := -II*g[111111]*(diff(w(x), x, x, x, x, x, x))-II*c[1111]*(diff(w(x), x, x, x, x))+A*g[113113]*(diff(w(x), x, x, x, x))-A*f[3113]*(diff(diff(Phi(x, z), z), x, x))-A*(c[1111]*(diff(u[0](x), x, x)+(1/2)*(diff((diff(w(x), x))^2, x)))+e[311]*(diff(diff(Phi(x, z), z), x))+q[311]*(diff(diff(psi(x, z), z), x)))*(diff(w(x), x))-A*(diff(w(x), x, x))*(c[1111]*(diff(u[0](x), x)+(1/2)*(diff(w(x), x))^2)+e[311]*(diff(Phi(x, z), z))+q[311]*(diff(psi(x, z), z))-beta[11]*`ΔT`);

-(diff(diff(diff(diff(diff(diff(w(x), x), x), x), x), x), x))-(diff(diff(diff(Phi(x, z), x), x), z))-(diff(diff(u[0](x), x), x)+(diff(w(x), x))*(diff(diff(w(x), x), x))+diff(diff(Phi(x, z), x), z)+diff(diff(psi(x, z), x), z))*(diff(w(x), x))-(diff(diff(w(x), x), x))*(diff(u[0](x), x)+(1/2)*(diff(w(x), x))^2+diff(Phi(x, z), z)+diff(psi(x, z), z)-1)

(4)

dsys := {B, J, R, S}; BCS := {D@@2*w(0) = 0, D@@2*w(L) = 0, Phi(x = 0) = 0, Phi(x = L) = 0, Phi(z = -(1/2)*h) = 0, Phi(z = (1/2)*h) = 0, psi(x = 0) = 0, psi(x = L) = 0, psi(z = -(1/2)*h) = 0, psi(z = (1/2)*h) = 0, w(x = 0) = 0, w(x = L) = 0, u[0](x = 0) = 0, u[0](x = L) = 0, (D(w))(0) = 0, (D(w))(L) = 0, (D(u[0]))(0) = 0, (D(u[0]))(L) = 0}

{D@@2*w(0) = 0, D@@2*w(L) = 0, Phi(x = 0) = 0, Phi(x = L) = 0, Phi(z = -(1/2)*h) = 0, Phi(z = (1/2)*h) = 0, psi(x = 0) = 0, psi(x = L) = 0, psi(z = -(1/2)*h) = 0, psi(z = (1/2)*h) = 0, w(x = 0) = 0, w(x = L) = 0, u[0](x = 0) = 0, u[0](x = L) = 0, (D(w))(0) = 0, (D(w))(L) = 0, (D(u[0]))(0) = 0, (D(u[0]))(L) = 0}

(5)

dsol5 := dsolve(dsys, numeric)

Error, (in dsolve/numeric/process_input) missing differential equations and initial or boundary conditions in the first argument: dsys

 

NULL

NULL

NULL

if former equations are not solvable , please help me for another way, in which at first two equation solve..in this way in equation [J and B] assume that q[311]=e[311]=0 and dsolve perform to find Φ and  ψ

after by finding Φ and  ψ is use for detemine w and u0

please see attached file below[bbb2_2.mw]

bbb2_2.mw

Download bbb2.mw

Dear All,

I am going to solve the following systems of ODEs but get the error: Newton iteration is not converging.
Could you please share your idea with me. In the case of AA=-0.2,0,0.2,0.4,...; I could get the solution.
Thank you in advance.


restart;
with(plots);
Pr := 2; Le := 2; nn := 2; Nb := .1; Nt := .1; QQ := .1; SS := .1; BB := .1; CC := .1; Ec := .1; MM := .2;AA:=-0.4;

Eq1 := diff(f(eta), `$`(eta, 3))+f(eta).(diff(f(eta), `$`(eta, 2)))-2.*nn/(nn+1).((diff(f(eta), eta))^2)-MM.(diff(f(eta), eta)) = 0; Eq2 := 1/Pr.(diff(theta(eta), `$`(eta, 2)))+f(eta).(diff(theta(eta), eta))-4.*nn/(nn+1).(diff(f(eta), eta)).theta(eta)+Nb.(diff(theta(eta), eta)).(diff(h(eta), eta))+Nt.((diff(theta(eta), eta))^2)+Ec.((diff(f(eta), `$`(eta, 2)))^2)-QQ.theta(eta) = 0;
Eq3 := diff(h(eta), `$`(eta, 2))+Le.f(eta).(diff(h(eta), eta))+Nt/Nb.(diff(theta(eta), `$`(eta, 2))) = 0;

bcs := f(0) = SS, (D(f))(0) = 1+AA.((D@@2)(f))(0), theta(0) = 1+BB.(D(theta))(0), phi(0) = 1+CC.(D(phi))(0), (D(f))(etainf) = 0, theta(etainf) = 0, phi(etainf) = 0

Error, (in dsolve/numeric/ComputeSolution) Newton iteration is not converging

I've been trying to numerically solve (and plot) this equation. Maple tells me that some matrix is singular - I have no idea, what I can do.

eq := diff(y(x), `$`(x, 3))-(diff(y(x), x))*(diff(y(x), x))+1 = 0;

cond := (D(y))(0) = 0, (D(y))(1) = 1, ((D@@2)(y))(0) = 0

de := dsolve({cond, eq}, y(x), numeric);

Error, (in dsolve/numeric/bvp) matrix is singular

I am unable to solve the attached optimal control problem,please any one who many help  me in guideing .tnx

restart:
unprotect('gamma');

L:=b[1]*c(t)+b[2]*i(t)+w[1]*(u[1])^2/2+w[2]*(u[2])^2/2+w[3]*(u[3])^2/2;
1 2 1 2 1 2
b[1] c(t) + b[2] i(t) + - w[1] u[1] + - w[2] u[2] + - w[3] u[3]
2 2 2
H:=L+lambda[1](t)*((1-p*Psi)*tau+phi* v + delta *r-lambda*(1-u[3])*s-u[1]*varphi*s -mu*s ) +lambda[2](t)*(p*Psi*tau + u[1]*vartheta*s -gamma*lambda* (1-u[3])*v-(mu+phi)*v ) +lambda[3](t)*( (1-u[3])*rho*lambda* (s +gamma*v)+(1-q)* u[2]*eta*i -(mu +beta +chi)*c ) +lambda[4](t)* ((1-rho)*(1-u[3])*lambda*( s +gamma*v) +chi*c - u[2]*eta*i - (mu +alpha )*i) +lambda[5](t)*( beta*c + u[2]*q*eta*i -(mu +delta)*r);
1 2 1 2 1 2
b[1] c(t) + b[2] i(t) + - w[1] u[1] + - w[2] u[2] + - w[3] u[3] + lambda[1](t
2 2 2

) ((1 - p Psi) tau + phi v + delta r - lambda (1 - u[3]) s - u[1] varphi s

- mu s) + lambda[2](t) (p Psi tau + u[1] vartheta s

- gamma lambda (1 - u[3]) v - (mu + phi) v) + lambda[3](t) ((1 - u[3]) rho

lambda (s + gamma v) + (1 - q) u[2] eta i - (mu + beta + chi) c) + lambda[4](t

) ((1 - rho) (1 - u[3]) lambda (s + gamma v) + chi c - u[2] eta i

- (mu + alpha) i) + lambda[5](t) (beta c + u[2] q eta i - (mu + delta) r)
du1:=diff(H,u[1]);

w[1] u[1] - lambda[1](t) varphi s + lambda[2](t) vartheta s
du2:=diff(H,u[2]);du3:=diff(H,u[3]);
w[2] u[2] + lambda[3](t) (1 - q) eta i - lambda[4](t) eta i

+ lambda[5](t) q eta i
w[3] u[3] + lambda[1](t) lambda s + lambda[2](t) gamma lambda v

- lambda[3](t) rho lambda (s + gamma v)

- lambda[4](t) (1 - rho) lambda (s + gamma v)

ddu1 := -A[1] u[1] + psi[1](t) beta x[1] x[3] - psi[2](t) beta x[1] x[3]

ddu2 := -A[2] u[2] - psi[3](t) k x[2]
sol_u1 := solve(du1, u[1]);
s(t) (lambda[1](t) varphi - lambda[2](t) vartheta)
--------------------------------------------------
w[1]
sol_u2 := solve(du2, u[2]);sol_u3 := solve(du3, u[3]);
eta i (-lambda[3](t) + lambda[3](t) q + lambda[4](t) - lambda[5](t) q)
----------------------------------------------------------------------
w[2]
1
---- (lambda (-lambda[1](t) s - lambda[2](t) gamma v + lambda[3](t) rho s
w[3]

+ lambda[3](t) rho gamma v + lambda[4](t) s + lambda[4](t) gamma v

- lambda[4](t) rho s - lambda[4](t) rho gamma v))
Dx2:=subs(u[1]= s*(lambda[1](t)*varphi-lambda[2](t)*vartheta)/w[1] ,u[2]= eta*i*(-lambda[3](t)+lambda[3](t)*q+lambda[4](t)-lambda[5](t)*q)/w[2], u[3]=-lambda*(lambda[1](t)*s+lambda[2](t)*gamma*v-lambda[3](t)*rho*s-lambda[3](t)*rho*gamma*v-lambda[4](t)*s-lambda[4](t)*gamma*v+lambda[4](t)*rho*s+lambda[4](t)*rho*gamma*v)/w[3] ,H );
2 2
s (lambda[1](t) varphi - lambda[2](t) vartheta)
b[1] c(t) + b[2] i(t) + -------------------------------------------------
2 w[1]

2 2 2
eta i (-lambda[3](t) + lambda[3](t) q + lambda[4](t) - lambda[5](t) q)
+ ------------------------------------------------------------------------- +
2 w[2]

1 / 2
------ \lambda (lambda[1](t) s + lambda[2](t) gamma v - lambda[3](t) rho s
2 w[3]

- lambda[3](t) rho gamma v - lambda[4](t) s - lambda[4](t) gamma v

/
\ |
+ lambda[4](t) rho s + lambda[4](t) rho gamma v)^2/ + lambda[1](t) |(1
\

/ 1
- p Psi) tau + phi v + delta r - lambda |1 + ---- (lambda (lambda[1](t) s
\ w[3]

+ lambda[2](t) gamma v - lambda[3](t) rho s - lambda[3](t) rho gamma v

- lambda[4](t) s - lambda[4](t) gamma v + lambda[4](t) rho s

\
+ lambda[4](t) rho gamma v))| s
/

2 \
s (lambda[1](t) varphi - lambda[2](t) vartheta) varphi |
- ------------------------------------------------------- - mu s| +
w[1] /

/
|
lambda[2](t) |p Psi tau
\

2
s (lambda[1](t) varphi - lambda[2](t) vartheta) vartheta /
+ --------------------------------------------------------- - gamma lambda |1 +
w[1] \

1
---- (lambda (lambda[1](t) s + lambda[2](t) gamma v - lambda[3](t) rho s
w[3]

- lambda[3](t) rho gamma v - lambda[4](t) s - lambda[4](t) gamma v

\
\ |
+ lambda[4](t) rho s + lambda[4](t) rho gamma v))| v - (mu + phi) v| +
/ /

// 1
lambda[3](t) ||1 + ---- (lambda (lambda[1](t) s + lambda[2](t) gamma v
\\ w[3]

- lambda[3](t) rho s - lambda[3](t) rho gamma v - lambda[4](t) s

\
- lambda[4](t) gamma v + lambda[4](t) rho s + lambda[4](t) rho gamma v))|
/

1 / 2 2
rho lambda (s + gamma v) + ---- \(1 - q) eta i (-lambda[3](t)
w[2]

\ \
+ lambda[3](t) q + lambda[4](t) - lambda[5](t) q)/ - (mu + beta + chi) c| +
/

/
| / 1
lambda[4](t) |(1 - rho) |1 + ---- (lambda (lambda[1](t) s
\ \ w[3]

+ lambda[2](t) gamma v - lambda[3](t) rho s - lambda[3](t) rho gamma v

- lambda[4](t) s - lambda[4](t) gamma v + lambda[4](t) rho s

\
+ lambda[4](t) rho gamma v))| lambda (s + gamma v) + chi c
/

2 2
eta i (-lambda[3](t) + lambda[3](t) q + lambda[4](t) - lambda[5](t) q)
- ------------------------------------------------------------------------
w[2]

\ /
| |
- (mu + alpha) i| + lambda[5](t) |beta c
/ \

+

2 2
eta i (-lambda[3](t) + lambda[3](t) q + lambda[4](t) - lambda[5](t) q) q
--------------------------------------------------------------------------
w[2]

\
|
- (mu + delta) r|
/
ode1:=diff(lambda[1](t),t)=-diff(H,s);ode2:=diff(lambda[2](t),t)=-diff(H,v);ode3:=diff(psi[3](t),t)=-diff(H,c);ode4:=diff(lambda[4](t),t)=-diff(H,i);ode5:=diff(lambda[5](t),t)=-diff(H,r);
d
--- lambda[1](t) = -lambda[1](t) (-lambda (1 - u[3]) - u[1] varphi - mu)
dt

- lambda[2](t) u[1] vartheta - lambda[3](t) (1 - u[3]) rho lambda

- lambda[4](t) (1 - rho) (1 - u[3]) lambda
d
--- lambda[2](t) = -lambda[1](t) phi
dt

- lambda[2](t) (-gamma lambda (1 - u[3]) - mu - phi)

- lambda[3](t) (1 - u[3]) rho lambda gamma

- lambda[4](t) (1 - rho) (1 - u[3]) lambda gamma
d
--- psi[3](t) = -lambda[3](t) (-mu - beta - chi) - lambda[4](t) chi
dt

- lambda[5](t) beta
d
--- lambda[4](t) = -lambda[3](t) (1 - q) u[2] eta
dt

- lambda[4](t) (-u[2] eta - mu - alpha) - lambda[5](t) u[2] q eta
d
--- lambda[5](t) = -lambda[1](t) delta - lambda[5](t) (-mu - delta)
dt
restart:
#Digits:=10:


unprotect('gamma');
lambda:=0.51:
mu:=0.002:
beta:=0.0115:
delta:=0.003:
alpha:=0.33:
chi:=0.00274:
k:=6.24:
gamma:=0.4:
rho:=0.338:;tau=1000:;Psi:=0.1:;p:=0.6:;phi:=0.001:;eta:=0.001124:q:=0.6:varphi:=0.9:;vatheta:=0.9:
b[1]:=2:;b[2]:=3:;w[1]:=4:;w[2]:=5:;w[3]:=6:
#u[1]:=s(t)*(lambda[1](t)*varphi-lambda[2](t)*vartheta)/w[1]:
#u[2]:=eta*i*(-lambda[3](t)+lambda[3](t)*q+lambda[4](t)-lambda[5](t)*q)/w[2]:;u[3]:=lambda*(-lambda[1](t)*s-lambda[2](t)*gamma*v+lambda[3](t)*rho*s+lambda[3](t)*rho*gamma*v+lambda[4](t)*s+lambda[4](t)*gamma*v-lambda[4](t)*rho*s-lambda[4](t)*rho*gamma*v)/w[3]:
ics := s(0)=8200, v(0)=2800,c(0)=1100,i(0)=1500,r(0)=200,lambda[1](20)=0,lambda[2](20)=0,lambda[3](20)=0,lambda[4](20)=0,lambda[5](20)=0:
ode1:=diff(s(t),t)=(1-p*Psi)*tau+phi* v(t) + delta *r(t)-lambda*(1-u[3])*s(t)-u[1]*varphi*s(t) -mu*s(t),
diff(v(t), t) =p*Psi*tau + u[1]*vartheta*s(t) -gamma*lambda* (1-u[3])*v(t)-(mu+phi)*v(t) ,
diff(c(t), t) =(1-u[3])*rho*lambda* (s(t) +gamma*v(t))+(1-q)* u[2]*eta*i(t) -(mu +beta +chi)*c(t),
diff(i(t), t) =(1-rho)*(1-u[3])*lambda*( s(t) +gamma*v(t)) +chi*c(t) - u[2]*eta*i(t) - (mu +alpha )*i(t),
diff(r(t), t) = beta*c(t) + u[2]*q*eta*i(t) -(mu +delta)*r(t),
diff(lambda[1](t), t) = -lambda[1](t)*(-lambda*(1-u[3])-u[1]*varphi-mu)-lambda[2](t)*u[1]*vartheta-lambda[3](t)*(1-u[3])*rho*lambda-lambda[4](t)*(1-rho)*(1-u[3])*lambda,diff(lambda[2](t),t)=-lambda[1](t)*phi-lambda[2](t)*(-gamma*lambda*(1-u[3])-mu-phi)-lambda[3](t)*(1-u[3])*rho*lambda*gamma-lambda[4](t)*(1-rho)*(1-u[3])*lambda*gamma,diff(lambda[3](t),t)=-lambda[3](t)*(-mu-beta-chi)-lambda[4](t)*chi-lambda[5](t)*beta,diff(lambda[4](t),t)=-lambda[3](t)*(1-q)*u[2]*eta-lambda[4](t)*(-u[2]*eta-mu-alpha)-lambda[5](t)*u[2]*q*eta,diff(lambda[5](t),t)=-lambda[1](t)*delta-lambda[5](t)*(-mu-delta);
d
--- s(t) = (1 - p Psi) tau + phi v(t) + delta r(t) - lambda (1 - u[3]) s(t)
dt

d
- u[1] varphi s(t) - mu s(t), --- v(t) = p Psi tau + u[1] vartheta s(t)
dt

d
- gamma lambda (1 - u[3]) v(t) - (mu + phi) v(t), --- c(t) = (1 - u[3]) rho lambda
dt

(s(t) + gamma v(t)) + (1 - q) u[2] eta - (mu + beta + chi) c(t), 0 = (1

- rho) (1 - u[3]) lambda (s(t) + gamma v(t)) + chi c(t) - u[2] eta - mu

d d
- alpha, --- r(t) = beta c(t) + u[2] q eta - (mu + delta) r(t), ---
dt dt

lambda[1](t) = -lambda[1](t) (-lambda (1 - u[3]) - u[1] varphi - mu)

- lambda[2](t) u[1] vartheta - lambda[3](t) (1 - u[3]) rho lambda

d
- lambda[4](t) (1 - rho) (1 - u[3]) lambda, --- lambda[2](t) =
dt
-lambda[1](t) phi - lambda[2](t) (-gamma lambda (1 - u[3]) - mu - phi)

- lambda[3](t) (1 - u[3]) rho lambda gamma

d
- lambda[4](t) (1 - rho) (1 - u[3]) lambda gamma, --- lambda[3](t) =
dt
d
-lambda[3](t) (-mu - beta - chi) - lambda[4](t) chi - lambda[5](t) beta, ---
dt

lambda[4](t) = -lambda[3](t) (1 - q) u[2] eta

- lambda[4](t) (-u[2] eta - mu - alpha) - lambda[5](t) u[2] q eta,

d
--- lambda[5](t) = -lambda[1](t) delta - lambda[5](t) (-mu - delta)
dt

sol := dsolve({c(0) = 0, i(0) = 0, r(0) = .1, s(0) = 0, v(0) = 0, diff(c(t), t) = (1-u[3])*rho*lambda*(s(t)+gamma*v(t))+(1-q)*u[2]*eta*i(t)-(mu+beta+chi)*c(t), diff(i(t), t) = (1-rho)*(1-u[3])*lambda*(s(t)+gamma*v(t))+chi*c(t)-u[2]*eta*i(t)-(mu+alpha)*i(t), diff(r(t), t) = beta*c(t)+u[2]*q*eta*i(t)-(mu+delta)*r(t), diff(s(t), t) = (1-p*Psi)*tau+phi*v(t)+delta*r(t)-lambda*(1-u[3])*s(t)-u[1]*varphi*s(t)-mu*s(t), diff(v(t), t) = p*Psi*tau+u[1]*vartheta*s(t)-gamma*lambda*(1-u[3])*v(t)-(mu+phi)*v(t), diff(lambda[1](t), t) = -lambda[1](t)*(-lambda*(1-u[3])-u[1]*varphi-mu)-lambda[2](t)*u[1]*vartheta-lambda[3](t)*(1-u[3])*rho*lambda-lambda[4](t)*(1-rho)*(1-u[3])*lambda, diff(lambda[2](t), t) = -lambda[1](t)*phi-lambda[2](t)*(-gamma*lambda*(1-u[3])-mu-phi)-lambda[3](t)*(1-u[3])*rho*lambda*gamma-lambda[4](t)*(1-rho)*(1-u[3])*lambda*gamma, diff(lambda[3](t), t) = -lambda[3](t)*(-mu-beta-chi)-lambda[4](t)*chi-lambda[5](t)*beta, diff(lambda[4](t), t) = -lambda[3](t)*(1-q)*u[2]*eta-lambda[4](t)*(-u[2]*eta-mu-alpha)-lambda[5](t)*u[2]*q*eta, diff(lambda[5](t), t) = -lambda[1](t)*delta-lambda[5](t)*(-mu-delta), lambda[1](20) = 0, lambda[2](20) = 0, lambda[3](20) = 0, lambda[4](20) = 0, lambda[5](20) = 0}, type = numeric);
Error, (in dsolve/numeric/process_input) invalid specification of initial conditions, got 1 = 0

sol:=dsolve([ode1,ics],numeric, method = bvp[midrich],maxmesh=500);

Error, (in dsolve/numeric/process_input) system must be entered as a set/list of expressions/equations

dsolve[':-interactive']({});
Error, `:=` unexpected
sol:=dsolve([ode1,ics],numeric, method = bvp[midrich],maxmesh=500);
Error, (in dsolve/numeric/process_input) system must be entered as a set/list of expressions/equations

eq1:=diff(s(t), t)=(1-p*Psi)*tau+phi* v(t) + delta *r(t)-lambda*(1-u[3])*s(t)-u[1]*varphi*s(t) -mu*s(t);
eq2:diff(v(t), t) =p*Psi*tau + u[1]*vartheta*s(t) -gamma*lambda* (1-u[3])*v(t)-(mu+phi)*v(t);
eq3:=diff(c(t), t) =(1-u[3])*rho*lambda* (s(t) +gamma*v(t))+(1-q)* u[2]*eta*i(t) -(mu +beta +chi)*c(t);
eq4:=diff(i(t), t) =(1-rho)*(1-u[3])*lambda*( s(t) +gamma*v(t)) +chi*c(t) - u[2]*eta*i(t) - (mu +alpha )*i(t);
eq5:=diff(r(t), t) = beta*c(t) + u[2]*q*eta*i(t) -(mu +delta)*r(t);

d
--- s(t) = (1 - p Psi) tau + phi v(t) + delta r(t) - lambda (1 - u[3]) s(t)
dt

- u[1] varphi s(t) - mu s(t)
d
--- v(t) = p Psi tau + u[1] vartheta s(t) - gamma lambda (1 - u[3]) v(t)
dt

- (mu + phi) v(t)
d
--- c(t) = (1 - u[3]) rho lambda (s(t) + gamma v(t)) + (1 - q) u[2] eta i(t)
dt

- (mu + beta + chi) c(t)
d
--- i(t) = (1 - rho) (1 - u[3]) lambda (s(t) + gamma v(t)) + chi c(t)
dt

- u[2] eta i(t) - (mu + alpha) i(t)
d
--- r(t) = beta c(t) + u[2] q eta i(t) - (mu + delta) r(t)
dt
eq6:=diff(Q(t),t)=b[1]*c(t)+b[2]*i(t)+w[1]*(u[1])^2/2+w[2]*(u[2])^2/2+w[3]*(u[3])^2/2;
d 1 2 1 2 1 2
--- Q(t) = b[1] c(t) + b[2] i(t) + - w[1] u[1] + - w[2] u[2] + - w[3] u[3]
dt 2 2 2
ics:=s(0)=8200, v(0)=2800,c(0)=1100,i(0)=1500,r(0)=200,Q(0)=6700;
s(0) = 8200, v(0) = 2800, c(0) = 1100, i(0) = 1500, r(0) = 200, Q(0) = 6700
sol0:=dsolve({eq1,eq2,eq3,eq4,eq5,eq6,ics},type=numeric,stiff=true,'parameters'=[u[1],u[2],u[3]],abserr=1e-15,relerr=1e-12,maxfun=0,range=0..50):
Error, (in dsolve/numeric/process_input) system must be entered as a set/list of expressions/equations
with(plots):
Q0:=6700;
6700
obj:=proc(u)
global sol0,Q0;
local ob1;
try
sol0('parameters'=[u[1],u[2],u[3]]):
ob1:=subs(sol0(20.),Q(t)):
catch :
ob1:=0;
end try;
#ob1:=subs(sol0(20.),Q(t));
if ob1>Q0 then Q0:=ob1;print(Q0,u);end;
ob1;
end proc;
proc(u) ... end;
obj([1,1,1]);
0
obj([3,2.5],4);
0
u0:=Vector(3,[0.,0.,0.],datatype=float[8]);
Vector[column](%id = 85973880)

Q0:=0;
Q0 := 0
with(Optimization);
[ImportMPS, Interactive, LPSolve, LSSolve, Maximize, Minimize, NLPSolve,

QPSolve]
sol2:=NLPSolve(3,obj,initialpoint=u0,method=nonlinearsimplex,maximize,evaluationlimit=100):
sol0('parameters'=[3.18125786060723, 2.36800986932868]);
sol0(parameters = [3.18125786060723, 2.36800986932868])
for i from 1 to 3 do odeplot(sol0,[t,x[i](t)],0..20,thickness=3,axes=boxed);od;
Error, (in plots/odeplot) input is not a valid dsolve/numeric solution

 

hello , 

how i can exract value from pdsolve ,i need to use dU(x,R)/dR 

thank you 

 

restart; with(plots)

n := 1/3;

1/3

(1)

Uu := (3*n+1)*(1-R^((n+1)/n))/(n+1);

-(3/2)*R^4+3/2

(2)

eq := Uu*(diff(theta(x, R), x))-4*(diff(R*(diff(theta(x, R), R)), R))/R;

(-(3/2)*R^4+3/2)*(diff(theta(x, R), x))-4*(diff(theta(x, R), R)+R*(diff(diff(theta(x, R), R), R)))/R

(3)

IBC := {theta(0, R) = 1, theta(x, 1) = 0, (D[2](theta))(x, 0) = 0};

{theta(0, R) = 1, theta(x, 1) = 0, (D[2](theta))(x, 0) = 0}

(4)

pds := pdsolve(eq, IBC, numeric);

module () local INFO; export plot, plot3d, animate, value, settings; option `Copyright (c) 2001 by Waterloo Maple Inc. All rights reserved.`; end module

(5)

U := subs(pds:-value(output = listprocedure), theta(x, R));

proc () local tv, xv, solnproc, stype, ndsol, vals; option `Copyright (c) 2001 by Waterloo Maple Inc. All rights reserved.`; Digits := trunc(evalhf(Digits)); solnproc := proc (tv, xv) local INFO, errest, nd, dvars, dary, daryt, daryx, vals, msg, i, j; option `Copyright (c) 2001 by Waterloo Maple Inc. All rights reserved.`; table( [( "soln_procedures" ) = array( 1 .. 1, [( 1 ) = (18446744074366926358)  ] ) ] ) INFO := table( [( "timestep" ) = 0.500000000000000e-1, ( "IBC" ) = b, ( "spaceidx" ) = 2, ( "fdepvars" ) = [theta(x, R)], ( "dependson" ) = [{1}], ( "eqnords" ) = [[1, 2]], ( "intspace" ) = Matrix(21, 1, {(1, 1) = .0, (2, 1) = .0, (3, 1) = .0, (4, 1) = .0, (5, 1) = .0, (6, 1) = .0, (7, 1) = .0, (8, 1) = .0, (9, 1) = .0, (10, 1) = .0, (11, 1) = .0, (12, 1) = .0, (13, 1) = .0, (14, 1) = .0, (15, 1) = .0, (16, 1) = .0, (17, 1) = .0, (18, 1) = .0, (19, 1) = .0, (20, 1) = .0, (21, 1) = .0}, datatype = float[8], order = C_order), ( "solvec2" ) = Vector(21, {(1) = .0, (2) = .0, (3) = .0, (4) = .0, (5) = .0, (6) = .0, (7) = .0, (8) = .0, (9) = .0, (10) = .0, (11) = .0, (12) = .0, (13) = .0, (14) = .0, (15) = .0, (16) = .0, (17) = .0, (18) = .0, (19) = .0, (20) = .0, (21) = .0}, datatype = float[8]), ( "allocspace" ) = 21, ( "solmat_ne" ) = 0, ( "depords" ) = [[1, 2]], ( "BCS", 1 ) = {[[1, 0, 1], b[1, 0, 1]], [[1, 1, 0], b[1, 1, 0]]}, ( "spacepts" ) = 21, ( "solvec3" ) = Vector(21, {(1) = .0, (2) = .0, (3) = .0, (4) = .0, (5) = .0, (6) = .0, (7) = .0, (8) = .0, (9) = .0, (10) = .0, (11) = .0, (12) = .0, (13) = .0, (14) = .0, (15) = .0, (16) = .0, (17) = .0, (18) = .0, (19) = .0, (20) = .0, (21) = .0}, datatype = float[8]), ( "autonomous" ) = true, ( "vectorproc" ) = proc (v, vp, vpp, t, x, k, h, n, vec) local _s1, _s2, _s3, _s4, _s5, _s6, xi; _s3 := 4*k; _s4 := -3*h^2; _s5 := 2*h*k; _s6 := 2*k*h^2; vec[1] := 0; vec[n] := 0; for xi from 2 to n-1 do _s1 := -vp[xi-1]+vp[xi+1]; _s2 := vp[xi-1]-2*vp[xi]+vp[xi+1]; vec[xi] := (_s4*vp[xi]*x[xi]^5+_s2*_s3*x[xi]-_s4*vp[xi]*x[xi]+_s1*_s5)/(_s6*x[xi]) end do end proc, ( "timeidx" ) = 1, ( "extrabcs" ) = [0], ( "pts", R ) = [0, 1], ( "solvec5" ) = 0, ( "timevar" ) = x, ( "t0" ) = 0, ( "solmat_v" ) = Vector(147, {(1) = .0, (2) = .0, (3) = .0, (4) = .0, (5) = .0, (6) = .0, (7) = .0, (8) = .0, (9) = .0, (10) = .0, (11) = .0, (12) = .0, (13) = .0, (14) = .0, (15) = .0, (16) = .0, (17) = .0, (18) = .0, (19) = .0, (20) = .0, (21) = .0, (22) = .0, (23) = .0, (24) = .0, (25) = .0, (26) = .0, (27) = .0, (28) = .0, (29) = .0, (30) = .0, (31) = .0, (32) = .0, (33) = .0, (34) = .0, (35) = .0, (36) = .0, (37) = .0, (38) = .0, (39) = .0, (40) = .0, (41) = .0, (42) = .0, (43) = .0, (44) = .0, (45) = .0, (46) = .0, (47) = .0, (48) = .0, (49) = .0, (50) = .0, (51) = .0, (52) = .0, (53) = .0, (54) = .0, (55) = .0, (56) = .0, (57) = .0, (58) = .0, (59) = .0, (60) = .0, (61) = .0, (62) = .0, (63) = .0, (64) = .0, (65) = .0, (66) = .0, (67) = .0, (68) = .0, (69) = .0, (70) = .0, (71) = .0, (72) = .0, (73) = .0, (74) = .0, (75) = .0, (76) = .0, (77) = .0, (78) = .0, (79) = .0, (80) = .0, (81) = .0, (82) = .0, (83) = .0, (84) = .0, (85) = .0, (86) = .0, (87) = .0, (88) = .0, (89) = .0, (90) = .0, (91) = .0, (92) = .0, (93) = .0, (94) = .0, (95) = .0, (96) = .0, (97) = .0, (98) = .0, (99) = .0, (100) = .0, (101) = .0, (102) = .0, (103) = .0, (104) = .0, (105) = .0, (106) = .0, (107) = .0, (108) = .0, (109) = .0, (110) = .0, (111) = .0, (112) = .0, (113) = .0, (114) = .0, (115) = .0, (116) = .0, (117) = .0, (118) = .0, (119) = .0, (120) = .0, (121) = .0, (122) = .0, (123) = .0, (124) = .0, (125) = .0, (126) = .0, (127) = .0, (128) = .0, (129) = .0, (130) = .0, (131) = .0, (132) = .0, (133) = .0, (134) = .0, (135) = .0, (136) = .0, (137) = .0, (138) = .0, (139) = .0, (140) = .0, (141) = .0, (142) = .0, (143) = .0, (144) = .0, (145) = .0, (146) = .0, (147) = .0}, datatype = float[8], order = C_order, attributes = [source_rtable = (Matrix(21, 7, {(1, 1) = .0, (1, 2) = .0, (1, 3) = .0, (1, 4) = .0, (1, 5) = .0, (1, 6) = .0, (1, 7) = .0, (2, 1) = .0, (2, 2) = .0, (2, 3) = .0, (2, 4) = .0, (2, 5) = .0, (2, 6) = .0, (2, 7) = .0, (3, 1) = .0, (3, 2) = .0, (3, 3) = .0, (3, 4) = .0, (3, 5) = .0, (3, 6) = .0, (3, 7) = .0, (4, 1) = .0, (4, 2) = .0, (4, 3) = .0, (4, 4) = .0, (4, 5) = .0, (4, 6) = .0, (4, 7) = .0, (5, 1) = .0, (5, 2) = .0, (5, 3) = .0, (5, 4) = .0, (5, 5) = .0, (5, 6) = .0, (5, 7) = .0, (6, 1) = .0, (6, 2) = .0, (6, 3) = .0, (6, 4) = .0, (6, 5) = .0, (6, 6) = .0, (6, 7) = .0, (7, 1) = .0, (7, 2) = .0, (7, 3) = .0, (7, 4) = .0, (7, 5) = .0, (7, 6) = .0, (7, 7) = .0, (8, 1) = .0, (8, 2) = .0, (8, 3) = .0, (8, 4) = .0, (8, 5) = .0, (8, 6) = .0, (8, 7) = .0, (9, 1) = .0, (9, 2) = .0, (9, 3) = .0, (9, 4) = .0, (9, 5) = .0, (9, 6) = .0, (9, 7) = .0, (10, 1) = .0, (10, 2) = .0, (10, 3) = .0, (10, 4) = .0, (10, 5) = .0, (10, 6) = .0, (10, 7) = .0, (11, 1) = .0, (11, 2) = .0, (11, 3) = .0, (11, 4) = .0, (11, 5) = .0, (11, 6) = .0, (11, 7) = .0, (12, 1) = .0, (12, 2) = .0, (12, 3) = .0, (12, 4) = .0, (12, 5) = .0, (12, 6) = .0, (12, 7) = .0, (13, 1) = .0, (13, 2) = .0, (13, 3) = .0, (13, 4) = .0, (13, 5) = .0, (13, 6) = .0, (13, 7) = .0, (14, 1) = .0, (14, 2) = .0, (14, 3) = .0, (14, 4) = .0, (14, 5) = .0, (14, 6) = .0, (14, 7) = .0, (15, 1) = .0, (15, 2) = .0, (15, 3) = .0, (15, 4) = .0, (15, 5) = .0, (15, 6) = .0, (15, 7) = .0, (16, 1) = .0, (16, 2) = .0, (16, 3) = .0, (16, 4) = .0, (16, 5) = .0, (16, 6) = .0, (16, 7) = .0, (17, 1) = .0, (17, 2) = .0, (17, 3) = .0, (17, 4) = .0, (17, 5) = .0, (17, 6) = .0, (17, 7) = .0, (18, 1) = .0, (18, 2) = .0, (18, 3) = .0, (18, 4) = .0, (18, 5) = .0, (18, 6) = .0, (18, 7) = .0, (19, 1) = .0, (19, 2) = .0, (19, 3) = .0, (19, 4) = .0, (19, 5) = .0, (19, 6) = .0, (19, 7) = .0, (20, 1) = .0, (20, 2) = .0, (20, 3) = .0, (20, 4) = .0, (20, 5) = .0, (20, 6) = .0, (20, 7) = .0, (21, 1) = .0, (21, 2) = .0, (21, 3) = .0, (21, 4) = .0, (21, 5) = .0, (21, 6) = .0, (21, 7) = .0}, datatype = float[8], order = C_order))]), ( "indepvars" ) = [x, R], ( "maxords" ) = [1, 2], ( "solvec1" ) = Vector(21, {(1) = .0, (2) = .0, (3) = .0, (4) = .0, (5) = .0, (6) = .0, (7) = .0, (8) = .0, (9) = .0, (10) = .0, (11) = .0, (12) = .0, (13) = .0, (14) = .0, (15) = .0, (16) = .0, (17) = .0, (18) = .0, (19) = .0, (20) = .0, (21) = .0}, datatype = float[8]), ( "startup_only" ) = false, ( "solvec4" ) = 0, ( "explicit" ) = false, ( "solmatrix" ) = Matrix(21, 7, {(1, 1) = .0, (1, 2) = .0, (1, 3) = .0, (1, 4) = .0, (1, 5) = .0, (1, 6) = .0, (1, 7) = .0, (2, 1) = .0, (2, 2) = .0, (2, 3) = .0, (2, 4) = .0, (2, 5) = .0, (2, 6) = .0, (2, 7) = .0, (3, 1) = .0, (3, 2) = .0, (3, 3) = .0, (3, 4) = .0, (3, 5) = .0, (3, 6) = .0, (3, 7) = .0, (4, 1) = .0, (4, 2) = .0, (4, 3) = .0, (4, 4) = .0, (4, 5) = .0, (4, 6) = .0, (4, 7) = .0, (5, 1) = .0, (5, 2) = .0, (5, 3) = .0, (5, 4) = .0, (5, 5) = .0, (5, 6) = .0, (5, 7) = .0, (6, 1) = .0, (6, 2) = .0, (6, 3) = .0, (6, 4) = .0, (6, 5) = .0, (6, 6) = .0, (6, 7) = .0, (7, 1) = .0, (7, 2) = .0, (7, 3) = .0, (7, 4) = .0, (7, 5) = .0, (7, 6) = .0, (7, 7) = .0, (8, 1) = .0, (8, 2) = .0, (8, 3) = .0, (8, 4) = .0, (8, 5) = .0, (8, 6) = .0, (8, 7) = .0, (9, 1) = .0, (9, 2) = .0, (9, 3) = .0, (9, 4) = .0, (9, 5) = .0, (9, 6) = .0, (9, 7) = .0, (10, 1) = .0, (10, 2) = .0, (10, 3) = .0, (10, 4) = .0, (10, 5) = .0, (10, 6) = .0, (10, 7) = .0, (11, 1) = .0, (11, 2) = .0, (11, 3) = .0, (11, 4) = .0, (11, 5) = .0, (11, 6) = .0, (11, 7) = .0, (12, 1) = .0, (12, 2) = .0, (12, 3) = .0, (12, 4) = .0, (12, 5) = .0, (12, 6) = .0, (12, 7) = .0, (13, 1) = .0, (13, 2) = .0, (13, 3) = .0, (13, 4) = .0, (13, 5) = .0, (13, 6) = .0, (13, 7) = .0, (14, 1) = .0, (14, 2) = .0, (14, 3) = .0, (14, 4) = .0, (14, 5) = .0, (14, 6) = .0, (14, 7) = .0, (15, 1) = .0, (15, 2) = .0, (15, 3) = .0, (15, 4) = .0, (15, 5) = .0, (15, 6) = .0, (15, 7) = .0, (16, 1) = .0, (16, 2) = .0, (16, 3) = .0, (16, 4) = .0, (16, 5) = .0, (16, 6) = .0, (16, 7) = .0, (17, 1) = .0, (17, 2) = .0, (17, 3) = .0, (17, 4) = .0, (17, 5) = .0, (17, 6) = .0, (17, 7) = .0, (18, 1) = .0, (18, 2) = .0, (18, 3) = .0, (18, 4) = .0, (18, 5) = .0, (18, 6) = .0, (18, 7) = .0, (19, 1) = .0, (19, 2) = .0, (19, 3) = .0, (19, 4) = .0, (19, 5) = .0, (19, 6) = .0, (19, 7) = .0, (20, 1) = .0, (20, 2) = .0, (20, 3) = .0, (20, 4) = .0, (20, 5) = .0, (20, 6) = .0, (20, 7) = .0, (21, 1) = .0, (21, 2) = .0, (21, 3) = .0, (21, 4) = .0, (21, 5) = .0, (21, 6) = .0, (21, 7) = .0}, datatype = float[8], order = C_order), ( "depvars" ) = [theta], ( "solmat_is" ) = 0, ( "adjusted" ) = false, ( "matrixhf" ) = true, ( "norigdepvars" ) = 1, ( "stages" ) = 1, ( "theta" ) = 1/2, ( "ICS" ) = [1], ( "multidep" ) = [false, false], ( "soltimes" ) = Vector(3, {(1) = .0, (2) = .0, (3) = .0}, datatype = float[8]), ( "depeqn" ) = [1], ( "method" ) = theta, ( "depshift" ) = [1], ( "depdords" ) = [[[1, 2]]], ( "matrixproc" ) = proc (v, vp, vpp, t, x, k, h, n, mat) local _s1, _s2, _s3, xi; _s1 := h^2; _s2 := -(3/2)/k; _s3 := (1/2)*(8*k+3*_s1)/(k*h^2); mat[3] := -(3/2)/h; mat[4] := 2/h; mat[5] := -(1/2)/h; mat[7*n-4] := 1; for xi from 2 to n-1 do mat[7*xi-4] := _s2*x[xi]^4+_s3; mat[7*xi-5] := (h-2*x[xi])/(_s1*x[xi]); mat[7*xi-3] := -(h+2*x[xi])/(_s1*x[xi]) end do end proc, ( "solution" ) = Array(1..3, 1..21, 1..1, {(1, 1, 1) = .0, (1, 2, 1) = .0, (1, 3, 1) = .0, (1, 4, 1) = .0, (1, 5, 1) = .0, (1, 6, 1) = .0, (1, 7, 1) = .0, (1, 8, 1) = .0, (1, 9, 1) = .0, (1, 10, 1) = .0, (1, 11, 1) = .0, (1, 12, 1) = .0, (1, 13, 1) = .0, (1, 14, 1) = .0, (1, 15, 1) = .0, (1, 16, 1) = .0, (1, 17, 1) = .0, (1, 18, 1) = .0, (1, 19, 1) = .0, (1, 20, 1) = .0, (1, 21, 1) = .0, (2, 1, 1) = .0, (2, 2, 1) = .0, (2, 3, 1) = .0, (2, 4, 1) = .0, (2, 5, 1) = .0, (2, 6, 1) = .0, (2, 7, 1) = .0, (2, 8, 1) = .0, (2, 9, 1) = .0, (2, 10, 1) = .0, (2, 11, 1) = .0, (2, 12, 1) = .0, (2, 13, 1) = .0, (2, 14, 1) = .0, (2, 15, 1) = .0, (2, 16, 1) = .0, (2, 17, 1) = .0, (2, 18, 1) = .0, (2, 19, 1) = .0, (2, 20, 1) = .0, (2, 21, 1) = .0, (3, 1, 1) = .0, (3, 2, 1) = .0, (3, 3, 1) = .0, (3, 4, 1) = .0, (3, 5, 1) = .0, (3, 6, 1) = .0, (3, 7, 1) = .0, (3, 8, 1) = .0, (3, 9, 1) = .0, (3, 10, 1) = .0, (3, 11, 1) = .0, (3, 12, 1) = .0, (3, 13, 1) = .0, (3, 14, 1) = .0, (3, 15, 1) = .0, (3, 16, 1) = .0, (3, 17, 1) = .0, (3, 18, 1) = .0, (3, 19, 1) = .0, (3, 20, 1) = .0, (3, 21, 1) = .0}, datatype = float[8], order = C_order), ( "totalwidth" ) = 7, ( "rightwidth" ) = 0, ( "solmat_i2" ) = 0, ( "minspcpoints" ) = 4, ( "erroraccum" ) = true, ( "eqndep" ) = [1], ( "errorest" ) = false, ( "banded" ) = true, ( "solspace" ) = Vector(21, {(1) = .0, (2) = .0, (3) = .0, (4) = .0, (5) = .0, (6) = .0, (7) = .0, (8) = .0, (9) = .0, (10) = .0, (11) = .0, (12) = .0, (13) = .0, (14) = .0, (15) = .0, (16) = .0, (17) = .0, (18) = .0, (19) = .0, (20) = .0, (21) = 1.0}, datatype = float[8]), ( "solmat_i1" ) = 0, ( "timeadaptive" ) = false, ( "spacestep" ) = 0.500000000000000e-1, ( "initialized" ) = false, ( "vectorhf" ) = true, ( "linear" ) = true, ( "spacevar" ) = R, ( "periodic" ) = false, ( "spaceadaptive" ) = false, ( "mixed" ) = false, ( "inputargs" ) = [(-(3/2)*R^4+3/2)*(diff(theta(x, R), x))-4*(diff(theta(x, R), R)+R*(diff(diff(theta(x, R), R), R)))/R, {theta(0, R) = 1, theta(x, 1) = 0, (D[2](theta))(x, 0) = 0}], ( "bandwidth" ) = [1, 3], ( "PDEs" ) = [(-(3/2)*R^4+3/2)*(diff(theta(x, R), x))-4*(diff(theta(x, R), R)+R*(diff(diff(theta(x, R), R), R)))/R], ( "leftwidth" ) = 1 ] ); if xv = "left" then return INFO["solspace"][1] elif xv = "right" then return INFO["solspace"][INFO["spacepts"]] elif tv = "start" then return INFO["t0"] elif not (type(tv, 'numeric') and type(xv, 'numeric')) then error "non-numeric input" end if; if xv < INFO["solspace"][1] or INFO["solspace"][INFO["spacepts"]] < xv then error "requested %1 value must be in the range %2..%3", INFO["spacevar"], INFO["solspace"][1], INFO["solspace"][INFO["spacepts"]] end if; dary := Vector(3, {(1) = .0, (2) = .0, (3) = .0}, datatype = float[8]); daryt := 0; daryx := 0; dvars := []; errest := false; nd := nops(INFO["depvars"]); if dary[nd+1] <> tv then try `pdsolve/numeric/evolve_solution`(INFO, tv) catch: msg := StringTools:-FormatMessage(lastexception[2 .. -1]); if tv < INFO["t0"] then error cat("unable to compute solution for %1<%2:
", msg), INFO["timevar"], INFO["failtime"] else error cat("unable to compute solution for %1>%2:
", msg), INFO["timevar"], INFO["failtime"] end if end try end if; if dary[nd+1] <> tv or dary[nd+2] <> xv then `pdsolve/interp2dto0d`(3, INFO["soltimes"], INFO["spacepts"], INFO["solspace"], nops(INFO["depvars"]), INFO["solution"], true, tv, xv, dary); if errest then `pdsolve/interp2dto0d`(3, INFO["soltimes"], INFO["spacepts"], INFO["err_t"], nops(INFO["depvars"]), INFO["solution"], true, tv, xv, daryt); `pdsolve/interp2dto0d`(3, INFO["soltimes"], INFO["spacepts"], INFO["err_x"], nops(INFO["depvars"]), INFO["solution"], true, tv, xv, daryx) end if end if; dary[nd+1] := tv; dary[nd+2] := xv; if dvars = [] then [seq(dary[i], i = 1 .. INFO["norigdepvars"])] else vals := NULL; for i to nops(dvars) do j := eval(dvars[i]); try if errest then vals := vals, evalhf(j(tv, xv, dary, daryt, daryx)) else vals := vals, evalhf(j(tv, xv, dary)) end if catch: userinfo(5, `pdsolve/numeric`, `evalhf failure`); try if errest then vals := vals, j(tv, xv, dary, daryt, daryx) else vals := vals, j(tv, xv, dary) end if catch: vals := vals, undefined end try end try end do; [vals] end if end proc; stype := "1st"; if nargs = 1 then if args[1] = "left" then return solnproc(0, "left") elif args[1] = "right" then return solnproc(0, "right") elif args[1] = "start" then return solnproc("start", 0) else error "too few arguments to solution procedure" end if elif nargs = 2 then if stype = "1st" then tv := evalf(args[1]); xv := evalf(args[2]) else tv := evalf(args[2]); xv := evalf(args[1]) end if; if not (type(tv, 'numeric') and type(xv, 'numeric')) then if procname <> unknown then return ('procname')(args[1 .. nargs]) else ndsol := pointto(solnproc("soln_procedures")[1]); return ('ndsol')(args[1 .. nargs]) end if end if else error "incorrect arguments to solution procedure" end if; vals := solnproc(tv, xv); vals[1] end proc

(6)

NULL

gg := U(x, 1):

NULL

thm := int(U(x, R)*Uu, R = 0 .. 1):

 

 

NULL

 

Download U(R)_numériqueg2.mw

I have a nonlinear function Q(a,b,c,d,x,y) and I'd like to get the optimum (x*,y*) for different values of (a,b,c,d). The usual sintax:

NLPSolve(Q(10, 1, 5, 2, x,y), x= 0 .. 50, y = 0 .. 50, initialpoint = {x = 2,y= .5}, assume = nonnegative) does not work when Q contains numerical integration, that is evalf (Int). I have no problem with the definite integral evalf(int). The problem is that most of the cases required numerical integration so I need the former expression.

I'd appreciate very much if someone could help me.

vz := 2*(-eta^2+1);

D_im := .22;

r0 := 1;

pde := diff(vz*Y(eta, z), z)-D_im*((diff(eta*(diff(Y(eta, z), eta)), eta))/eta+diff(Y(eta, z), `$`(z, 2)))/r0 = 0;

pde := expand(%);

ibc := [Y(1, z) = 0, (D[1](Y))(0, z) = 0, Y(eta, 0) = 1, (D[2](Y))(eta, 0) = 0];

sol := pdsolve(pde, ibc, numeric, time = z, range = 0 .. 1);

pds := sol:-value(z = 0, output = listprocedure);

sol:-plot(z = 0.1e-3, numpoints = 50, color = blue, view = 0 .. 1)

So I was trying to solve this conservation equation for the radial coordinate eta and the z coordinate being treated as time. The flow is in z direction. Now unfortunately it is diverging. Not sure why though. What am I doing wrong?

Hello, My problem is as following:

 

I have tried 2 options for solving the problem below, trying to plot the behaviour of a system to a predetermined function.

First I tried to use dsolve as usual:

restart; with(plots); C := setcolors(); with(LinearAlgebra);
eq1 := Force = Mass*(diff(y(t), `$`(t, 2)));
formula1 := 2.6*BodyWeight*abs(sin(4*Pi*t));
2.6 BodyWeight |sin(4 Pi t)|
BodyWeight := 80*9.81;
plot(formula1, t = 0 .. 2);


eq2 := formula1-SpringConstant*(diff(y(t), t)) = Mass*(diff(y(t), `$`(t, 2)));
/ d \ / d /
2040.480 |sin(4 Pi t)| - SpringConstant |--- y(t)| = Mass |--- |
\ dt / \ dt \

d \\
--- y(t)||
dt //
Mass := .200;
Springt := 200;
200
SpringConstant := Youngsmodulus*Surface/DeltaLength;
DeltaLength := 0.2e-1-y(t);
Surface := .15;
Youngsmodulus := 1600*10^6-20*t^2;
eq2;
/ 2 \ / d \
0.15 \-20 t + 1600000000/ |--- y(t)|
\ dt /
2040.480 |sin(4 Pi t)| - ------------------------------------- =
0.02 - y(t)

/ d / d \\
0.200 |--- |--- y(t)||
\ dt \ dt //

incs := y(0) = 0, (D(y))(0) = 0;
eq4 := dsolve({eq2, incs});
Warning: System is inconsistent

 

Second, I tried using a numerical solving, with maxfun.


restart; with(plots); C := setcolors(); with(LinearAlgebra);
eq1 := Force = Mass*(diff(y(t), `$`(t, 2)));
formula1 := 2.6*BodyWeight*abs(sin(4*Pi*t));
2.6 BodyWeight |sin(4 Pi t)|
BodyWeight := 80*9.81;
plot(formula1, t = 0 .. 2);


eq2 := formula1-SpringConstant*(diff(y(t), t)) = Mass*(diff(y(t), `$`(t, 2)));
/ d \ / d /
2040.480 |sin(4 Pi t)| - SpringConstant |--- y(t)| = Mass |--- |
\ dt / \ dt \

d \\
--- y(t)||
dt //
Mass := .200;
Springt := 200;
200
SpringConstant := Youngsmodulus*Surface/DeltaLength;
DeltaLength := 0.2e-1-y(t);
Surface := .15;
Youngsmodulus := 1600*10^6-20*t^2;
eq2;
/ 2 \ / d \
0.15 \-20 t + 1600000000/ |--- y(t)|
\ dt /
2040.480 |sin(4 Pi t)| - ------------------------------------- =
0.02 - y(t)

/ d / d \\
0.200 |--- |--- y(t)||
\ dt \ dt //

incs := y(0) = 0, (D(y))(0) = 0;
eq4 := dsolve({eq2, incs}, y(t), type = numeric, output = listprocedure, maxfun = 10^7);
[
[t = proc(t) ... end;, y(t) = proc(t) ... end;,
[

d ]
--- y(t) = proc(t) ... end;]
dt ]

test := rhs(eq4[2]);
proc(t) ... end;

This one does plot, but no further than 0.2*10^-6. I have tried compiling the data, but this has not worked yet.

 

Does anyone know a way to work around such a problem. Is it possible to plot the equation using a for loop? If yes, how?

 

 

1 2 3 4 5 6 7 Last Page 1 of 22