Rouben Rostamian

MaplePrimes Activity


These are answers submitted by Rouben Rostamian

You have terminated almost every command in your worksheet with a colon, in effect depriving yourself from the benefit of Maple's feedback.

Change the colons to semicolons to see what Maple thinks of your input.  When you do that, you will see that the fourth equatiuon in Eqns is not what you would have expected.  Examining the cause, you will find out that you have entered vyH where you should have vyH(t).

Carl has pointed out several syntax errors and other issues with your equations.

I believe there are more problems.  If the equations are to model an oscillator, the x on the right-hand side of the first equation should be z, and the y on the right-hand side of the second equations should be x.  So we change the right-hand sides to:

f := (x,z) -> alpha*(1 - (x^2/a^2 + z^2/b^2))*z + w*a/b*z;
g := (x,z) -> beta*(1-(x^2/a^2 + z^2/b^2))*x + w*b/a*x;

Then the differential equations are:

de1 := diff(x(t),t) = f(x(t),z(t));
de2 := diff(z(t),t) = g(x(t), z(t));

We introduce the parametes

params := alpha=1, beta=1, a=0.4, b=0.2, w=1;

and evaluate the equations accordingly:

sys := eval([de1, de2], [params]);

To get an idea of the region of interest in the x-z plane, we find the system's equilibria:

RealDomain[solve](
    eval([f(x,z)=0, g(x,z)=0], [params])
);

We see that in addition to the origin (0,0), there are four equilibria at x=0, z=+/- 0.35, and x=+/-0.49, z=0.  Thus we set

xmax := 0.6;  zmax := 0.5; tmax := 8;

and generate a reasonable set of initial data:

 ic1 := seq([x(0)=h,z(0)=0], h=-xmax..xmax, 0.1),
      seq([x(0)=h, z(0)=-zmax], h=-xmax..0, 0.05),
      seq([x(0)=h, z(0)=zmax], h=0..xmax, 0.05);

Now we call DEplot to produce a phase portrait:

DEtools[DEplot](sys, [x(t),z(t)], t= 0..tmax, [ic1],
  linecolor=black, thickness=1, stepsize=tmax/(3*48),
  x(t)=-xmax..xmax, z(t)=-zmax..zmax, arrows=none);

Lo and behold, the Starbucks siren emerges!

We may want to add extra orbits to cover the chest area:

ic2 := seq([x(0)=h,z(0)=0], h=-xmax..xmax, 0.1),
      seq([x(0)=h, z(0)=-zmax], h=-xmax..0, 0.05),
      seq([x(0)=h, z(0)=zmax], h=0..xmax, 0.05),
      seq([x(0)=0, z(0)=h], h=-0.34..0.34, 0.04);

DEtools[DEplot](sys, [x(t),z(t)], t= 0..tmax, [ic2],
  linecolor=black, thickness=1, stepsize=tmax/(3*48),
  x(t)=-xmax..xmax, z(t)=-zmax..zmax, arrows=none);

 

I think you are asking for a surface of revolution.  A tubeplot is something entirely different.

A surface of revolution may be ploted in Maple in a very many different ways, the simplest of which is with the help of Student[Calculus1] package, like this:

restart;
with(Student[Calculus1]):
SurfaceOfRevolution(x^2+x-1, x=0..1, axis=vertical, output=plot, axes=boxed);  p1 := %:

We see that the vertical axis goes up to z=1.  You want to continue this straight up, to height z=5.  The graph of that extra part is an ordinary cylinder, which we produce through

plot3d(1, t=0..2*Pi, z=1..5, coords=cylindrical, color="Maroon", style=surface);  p2 := %:

Now p1 and p2 represent thos two plots. We put them together with the help of plots:-display command:

plots:-display([p1,p2], scaling=constrained);

For some reason which I don't understand, this does not do the expected thing.  It shows only the p1 graph  This may be a bug in Maple, or perhaps something else that I am doing wrong.  I will let others comment on that.

Added later:

I just noticed that if we specify a view=... option to the latter command, it displays the composite picture properly:

plots:-display([p1,p2], scaling=constrained, view=-1..6);

The reason that it does not work without the view is that the SurfaceOfRevolution command sets the view option to its own extents, which is view=-1..1 in this case.  That option caries over to the plots:-display command, and therefore the result is that the composite surface is truncated at the z=1 level.  We specify the view=-1..6 option explicitly to override the one set by SurfaceOfRevolution.

The way I see it, this is sort of a minor bug, but perhaps not worth making a big fuss about.

 

You don't need DEtools for this, but you need plots:

restart;
with(plots):
w := 1+t/100;
sys := {diff(x(t), t) = y(t), diff(y(t), t) = -w^2*x(t), x(0) = -1, y(0) = 2};
dsn := dsolve(sys, numeric);
odeplot(dsn, [[t,x(t)],[t,y(t)]], t=0..10);
E := 1/2*(y(t)^2 + w^2*x(t)^2);
odeplot(dsn, [t, E], t=0..10);
J := E/w;  # I is predefined in Maple; use J
odeplot(dsn, [t, J], t=0..10);

I have never felt the need to use animate() whose syntax I find rather cluncky.  I just prepare the sequence of the individual frames first, then use display() to animate it, as in:

with(plots):

frames := seq(
  dualaxisplot(plot(sin(x+A), x=0..5), plot(cos(x+A), x=0..5)),
  A=0..5
):

display([frames], insequence=true);

I couldn't execute your worksheet in Maple 2015.  It ran into errors.

But just looking through it I formed an idea of what it is that you want to do and produced this worksheet: mw.mw

The end result is:

As to your question about the area of the spherical triangle, there are very many neat formulas for it.  See the Wikipedia article:

    https://en.wikipedia.org/wiki/Spherical_trigonometry

where you wil find formulas for the area in terms of the angles a,b,c (these are the angle between the vectors shown in the graphics above) or the angles A,B,C (these are the actual angles at the triangle's vertices measured on the sphere's surface).  Both sets of angles may be calculated easily in terms of dot products and cross products from the given data.

 

If I am not mistaken, the following worksheet does what you are asking.

mw.mw

Here is the graph of the function which is expressed as an integral in your post:

Here is one of many different ways of writing a loop:

for V from 1 to 30 do
   ....
end do;

 

convert(abs(1-b)+abs(1+b), piecewise);

D[1,1](w) means the second derivative of w(x,y) with respect to x.  Then D[1,1](w)(x,0) means the second derivative in the x direction along the bottom edge.  Surely you don't mean that.  If your intent is to specify zero bending moment along the bottom edge, you want D[2,2](w)(x,0)=0.

General advice: If you are having difficulty with a complex problem, see if you can do a simpler problem first.

So forget about your F_T for the moment.  Do you know how to work with your Esolve()?

For that, examine what Esolve() does, and what its parameters are.

We see that Esolve() is designed to solve a system of differential equations of the form dy/dx = f(y,x) for the unknown vector-valued function y(x).

The intial condition is y(x0) = y0.

The solution is produced by taking N steps of stepsize h.

You are interested in solving a problem of the type dv/dt = F(v,t).  Thus, your v and t correspond to Esolve()'s y and x, respectively.

Let's give it a try with a function f which takes a vector <v[1], v[2]> as argument and returns the vector <-v[2],v[1]>.  Thus, in Maple we set:

f := v -> <-v[2], v[1]>;

For the initial value pick v(0) = <1,0>.  Now solve with stepsize h = 0.1 and take 25 steps:

sol := Esolve(f, 0.1, 0, <1,0>, 25);

Plot the solution, that is, v(t) versus t:

plots:-pointplot([seq(sol[i], i=0..25)], scaling=constrained, color=red, connect=true);

See if you can take it up from there.

This solution is obtained by eyeballing the figure, not through any algorithm:

Replace 1 with 1.0001 in the square root:

plots:-shadebetween(2,1/sqrt(1.0001-x^2), x=0..sqrt(3)/2);

That looks like a bug to me.  I suggest that you report it.

In the meanwhile, you may get the desired plot by "cheating" as follows.  Exchange the x and y axes temporarilty, and plot whatever you want "sideways":

plot(
  [[2, x, x=0..sqrt(3)/2],
   [1/sqrt(1-x^2), x, x=0..sqrt(3)/2]
  ], color=[red,blue], thickness=3, filled=[color=yellow]);    p := %:

Then apply plottools:-transform() to restore the x and y axes to their standard positions:

plottools:-transform((x,y)->[y,x])(p);

Producing a plot of that motion is quite trivial in Maple, and it is quite likely that that several people will show you how.  However, since you are taking a course in advanced dynamics, I think you will benefit from figuring this out in your head.  I suggest that you begin with simpler problems.

  1. What does the motion

    r(t) = cos(ωt) i + sin(ωt) j

    look like?  If you don't see that, then I am afraid that this is not quite the right time for advanced dynamics.

  2. What does the motion

    r(t) = (a cos(ωt))i + (b sin(ωt))j

    look like?  This is only slightly different from the previous one.

  3. Now you should be ready to answer your original question.  None of these requires help from Maple.
First 50 51 52 53 54 55 56 Page 52 of 58