Rouben Rostamian

MaplePrimes Activity

These are replies submitted by Rouben Rostamian

@Rakshak I don't know how to apply the reduction automatically to every expression.  Sorry.

I don't have Maplesim, so I cannot view your worksheet.  Did Maplesim calculate those linearized equations or did you?  In any case, those equations are not correct.  Here are the linearized equations of the double-pendulum as calculated in Maple:

phi and psi: the angles of the pendulum's links with resepct to the vertical
a and b: the lengths of the pendulum's links

m__1 and m__2: the masses

g: gravitational acceleration


de1 := -m[2]*a*b*diff(phi(t),t,t) - g*m[2]*b*psi(t)
                 - m[2]*b^2*diff(psi(t),t,t) = 0;

-m[2]*a*b*(diff(diff(phi(t), t), t))-g*m[2]*b*psi(t)-m[2]*b^2*(diff(diff(psi(t), t), t)) = 0

de2 := -(a*(m[1] + m[2])*diff(phi(t),t,t) + diff(psi(t),t,t)*b*m[2]
                + g*phi(t)*(m[1] + m[2]))*a;

-(a*(m[1]+m[2])*(diff(diff(phi(t), t), t))+(diff(diff(psi(t), t), t))*b*m[2]+g*phi(t)*(m[1]+m[2]))*a

or equivalently:

sys := {
        diff(u(t),t,t) = -g*(m[1] + m[2])*u(t)/(a*m[1]) + g*m[2]*v(t)/(a*m[1]),
        diff(v(t),t,t) = g*(m[1] + m[2])*u(t)/(b*m[1]) - g*(m[1] + m[2])*v(t)/(b*m[1])

{diff(diff(u(t), t), t) = -g*(m[1]+m[2])*u(t)/(a*m[1])+g*v(t)*m[2]/(a*m[1]), diff(diff(v(t), t), t) = g*(m[1]+m[2])*u(t)/(b*m[1])-g*(m[1]+m[2])*v(t)/(b*m[1])}

where I have written u(t) and v(t) for the angles to simplify the typing.

@tomleslie You are bringing up multiple and mutually disjoint issues.

  1. Mathematical modeling issue.  You note that a more accurate model will account for a thermal impedance at a contact area between two bodies of differing temperatures.  Okay, one may take that into account, but then we would be solving a different problem, not the one that the OP has asked.

    Even if we account for the impedance, the model may still be criticized as "physically inaccurate" since it ignores that everything is made of atoms and the model is ignoring those. Where does one stop?
  2. Theoretical considerations.  Within the context of the original model, we have well-established existence and uniqueness theorems of OP's classical initial/boundary value problem.  A unique solution exists even when boundary conditions are discontinuous. In fact, the discontinuities need not be only at the "corners".  A discontinuous initial condition
    u(y,0) = piecewise(y< N/2, 0, 1);   #  0 < y < N

    is perfectly acceptable as well.  In all cases we have an infinitely smooth solution away from the boundaries, but the solution may have discontinuities at the boundaries.  Nothing wrong with that.

  3. Numerical solution. As noted in the previous item, the problem has a well-defined unique solution.  Whether a numerical algorithm can capture that unique solution is a different issue. As you have noted, tweaking the parameters yields various candidates for the solution, each of which is an attempt, by the designer of the numerical algorithm, to capture the problem's one and only true solution. The failure of a numerical algorithm should not be construed as a shortcoming of the mathematical model as you seem to be implying.


@tomleslie The mismatch between the boundary values at (0,0) is not a problem.  In fact, that's a very common occurrence and is harmless as far as solving the PDE is concerned.  All it says is that the solution is discontinuous at that point on the boundary.

Consider the following scenario which leads to such a discontinuity.  Take a metal bar which is initially at temperature u=0 throughout.  At time t=0 bring one end of the bar in contact with the wall of a hot (say 300 degrees) oven.  The temperature in the bar will evolve according to the heat equation u_t = u_xx.  The initial condition is u(x,0)=0.  The boundary condition at x=0 is u(0,t)=300.  The temperature has a discontinuity at (0,0) but the bar doesn't mind.


@acer Very clever idea.  Vote up!

@acer That's a good catch!  I should have thought of that.  Thanks for pointing it out.

@fkohlhepp What you have posted is still incomplete since the data file NACA0012 CDalpha.csv is missing.

But all those details should be unnecessary.  Try entering


If MapleFlow works like Maple, that should print out the equation in all its glory.  Just copy and paste that single equation.  Nothing else is needed.

@fkohlhepp As Axel Vogt has suggested, it would be best if you could just post your eqn2 in a plain text Maple notation.  The details of how that equation is arrived at would be immaterial.  Then people can copy and paste that equation into Maple and suggest good ways to solve it.

In the absence of that equation, I have made up something which is remotely similar.  Maple's fsolve() indeed has difficulty solving the equation—I don't know why—but the NewtIt() proc which I posted earlier solves the equation right away.


NewtIt := proc(f, x0, {eps:=1e-5, nmax:=10})
        local x, xnew;
        x := evalf(x0);
        xnew := x - f(x)/D(f)(x);
        if is(abs(xnew - x) < eps) then
                return xnew;
        elif nmax = 0 then
                error "number of iterations exceeded nmax";
                return `procname`(f, xnew, :-eps=eps, :-nmax=nmax-1);
        end if;
end proc:

f := u -> int(cos(1+arctan(u,r)), r=1..2, numeric) - u^2;

proc (u) options operator, arrow; int(cos(1+arctan(u, r)), r = 1 .. 2, numeric)-u^2 end proc

plot(f(u), u=0..1);

fsolve() fails; I don't know why:

fsolve(f(u)=0, u=0..1);

Error, (in fsolve) cannot determine if this expression is true or false: 1.*abs(Im(r))+1.*abs(Re(r)) <= 0.

but NewtIt() works fine:

NewtIt(f, 0.2);




@fkohlhepp Regarding "its already essentially here in another post": That's not good enough. Your problem lies in the definition of the function Thust which you have not shown. The more information you provide, the better help you will get, and perhaps it may turn out that your "home-built" iteration scheme is not needed at all.


@Mohamed Mohsin If you don't know how to solve a problem, try to find an easier problem that you can solve.  If you don't know how to solve that easier problem, find even an easier one, and so on, until you arrive at a point where you are able to do something.

For instance, in your case, do you know how to solve

D^γ S(t) = S(t)

and if not, do you know how to solve

D S(t) = S(t)

and so on.

You should also ask yourself:

  • Do I know anything about fractional differential equations?
  • Do I know anything about the homotopy perturbation method?  

Answers to these questions have nothing to do with Maple. Begin with investigating these mathematical prerequisites of your master's thesis.  Coding and calculating with Maple will come only after you have positively answered those questions.


I assume that you have done more for your master's thesis than what you have shown here.  Explain what you have so far, and where it is that you need help with Maple.

Otherwise it looks like you are asking for someone else to do your master's thesis for you.

@ijuptilk I am puzzled by your statement that you have an analytic solution for the problem since your initial condition v(z,0)=0 is inconsistent with the rest of the system.  Examine what you think is an analytical solution.  It's either wrong or it violates v(z,0)=0.

The reason that v(z,0)=0 is inconsistent is that there is no t derivative of v(z,t) in your system of PDEs, so you cannot specify an initial condition on v at t=0.

If you remove the v(z,0)=0 condition, you can obtain a symbolic solution to the problem without much difficulty.

@ijuptilk As it has been pointed out to you, sin(pi/2) is 1, so saying theta(z,0) = thetab*sin(pi/2)*z is the same as saying theta(z,0) = thetab*z.  Do you see that?

I suspect that what you really want to say is theta(z,0) = thetab*sin(pi/2*z).  Look closely at what you write. Mathematics is a very precise language.

Furthermore, are you sure you want sin(pi/2*z)?  Changing it to sin(pi/d*z) makes it compatible with the boundary conditions.

Finally, since you know the symbolic solution, why are you bothering with a numeric one?  Can you give the reason?

@ecterrab Thank you for your customary lightningly quick reaction to issues brought up here in this forum.  I will either install the Physics Update or just do `-`(1, p.q) for now.  Both methods will enable me to move forward.

Thank you again!


@ecterrab Your explanation regarding type versus Identiry makes good sense.  With that hint I have been able to move over a good deal of my (rather extensive) old code from the LinearAlgebra vectors to Physics vectors.

I have run into a puzzle which I would appreciate if you could illuminate for me.  Consider the following code:




my_module := module()
        uses Physics:-Vectors;
        export my_proc;
        my_proc := proc(p::PhysicsVectors, q::PhysicsVectors)
                1 - p.q;
        end proc:

end module:



Error, (in my_proc) invalid input: Physics:-Vectors:-- uses a 2nd argument, b, which is missing



I don't see why Maple complains.  I can get around that by replacing the line 1 - p.q with `-`(1, p.q) but that shouldn't be necessary, or should it?



3 4 5 6 7 8 9 Last Page 5 of 84