Rouben Rostamian

MaplePrimes Activity


These are answers submitted by Rouben Rostamian

Here is the desired volume:

See the attached worksheet, volume-of-a-spherical-wedge-ver2.mw, for the derivation.

I have replaced the H1 and H2 with a and b to simplify the presentation.

You have

Eq:= (G-0.043)(Lc/2)+(-G)(Ls/2)+G*Lcon-0.043*Lc=14;

That should be

Eq:= (G-0.043)*(Lc/2)+(-G)*(Ls/2)+G*Lcon-0.043*Lc=14;

That will fix the problem and you will get a solution from fsolve().

However, your equation has two solutions.  There are at least two ways to get the second solution:

  1.  You can tell fsolve() the interval to search for the second solution.  You can get a rough idea of where that solution lies by plotting the graph of Eq - 14.
  2. Since your equation is particularly simple, you may find the two solutions at once with solve(), as in solve(Eq), instead of fsolve(Eq).

A side comment: The subject line of your inquiry refers to simultaneous equations.  But Eq is a single equation, there is nothing simultaneous about it.

 

I see no reason for such an error message in the code that you have shown.  When I paste that code to a Maple worksheet, it executes without a flaw.  Perhaps there is something else in your worksheet that is interfering.  Consider uploading your troublesome worksheet for comments.  Use the large green up-arrow in the reply screen to upload

Let

f := (x,t) ->  x + sin(x) +  ln(t);

Then consider f(x,t) = 0 as defining x(t) and a function of t.  Note that x(1)=0.

We have f(x(t),t) = 0,  therefore the derivative of f(x(t),t) with respect to t is zero.  This gives us a differential equation, as you have noted:

diff(f(x(t),t),t);
solve(%, diff(x(t),t));
de := diff(x(t),t) = %;

Solve the differential equation

dsol := dsolve([de, x(1)=0], x(t), numeric);

Now define

F := x -> x + x^2;

and plot F(x(t)):

plots:-odeplot(dsol, [t, F(x(t))], t=0..10);

Here is the solution to Case (vi) in the picture that you have posted.

restart;
doit := proc(a,b,c)  # Case (vi): assumes a < b < c
   plot([seq([x + a*t, t, t=0..4], x = -6..-1, 0.3)], color="Red"),
   plot([seq([x + b*t, t, t=0..4], x = -1..1, 0.3)], color="Green"),
   plot([seq([x + c*t, t, t=0..4], x = 1..6, 0.3)], color="Blue"),
   plot([seq([-1 + s*t, t, t=0..4], s=a..b, 0.15)], color="Magenta"),
   plot([seq([ 1 + s*t, t, t=0..4], s=b..c, 0.15)], color="Cyan")
end proc:
plots:-display([doit(0,1,2)], axes=boxed, scaling=constrained, view=[-4..4, 0..3], labels=[x, t]);
 

The characteristic lines for Case (vi) with a=0, b=1, c=2

Another sample:

plots:-display([doit(-1,0.5,2)], axes=boxed, scaling=constrained, view=[-4..4, 0..3], labels=[x, t]);

The characteristic lines for Case (vi) with a=-1, b=0.5, c=2

 

Your code is correct other than

h=evalf((xn-x0)/n);

should be

h:=evalf((xn-x0)/n);

Aside: To verify the correctnest of your integrator, try it on a cubic polynomial.  Simpson 3/8 should produce the exact answer. (Yours does.)  Remove the evalf()s to see that clearly.  Thus:

expand(simp38(x -> x^3, a, b, 9));

results in 

b^4/4 - a^4/4

as expected.

 

 

If you open the file in Maple, you will get everything that the user has typed.  You don't need a special command for that.  But perhaps I have misunderstood your question.

Your complexplot3d() is not producing a plot because complexpplot3d() is defined in the plots package, but you are loading the plots package too late.  Move the with(plots) line further up in the worksheet.

The comlexplot() at the worksheet's end does produce a graph.  To help you see the graph, change the plotting command to

complexplot({a, b}, numpoints = 100, color = red, thickness = 3, filled = true, title = "Stability Region");

What you have written in your message focuses on Maple.  That's not quite right.  I suggest that you set Maple aside for now, and make an effort to understand what Euler's method is supposed to do.  Remember: Euler lived about 250 years ago and he did not have Maple.

Only after you have understood what the method does, attempt to ask Maple to do that thing for you.

The answer is approximately 1/2015!.  Here is a sketch.

Let's write p(x) for mul(x+j, j = 0 .. 2015).  We are looking for the positive root of p(x)=1.  We have:

p(x) = (x^2015 + ... + 2015!) . x

We see that p'(0) = 2015!, thererefore the graph of p rises quite steeply(!) from p(0)=0 to p(t)=1, where t is the root you are looking for.  Approximating the graph by a straight line, we get t = 1/2015!  We note that

> 1.0/2015!;
                                    -5786
                      8.667800428 10     

which agrees withKitonum's alternative computation.

I don't see why you would want to modify the standard Laplace transform code since you can reduce your integrals to the usual Laplace transform through a change of variables.

The integral int(f(u*t)*exp(-s*t), t=0..infinity) after the change of variables u*t=tau takes the form

(1/u)*int(f(tau)*exp(-s/u*tau), tau=0..infinity),

which is the same thing as

(1/u)*laplace(f(tau), tau, s/u).

Similarly, the integral int(f(u*t)*exp(-t), t=0..infinity) changes over to

(1/u)*laplace(f(tau),tau,1/u).

Whether you want to to do these changes of variables or not depend on your ultimate goal.  In practice, both of those integrals may be computed without resort to any changes of variables:

int(f(u*t)*exp(-s*t), t=0..infinity) = laplace(f(u*t), t, s);

int(f(u*t)*exp(-t), t=0..infinity) = laplace(f(u*t), t, 1);

I don't quite understand what it is that you want to do.  The graphs you have shown are those of y versus psi.  You don't need the velocity V or the DEtools package to produce those.  You will obtain the first of your three diagrams if you replace your plot() command with

R := 0.2;
P := seq(
       [psi, y, y=0..1],
     J in As);
plot([P],
  color=[green,red,blue],
  linestyle=[solid,dash,dot],
  thickness=3,
  axes=boxed,
  labels=['psi','y'],
  legend=[seq(J=j, j in As)],
  title=sprintf("y vs psi at R = %g", R)
 );

Change R=2.0 to 5.0 or 6.5 to obtain the other two diagrams.  By the way, the three diagrams that you have shown are identical.   They would be different  if they were done correctly.

 

Introduce velocities as new variables and convert your second order differential equations to first order.  Then you will have a system of four first order differential equations to solve, and events on velocities may be specified in a straightforward way.  See the attached worksheet.  Here is what we get:

Worksheet here: mw.mw

This plots p versus w:

w := 1560*Pi*(9*z-(1/3)*(z-2.7)^3*((1/3)*2.7^3));

plots[odeplot](Sol, [w, p(z)], z = 0 .. 26, color = green);

As acer has noted, this may be a bug.  Or perhaps this is due to a not-so-well documented requirement in Eigenvectors().  I keep a copy of the old Maple 11 for testing purposes.  Running your problem in Maple 11 gives an error message which provides a clue:

> LinearAlgebra:-Eigenvectors(A);

Error, (in LinearAlgebra:-LA_Main:-Eigenvectors) expecting either Matrices of rationals, rational functions, radical functions, algebraic numbers, or algebraic functions, or Matrices of complex(numeric) values

"Expecting rationals"?  So, we change your 0.991 and 0.997 to 991/1000 and 997/1000, and voila, we get an answer! This works  equally well on Maple 11 and Maple 18.

 

   

 

First 52 53 54 55 56 57 58 Page 54 of 58