Rouben Rostamian

MaplePrimes Activity


These are replies submitted by Rouben Rostamian

Here are a couple of things to think about:

  1. Your data set consists of 42 points.  Simpon's rule of integration requires an odd number of data points. Why are you focused on Simpson's rule in this case?
  2. The factor sin(x) looks odd.  The "angle" x goes from 1 to 42.  What are the units of measurement here? Does sin(42) mean the sine of 42 degrees, 42 grads, or 42 radians?  Unless you qualify this, Maple will assume that you mean 42 radians.  Is that what you intend?

 

@mary120 The "trick", if it can be called so, is to plot the system's vector field first, without worrying about the solution curves.  Thus, you will do:
DEtools:-DEplot({de1,de2}, [x(t),y(t)], t=0..1, x=-4..4, y=-4..4);
and get:

From that picture you can see that if you start anywhere on the positive x axis, let's say [x(0)=3, y(0)=0], as you have done, the orbit (that's the correct word for what you have called "solution curves") takes off from that point and moves up and to the right until it exits the right edge of the plot.  Something that may not be immediately obvious to you is that if you start at that same initial point and move backward in time, that is, go against the arrows, you will move down and to the right until you exit the right edge of the plot.  Be sure to understand what I am saying here.

Therefore, to sketch the orbit starting at [x(0)=3, y(0)=0], going both backward and forward in time (note the t=−1..1), we do

DEtools:-DEplot({de1,de2}, [x(t),y(t)], t=-1..1, [[x(0)=3, y(0)=0]], x=-4..4, y=-4..4, stepsize=0.01);

Once you see this, you can add a sequence of initial conditions spread along the x axis through the seq() constructions which I showed in my previous answer.  To complete the diagram, you will then add initial conditions along the y axis in the same way.  Then you will adjust the line thickness and color, again as I did in my previous answer.

If you want to ask further questions, include your worksheet to show what you have done.

 

 

 

@dharr I agree with your observation regarding user-friendliness of error messages.  As much as it is tempting to perform type-checking at the level of parsing a proc's args, that can lead to cryptic error messages as you have noted.

Customized parsing of the args adds a few extra lines of code but that can be worth the effort in a proc which is meant to be called at the top-level by other users.  And as you have noted, that "other user" can be you!

That said, Carl in his responses has pointed out that there are cases where type-checking at the argument level can be justified.  That's something to have in mind when writing a general-purpose code that's expected to have a lifetime of more that a few days.

 

 

@Carl Love That's a great overview of Maple's algebra of types, and now I clearly see the logic behind the solution that you had offered to my question.  It would have been difficult for me to deduce all that by reading the help pages.

Among other things, I was surprised to learn that {f}(x) is equivalent to {f(x)}.  Very curious.

I am saving your mini-tutorial for future reference.

 

@Carl Love That's excellent but a bit over my head.  But that's a good thing—it motivates me to dig into the manuals to figure it out.

@acer That's excellent.  Thanks!

@abdulganiy Your PDE is nonlinear, and consequently the system of ODEs produced by the Method of Lines is also nonlinear.  There is no such a thing as a "matrix representation" for a nonlinear system.

 

@Hnx You have 

hz3 := subs(phi2 = hz2, hz1);

That should be

hz3 := subs(phi2(t) = hz2, hz1);

Additionally, the assumptions on your first line don't seem to be essential but they result in clutter in the displayed formulas.  I suggest removing them.  Here is what we get (in the worksheet mode with 1D input, which is my preference):

#Trebuchet, Phase I, 2020-05-27 Ki

restart;

hz1 := g*(r3*m3*cos(phi1(t)) - r1*m2*sin(phi2(t))/sin(phi2(t) - phi1(t)))/theta3;

g*(r3*m3*cos(phi1(t))+r1*m2*sin(phi2(t))/sin(-phi2(t)+phi1(t)))/theta3

hz2 := phi1(t) + arcsin((h + r1*sin(phi1(t)))/r2);

phi1(t)+arcsin((h+r1*sin(phi1(t)))/r2)

hz3 := subs(phi2(t) = hz2, hz1);

g*(r3*m3*cos(phi1(t))+r1*m2*sin(phi1(t)+arcsin((h+r1*sin(phi1(t)))/r2))/sin(-arcsin((h+r1*sin(phi1(t)))/r2)))/theta3

deq := diff(phi1(t), t, t) - hz3 = 0;

diff(diff(phi1(t), t), t)-g*(r3*m3*cos(phi1(t))-r1*m2*sin(phi1(t)+arcsin((h+r1*sin(phi1(t)))/r2))*r2/(h+r1*sin(phi1(t))))/theta3 = 0


 

Download mw.mw

 

@Earl I don't see how that can be done.

@acer Yes, I too have noticed that statement regarding plot3d,viewpoint in the help page, but I haven't run into problems in the few instances where I have combined animation with viewpoint.

As to the number of viewpoint's frames, it appears that it is automatically inherited from the underlying animation.  The animation in this post consists of 181 frames.  In the last line, changing the viewpoint=circleleft to viewpoint=[circleleft, frames=100] (or any other number instead of 100) does not seem to affect the animation.  That's good if it's actually the case.

On an unrelated note, I want to point out that the polyhedron proc in the original post can be shortened quite a bit, like this:

polyhedron := proc(V::list, F::list)
  local f;
  plots:-display([seq(plottools:-polygon(V[f]), f in F)]);
end proc:

 

@Joe Riel You are right.  I don't know why I thought I was getting zero.

@janhardo I suggest that you post your attempts at solving those exercises, otherwise it is not clear what level of help will be beneficial to you.

 

@Kitonum That's very good.  Vote up!  But do you know why that doesn't work when we replace posint with integer?

int(cos(m*t)*cos(n*t), t=-Pi..Pi, allsolutions) assuming integer;

This one returns 0 in Maple 2020.

@Thomas Richard Thanks for looking into this.

@Hamza Brhl 

The wave equation u__tt = u__xx (where subscripts indicate derivatives)

may be converted to a first-order system by introducing the v(x, t)

and w(x, t) variables as follows:

v(x, t) = u__t(x, t)

w(x, t) = u__x(x, t)

Then we note that

v__t = u__tt and u__tt = u__xx and u__xx = w__x

w__t = u__tx and u__tx = v__x

In conclusion, the wave equation u__tt = `u__xx `is equivalent to the system
u__t = v*v__t and v*v__t = w__x*w__t and w__x*w__t = v__x

We apply Maple's pdsolve to that system and get our hands

on the values of v(x, t)(the solution), v(x, t) (the velocity), and

w(x, t)(the slope).

restart;

pde := diff(u(x,t),t) = v(x,t),
       diff(v(x,t),t) = diff(w(x,t),x),
       diff(w(x,t),t) = diff(v(x,t),x);

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

Let's say we wish to solve the wave equation on the interval 0 < x and x < 1 subject

to the initial conditions u(x, 0) = f(x), u__t(x, 0) = g(x).  In terms of the variables

u, v, w this translates to u(x, 0) = f(x), v(x, 0) = g(x), "w(x,0)=f '(x)". For

illustration, I take f(x) = sin*Pi*x, g(x) = 0, and zero for the boundary conditions.

ic := u(x,0) = sin(Pi*x), v(x,0) = 0, w(x,0) = Pi*cos(Pi*x);

u(x, 0) = sin(Pi*x), v(x, 0) = 0, w(x, 0) = Pi*cos(Pi*x)

bc := u(0,t)=0, u(1,t)=0;

u(0, t) = 0, u(1, t) = 0

pdsol := pdsolve({pde}, {bc,ic}, numeric);

_m140058886483840

plots:-display(
  pdsol:-animate(u, t=0..2, color="Red"),
  pdsol:-animate(v, t=0..2, color="Green"),
  pdsol:-animate(w, t=0..2, color="Blue"));

 

Download mw.mw

 

First 37 38 39 40 41 42 43 Last Page 39 of 99