Items tagged with numeric numeric Tagged Items Feed

This is the code


Let me explain it.

I am sure that the mistakes must be about the expresstion of the I1(t) and I2(t). Actually if you delete I1(t) and I2(t) , the whole code works and get the picture at the bottom. 

What I want is to put the expresstion of the I1(t) and I2(t) into 'sol:=dsolve...' and 'plots...' to get the picture of I1(t) and I2(t) with respect to t. Before the t* which subject to Phi(t*)=0 (The blue line in the picture at the bottom is Phi) I want I1(t) and after t* I want I2(t).

I1(t) = (int(sqrt(2*(H(t)+omega*cos(q(t)))), q(t) = q(t)-2*Pi .. q(t), numeric))/Pi.    what I want of this experesstion is to get  'int(sqrt(2*(H(t)+omega*cos(q(t)))' from  'q(t)-2*Pi' to 'q(t)' by numeric method.This q(t) is the solution of the ODE sys.

For example(the number I used is not true,just for example) , at the point t=20, q(t)=30-2*Pi.

so I1(t)= (int(sqrt(2*(H(t)+omega*cos(x))), x = 30-2*Pi .. 30, numeric))/Pi.The I2(t) I want is similar to I1(t).


How can I solve it?

when the variable is in limits of integration?

there exist an error : 

Error, (in int) unable to compute a numeric answer for symbolic limits, q = -3.141592654+arccos(11./(1.+.1000000000*t)) .. 3.141592654-1.*arccos(11./(1.+.1000000000*t))

How to make this code work?

I have in mind


My unsuccessful trials are

evalf(Int(cos(log(x)/x)/x, x = 0 .. 1, method = _d01akc));

Error, (in evalf/int) NE_QUAD_MAX_SUBDIV:
  The maximum number of subdivisions has
  been reached: max_num_subint = 500


evalf(Int(cos(log(x)/x)/x, x = 0 .. 1, method = _d01akc, methodoptions = [maxintervals = 2000]));

Error, (in evalf/int) NE_QUAD_MAX_SUBDIV:
  The maximum number of subdivisions has
  been reached: max_num_subint = 2000


evalf(Int(cos(log(x)/x)/x, x = 0 .. 1, method = _d01ajc));

Error, (in evalf/int) NE_QUAD_NO_CONV:
  The integral is probably divergent or slowly convergent.

Mathematica cracks it without any options, producing 0.323367.



with(DETools, diff_table);

kB := 0.138064852e-22;

R := 287.058;

T := 293;

p := 101325;

rho := 0.1e-2*p/(R*T);

vr := diff_table(v_r(r, z));

vz := diff_table(v_z(r, z));

eq_r := 0 = 0;

eq_p := (vr[z]-vz[r])*vr[] = (vr[]*(vr[r, z]-vz[r, r])+vz[]*(vr[z, z]-vz[z, r]))*r;

eq_z := 0 = 0;

eq_m := r*vr[r]+r*vz[z]+vr[] = 0;

pde := {eq_m, eq_p};

IBC := {v_r(1, z) = 0, v_r(r, 0) = 0, v_z(1, z) = 0, v_z(r, 0) = r^2-1};

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


what am I doing wrong?

it's telling me: Error, (in pdsolve/numeric/par_hyp) Incorrect number of boundary conditions, expected 3, got 2
but i did just as in the example :-/

and plot  function I? This I is the area which I wrote at the paper.

Could you give me the code which can be used to solve the ODE by numerical method and plot I with respect to t?

I think I have write down everything clearly but if you feel confused please ask me.

I am eager to know the code. Thanks very much!

Perhaps the question is trivial, but I could not find the solution.

I am solving numerically an ODE (e.g., a simple harmonic oscillator) with the righthand side that contains a random part. For example it is

eq:=diff(a(t),t,t) + a(t) = (1+0.01*R(t))*cos(5*t);

where R(t) is a random function of t. How can I make such a function?

The naive attempt: 

eq:=diff(a(t),t,t) + a(t) = 3.*(1+0.1*rand()/1000000000000.)*cos(5.*t);

gave me a fixed (while random) value, e.g.


But, I need this coefficeint to be random each time step.

Any suggestions are very welcome!




Are there any commands in maple that will help me find a suitable function that approximates the numerical solution of:

  PDE := diff(v(x, t), t) = diff(v(x, t), x, x);
  JACOBIINTEGRAL := int(JacobiTheta3(0, exp(-Pi^2*s))*v(1, t-s)^4, s = 0 .. t);
  IBC:= D[1](v)(0,t)=0,
        D[1](v)(1,t)=-0.000065*v(1, t)^4,
# For x=0..1, t=0..1, the solution varies only very slowly
# so I have increased the timestep/spacestep, just to speed
# up results generation for diagnostic purposes
  pds := pdsolve( PDE, [IBC], numeric, time = t, range = 0 .. 1,
                  spacestep = 0.1e-1, timestep = 0.1e-1,

diff(v(x, t), t) = diff(diff(v(x, t), x), x)


int(JacobiTheta3(0, exp(-Pi^2*s))*v(1, t-s)^4, s = 0 .. t)


(D[1](v))(0, t) = 0, (D[1](v))(1, t) = -0.65e-4*v(1, t)^4, v(x, 0) = 1




# Plot the solution over the ranges x=0..1,
# time=0..1. Not a lot happens!
  pds:-plot(x=1, t=0..1);


# Plot the estimated error over the ranges x=0..1,
# time=0..1
  pds:-plot( err(v(x,t)), x=1,t=0..1);


# Get some numerical solution values
  pVal:=pds:-value(v(x,t), output=procedurelist):
  for k from 0 by 0.1 to 1 do
      pVal(1, k)[2], pVal(1, k)[3];


t = 0., v(x, t) = Float(undefined)


t = .1, v(x, t) = .999977377613528229


t = .2, v(x, t) = .999967577518313666


t = .3, v(x, t) = .999959874331053822


t = .4, v(x, t) = .999952927885405241


t = .5, v(x, t) = .999946262964885979


t = .6, v(x, t) = .999939702966688881


t = .7, v(x, t) = .999933182128311282


t = .8, v(x, t) = .999926675964661227


t = .9, v(x, t) = .999920175361791563


t = 1.0, v(x, t) = .999913676928735229





Download PDEprob2_(2).mw


I am refering to the first graph, is there a way in maple to find an explicit suitable approximating function?

I.e, I want the function to have the same first graph obviously, it seems like addition of exponent and a line function, I tried plotting exp(-t)-0.3*t, it doesn't look like it approximates it very well. Any suggestion on how to implement this task in maple?



 of solution to plot a function?


sys := {diff(x(t), t) = y(t), diff(y(t), t) = -(1+(1/100)*t)^2*x(t), x(0) = -1, y(0) = 2};

sys := {diff(x(t), t) = y(t), diff(y(t), t) = -(1+(1/100)*t)^2*x(t), x(0) = -1, y(0) = 2}

dsn := dsolve(sys, numeric);

dsn:=proc(x_rkf45) ... end proc


My question is how to use dsn to plot functions E(t) and I(t) with respect to t in a given interval [t1,t2].

E(t)=0.5*(diff(x(t),t)^2+w^2*x(t)^2), w=1+t/100.  I(t)=E(t)/w.


Many thanks to everyone helping me : )

Hello altogether,

I want to plot the numerical result of an ODE, which seems to be pretty simple at first sight, but the difficulty is that the boundaries are depending on the solution.

The following pseudo-code describes what I want to have, but it doesn't work. This code fills the RAM pretty fast and you will have to kill the process.

Is it possible to calculate a solution to this problem numerically (or even analytically) and if yes, how?

Since I am new here, I am sorry for any bad-to-read maple code or any noob errors I have made. I would be very thankful, for any response and help.



Dear Collagues

I wrote a code to solve a system of ode numerically. When the equations have been solved, I want to find the value of MTE which is in integral form. I check and see that the value of MTE is calculated wrong. I dont know why?

I calculated the values of MTE twice (another variable is masst below res1). Two different answers!!

I would be most grateful if you could help me.

Thank you

Here is my code







eq1:=(a[mu1]*(diff(phi(eta), eta))+2*b[mu1]*phi(eta)*(diff(phi(eta), eta)))*(diff(u(eta), eta))+(1+a[mu1]*phi(eta)+b[mu1]*phi(eta)^2)*(diff(u(eta), eta, eta))+rhost(eta)-Ha^2*u(eta);




res1 := dsolve(subs(NBT=1,gama1=0.3,phiv=0.06,{eq1,eq2,eq3,u(0)=lambda*D(u)(0),D(u)(1)=0,T(0)=0,phi(1)=phiv,T(1)=1}), numeric,method=bvp[midrich],output=listprocedure,continuation=cc):

masst:=evalf(int((1-G0(eta)+G0(eta)*rhop/rhobf)*G1(eta), eta = 0..1));

for iNBT from 1 to 3 do

for iphi from 1 to 3 do

for igamma from 1 to 3 do

res := dsolve(subs(NBT=ANBT(iNBT),gama1=Agama1(igamma),phiv=Aphiv(iphi),{eq1,eq2,eq3,u(0)=lambda*D(u)(0),D(u)(1)=0,T(0)=0,phi(1)=phiv,T(1)=1}), numeric,method=bvp[midrich],output=listprocedure,continuation=cc):
MTE[iNBT,iphi,igamma]:=evalf(Int((1-F0(eta)+F0(eta)*rhop/rhobf)*F1(eta), eta = 0..1)):

end do;

end do;

end do;

for iNBT from 1 to 3 do

for iphi from 1 to 3 do

for igamma from 1 to 3 do

print (convert(HTC[iNBT,iphi,igamma],float,9));

end do;

end do;

end do;
for iNBT from 1 to 3 do

for iphi from 1 to 3 do

for igamma from 1 to 3 do

print (convert(MTE[iNBT,iphi,igamma],float,9));

end do;

end do;

end do;



Hi everyone.

I'm trying to solve the following PDE


but I'm getting this error:

Error, (in pdsolve/numeric/process_IBCs) initial/boundary conditions can only contain derivatives which are normal to the boundary, got (D[1, 1](w))(x, 0)

The PDE represents the bending of a thin plate.

See File:

Hi, I am working on an assignment and I am having difficulties with this Eurler method since I am working with vector.

The assignment:

If the aircraft's engine transmits power P to the air, the force of the propellar, Fp, satisfies P=Fpv. Since the force Fp acts in the direction of motion, it can be written in vector form as: (what I have so far)


v:=Vector(2,1,[vx,vy]); velocity vector

nu:=v->Norm(Vector(v),2); (length of v)

Fd:=-1/2*C*rho*A*nu*v # air resistance

Fl:=1/2*Cl*rho*A*nu*v # liftforce in y direction here v=[-vy,vx]

Fp:=(P/nu^2)*v^2 # power from aircraft

All this adds up to:

F_T:=v->Vector(evalm(Fp/nu(v)^2*v+Fd*nu(v)*v+Fl*nu(v)*zip((x,y)->x*y,v([2,1]),Vector([-1,1]))+Vector(2,1,[0,-m*g]))); # I have tried it and it works


so up til here I am good. It is this next part that I can't seam to get.

Write a routine which uses Euler's method to solve the equation of motion numerically for velocity v(t):


Test the program by trying a plane with intitial velocity [v0,0] such that it maintains this speed when the engines power is P=P0

We have been given an Esolve we can use:

Esolve:=proc(f::procedure,h,x0,y0,N::integer) # Calculates y(x0+n*h) for n=1..N given y0=y(x0) and f(x,y)

local n,y,x;

y[0]:=y0; # Start Value

for n from 1 to N do
y[n]:= y[n-1]+h*f(y[n-1],x); # Euler's formula
x:=x+h; # Next x
end do;


end proc;


If anyone could help me understand how it works and how to get this to work with vectors I would be very gratefull.

Hi !

I am trying to solve a pde with initial/boundary conditions, in the numeric mode. It works very well when I provide 3 numerical conditions : 1 initial condition (uniform distribution at first), and two boundary conditions.

Now I want to switch the initial condition to a function of r. It's a polynome I obtained after interpolation of another result. I tested it and the function interp20(r) works. But the pdsolve doesn't seem to evaluate the function, when it comes to start the procedure and pdsolve doesn't return a module as it should, although it doesn't send any message error.

You can see the code following, with the error messages I get.


Is it strictly impossible to use a function as an initial condition ?
Is it just a problem I could solve by converting the function or its result to something else ?
Is float the right type of input ? For example, if I write a:=283.15, is it a float ?
Have you already had similar problems ? How did you solve them ? Where could I find working code examples on this ?
Could I solve this problem with non-uniform initial condition with another Maple function ?

Then you very much for your help !



Bad title, I know. 

Suppose I have solved an ODE numerically, say for a function f(x). Now, I want to solve another ODE numerically, given by

diff(g(x),x) = f(x)*g(x)

I want to solve for g(x), and I only know f(x) from the numerical solution I obtained in the first ODE. How do I implement this in Maple?

I attached a worksheet with my naive attempt with this simple example. Basically I do not know how to set up the differential equation which I would then input to dsolve.


hi i a have a problem with program which attached..please help me

thanks alot

dsys3 := {(diff(w(x), x, x, x, x))/x^2 = .4/(x^3*(1-w(x))^3)+40*(1+.65*(1-w(x)))/(x*(1-w(x))^2), w(0) = 0, w(1) = 0, ((D@@1)(w))(0) = 0, ((D@@1)(w))(1) = 0}

dsol5 := dsolve(dsys3, numeric, output = array([.5]))

Error, (in dsolve/numeric/bvp) system is singular at left endpoint, use midpoint method instead

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