Carl Love

Carl Love

28055 Reputation

25 Badges

12 years, 357 days
Himself
Wayland, Massachusetts, United States
My name was formerly Carl Devore.

MaplePrimes Activity


These are answers submitted by Carl Love

When the number variables being solved for is less than the number of equations you should generally use eliminate instead of solve.

eliminate(sys, drivers);

Note that eliminate always returns pairs of sets: The first set is the elimination equations; the second set is the original equations with the appropriate variables eliminated via the substitutions from the first set.

In the pdsolve command, use the options timestep and spacestep. For example,

pds := pdsolve(pde, IBC, numeric, timestep=.01, spacestep= 1);

I cannot vouch for the accuracy of the plots that you get with this, but you will get the plots instead of error messages. I tried smaller values, but the left sides of the plots still have that sharp drop. But note that the boundary condition u(t,-50) = 0 necessitates that the left sides go to 0, so maybe those plots are accurate.

Here's how to make the input the differential equation itself. The procedure preprocessODE takes the ODE (given in the standard Maple form accepted by dsolve) and returns the function f that is passed to eul. I made a slight modification to eul so that it calls preprocessODE. (I am not sure if this will work in Maple 5; the error command syntax might be different. Let me know.)

preprocessODE:= proc(eq::`=`)
local d, y, x, f;
    d:= indets(eq, specfunc(anything, diff));
    if nops(d) = 0 then error "No derivative found." end if;
    if nops(d) > 1 then error "More than 1 derivative found." end if;
    d:= d[];
    y:= op([1,0], d);  x:= op([1,1], d);
    f:= {solve(eq, d)};
    if nops(f) <> 1 then error "Cannot isolate derivative." end if;
    unapply(subs(y(x)= y, f[]), [x,y])
end proc:
    
eul:=proc(ODE::`=`,h,x0,y0,xn)
  local f,no_points,x1,y1,i;
  f:= preprocessODE(ODE);
  no_points:=round(evalf((xn-x0)/h));
  x1:=x0;
  y1:=y0;
 
  for i from 1 to no_points do
      y1:=y1+evalf(h*f(x1,y1));
      x1:=x1+h
  end do;
  y1
end proc:
       

Example of use:

eul(diff(y(x), x) = x^2*y(x)^3, .01, 0, 1, 1);

1.68181553277060

You cannot use indexed variables where the index or the base are the same as the variable of integration. So you need to change lambda[v], lambda[t], and t[c]. I changed them to L__v, L__t, and t__c  (which print as subscripted variables in Maple 17 & 18). Also, you need to assume that the lambdas are positive.

int(
     int(
          L__v*L__t*exp(-L__v*v-L__t*t),
          v = (1/2)*(q[p]+q[p]*t__c*t+2*S[di]*h*t)/(h*t) .. infinity
     ), t = 0 .. infinity
) assuming L__v > 0, L__t > 0;

 

 

The command to plot the cylinder (which is a surface, not a curve) is plot3d, not spacecurve. You may keep the rest of the command as is. However, I do recommend that you omit the options numpoints= 4000, scaling= constrained, color= blue, and thickness= 2 at least until you get the actual plot. The axis is a curve, and you may use spacecurve to plot that (and the thickness and color options would be more a appropriate for the axis).

Please attach your file as Maple worksheet instead of a pdf so that I can work with it. I believe that option is working again. (I've seen a few worksheets posted today.)

Here's an example using fsolve and dsolve's option output= listprocedure:

Sol:= dsolve(
     {
          diff(a[1](t),t$2) = -a[1](t), diff(a[2](t),t$2) = -a[2](t),
          a[1](0)=1, a[2](0)=0, D(a[1])(0)=0, D(a[2])(0)=1
     }, numeric, output= listprocedure
):
A[1]:= eval(a[1](t), Sol): #extract procedure for a[1]
A[2]:= eval(a[2](t), Sol):
fsolve(A[1], 1..infinity);
                       
10.9955744378927
fsolve(A[2], 1..infinity);
                       
3.14159269513977

It is also possible to locate the zeros with the events option to dsolve. I'll leave it to someone else to give an example of that.

Edit: Square brackets after D changed to parentheses.

Is this what you're looking for?

solve({0<x, x<1, x^2+y^2 - 1 > 0}, x, parametric);

If that's not what you want, could you give a specific solved example of what you want?

Your whole procedure could be replaced by the following:

ff:= (N::posint)-> Matrix([[-2 $ 2*N+1], [0 $ 2*N+1] $ 2*N-1, [-3 $ 2*N+1]]);

I don't know if correcting this will solve your problem, but you have an error in your second boundary condition. You have u(y,0) = 0, but it should be u(0,y) = 0.

You need to do

restart;

to clear the previous values of the variables. After that, Maple still will not be able to complete the definite integration (even with assumptions n::posint, h > 0); I don't know why. But you can do this:

restart:
eta:= 1000:  B:= 5/2:  #Convert 2.5 to an exact value.
F:= int(B*eta^(-B)*t^(B-1)*exp(-(t/eta)^B)*(t-n*h), t);
eval(F, t= (n+1)*h) - eval(F, t= n*h); #Apply Fund Thm of Calc.

Note that I assumed that you meant for the lower limit of integration to be n*h rather than n.

 

The command plots:-display is used to combine multiple plots into a single plot. Plots can be assigned to variables, as in

P1:= plot([Y||(1..3)], 0..10):
P2:= plot([YP||(1..3)], 0..10):

plots:-display([P1,P2]);

The command plottools:-arrow can be used to make a plot of an arrow, which can then be combined with other plots using display.

See ?plottools,arrow and ?plots,display .

Add option method= modifiednewton. Note that the Optimization commands only find local extrema.The Pi has nothing to do with your problem.

Assuming that your count of boundary conditions is correct, a common cause for this error is having a parameter whose numerical value has not been specified. I can't diagnose it further unless you post the code.

restart:
f:= 3-x+2*x^2:  g:= 4+k*x-x^2:
solve(int(f*g, x= -1..1), k);

If you want to get fancy with the inner product notation, you could do

restart:
f:= x-> 3-x+2*x^2:  g:= x-> 4+k*x-x^2:
local `<,>`:= (f,g)-> int(f*g, -1..1):
sqrt(<f,f>);

solve(<f,g>, k);

 

The following does not prove the divergence, but the wildly changing values makes convergence highly suspect:

restart:
Digits:= 15:
f:= z-> exp(-z^2*sin(z)^2):
F:= d-> int(f(z), z= 0..infinity, numeric, epsilon= 10.^(-d)):
'F(d)' $ d= 6..9;



These wildly changing values also make convergence highly suspect:

restart:
Digits:= 10:
f:= z-> exp(-z^2*sin(z)^2):
int(f(z), z= 0..infinity, numeric);

2.835068335

Digits:= 15:
int(f(z), z= 0..infinity, numeric);

2.68030993194038

First 308 309 310 311 312 313 314 Last Page 310 of 395