Rouben Rostamian

MaplePrimes Activity


These are replies submitted by Rouben Rostamian

@arshl Solving your modified problem leads to an integral of the type
which is not expressible in terms of elementary functions. This means that your problem has no symbolic solution.  If this problem arises in a real application, then you may seek a numerical solution.  But for that you need to supply the numerical values of all the parameters.

 

 

@Carl Love Your approach requires the initial condition to be twice differentiable.  We can see that by replacing the initial condition u(x,0)=x*x(x-1) by u(x,0)=g(x) and then observing the resulting expression for bounds.

Thus, for instance, this approach fails with the perfectly fine initial condition u(x,0)=Heaviside(x-1/2).

If, however, the twice differentiability condition on the initial condition is met, then your approach will produce the correct solution.

That said, let me reiterate that the original statement of the PDE is flawed.  The sign of the second derivative term is wrong.  To see the problem, in your solution change
psol:= pdsolve(pde, bounds, numeric):
to
psol:= pdsolve(pde, bounds, numeric, spacestep=s);
and watch what happens to the amplitude of the solution with s=0.01 and s=0.001.

Question: The use of 'union=' in your solution is new to me. Where is it documented?

Added later:
I see that C-style assignments such as
x := 12;
x += 1;
x *= 5;

also work.  These have been on my wish-list for a long time!

@abbastalebi Here is how:

plots:-odeplot(dsol, [eta, diff(f(eta),eta)]);      # first derivative
plots:-odeplot(dsol, [eta, diff(f(eta),eta,eta)]);  # second derivative

@Mariusz Iwaniuk That's excellent. A thumb up!  I wonder whether it's possible to recover the original integral through differentiating the answer.

@tomleslie No, I didn't mean to.  I played around with the equation a bit and, among other things, changed the sign of the exponent, but then I forgot to restore it.  Mea culpa!

That said, the function exp(t*x) grows astronomically large very quickly.  I can't think of any realistic problem leading to such a term.  I like the problem better with the minus in the exponent.

@kambiz1199 A complete answer to your question depends on the details of your worksheet that you have not shown, so people are just guessing what you may have in there.

Here is a complete worksheet where for simplicity I have taken the N matrix to be 2x4, and have taken some random expressions for N[1] and N[2].  If this works for you, then you should be able to extend it as needed.

restart;

interface(rtablesize=20):

N[1] :=  x^2 + y^2;
N[2] := x^3*y + 7*y^3;

x^2+y^2

x^3*y+7*y^3

Nmat := < seq('N[i],0', i=1..2);  
          seq('0,N[i]', i=1..2) >;

Matrix(2, 4, {(1, 1) = x^2+y^2, (1, 2) = 0, (1, 3) = x^3*y+7*y^3, (1, 4) = 0, (2, 1) = 0, (2, 2) = x^2+y^2, (2, 3) = 0, (2, 4) = x^3*y+7*y^3})

Z := < D[1], 0; 0, D[2]; D[1], D[2] >;

Matrix(3, 2, {(1, 1) = D[1], (1, 2) = 0, (2, 1) = 0, (2, 2) = D[2], (3, 1) = D[1], (3, 2) = D[2]})

Matrix(3,4, (i,j)-> add(Z[i,k](unapply(Nmat[k,j], [x,y]))(x,y),k=1..2));

Matrix(3, 4, {(1, 1) = 2*x, (1, 2) = 0, (1, 3) = 3*x^2*y, (1, 4) = 0, (2, 1) = 0, (2, 2) = 2*y, (2, 3) = 0, (2, 4) = x^3+21*y^2, (3, 1) = 2*x, (3, 2) = 2*y, (3, 3) = 3*x^2*y, (3, 4) = x^3+21*y^2})


 

Download mw2.mw

@arshl  We apply dsolve to find the solution v[3](t) of the differential equation.  The solution depends on many parameters such as p1, R, k11, a11, ro, cp2, etc.  There is nothing special about p1.  Maple is quite happy solving differential equations involving parameters.  Here is a simple example:

de := diff(x(t),t,t) + omega^2 * x(t) = 0;
dsolve(de, x(t));

The solution x is a function of t, so we write it as x(t).  It also depends on the parameter omega, but we don't write x(t,omega).

 

@tomleslie Thanks for pointing this out.  I should have noticed that the integral symbol returned by dsolve is colored gray, indicating an inert Int.

 

@arshl In my previous reply I said that all "pi" should be "Pi".  You haven't made that change yet.  You should.

Looking further down your code, I don't understand the meaning of 
m1(p1, my(t)) := (1/2)*b11*ro*pi*R^2*(my(t)*cp2/(R*p1)-cp3)/p1;
What are your attempting to do there?

Perhaps it would be easier if you explain just the mathematics of your problem (no Maple).  Then someone may suggest ways to convert that mathematics into Maple.

Your say you are applying dsolve, but you haven't shown the code that produces the result that you have shown.  Showing only a part of your worksheet is not very informative.  Post your entire worksheet as an attachment.  How to do that? In the window in which you edit your message before posting, note the big green arrow near the top edge.  Click on that arrow to attach your worksheet.

But...

Before you do that, change all occurrences of "pi" in your worksheet to "Pi".  They are not the same thing.

 

@escorpsy

restart;

We have seen that the integral

exp(-h*Pi*v*(x+y)/beta) / cosh(Pi*v*(y-x)/beta)^(4+h):
Int(%, y=-infinity..z, x=-infinity..w);

Int(exp(-h*Pi*v*(x+y)/beta)/cosh(Pi*v*(y-x)/beta)^(4+h), y = -infinity .. z, x = -infinity .. w)

may be evaluated in terms of elementary functions when h is a given
negative integer.  In your latest worksheet you ask whether it is possible
to obtain the answer in terms of an unspecified h.

 

That is doubtful.  Consider the much simpler case

int(1/cosh(x)^h, x);

int(1/cosh(x)^h, x)

showing that Maple is unable to obtain an explicit formula for general h.

 

Specifying h as a non-integer does not help either:

int(1/cosh(x)^(7/3), x);

int(1/cosh(x)^(7/3), x)

But when h is an integer, we get an explicit formula:

int(1/cosh(x)^h, x) assuming h::posint;

(Sum((Product(1-1/(-h+2*j+1), j = 0 .. i))*sech(x)^(h-2*i-2)/h, i = 0 .. -ceil(-(1/2)*h)-1))*tanh(x)+2*(Product(1+1/(-h+2*j-1), j = 1 .. -ceil(-(1/2)*h)))*arctan(exp(x))

Unfortunately the result is too complex for evaluating the integral on -infinity, 0.

int(1/cosh(x)^h, x=-infinity..0) assuming h::posint;

limit(-(Sum((Product((h-2*j)/(h-2*j-1), j = 0 .. i))*sech(x)^(h-2*i-2)/h, i = 0 .. -ceil(-(1/2)*h)-1))*tanh(x)-2*(Product((h-2*j)/(h-2*j+1), j = 1 .. -ceil(-(1/2)*h)))*arctan(exp(x))+(1/2)*(Product((h-2*j)/(h-2*j+1), j = 1 .. -ceil(-(1/2)*h)))*Pi, x = -infinity)

Your actual integral is more complex than the simple illustration above,

therefore I wouldn't hold out hope for an explicit formula for its value.

 

Download mw.mw

@mmcdara I pointed out the ambiguity of the definitions of c_1 and c_2, but farazhedayati said that he has reasons for treating xbar as an independent variable, so I just implemented what he is asking.  It's up to him to make sense of the result.

 

@farazhedayati To add the c_2 values, we do

add(c_2(my_x, my_xbar, i), i=1..n);

and we get 1.479074043.

The definition
    c1:=diff(b, xbar);
is ambiguous.  To illustrate that in a very simple case, consider the function f = m + 2*x + 2*y, where m is the mean of x and y.  You ask, what is the derivative of f with respect to m?

Answer #1:  Clearly df/dm = 1.

Answer #2:  By the definition of m we have m = (x+y)/2, therefore x+y=2*m, that is, 2*x+2*y=4*m.  Consequently, f=m+4*m=5*m, and therefore df/dm=5.

Which answer would you accept?

 

This alternative version provides a better demo.

Download stretch-octahedron2.mw

First 48 49 50 51 52 53 54 Last Page 50 of 99