ecterrab

13431 Reputation

24 Badges

19 years, 358 days

MaplePrimes Activity


These are replies submitted by ecterrab

@nm

I already uploaded v.264, then installed it in my copy of Maple, it worked fine. Version 264 has several refinements that make the code find a solution for some more difficult examples. These are two of them:
 

Articolo's textbook example 7.9.1

pde[1] := diff(u(r, theta, t), t, t) = (diff(u(r, theta, t), r)+r*(diff(u(r, theta, t), r, r))+(diff(u(r, theta, t), theta, theta))/r)/(4*r)

iv[1] := u(1, theta, t) = 0, u(r, 0, t) = 0, u(r, (1/2)*Pi, t) = 0, u(r, theta, 0) = (-r^4+r^2)*sin(2*theta), (D[3](u))(r, theta, 0) = 0

pdsolve([pde[1], iv[1]])

u(r, theta, t) = `casesplit/ans`(Sum(-(1/6)*sin(2*theta)*(BesselJ(0, lambda[n]*r)*lambda[n]*r-2*BesselJ(1, lambda[n]*r))*cos((1/2)*lambda[n]*t)*hypergeom([], [5], -(1/4)*lambda[n]^2)/((hypergeom([2, 5/2], [3, 3, 4], -lambda[n]^2)*lambda[n]-hypergeom([3/2, 2, 2], [1, 3, 3, 3], -lambda[n]^2)*lambda[n])*r), n = 1 .. infinity), {And(lambda[n] = BesselJZeros(2, n), 0 <= lambda[n])})

(1)


Articolo's Exercise 7.15, with 6 BC/IV, two for each variable

pde[2] := diff(u(x, y, t), t, t) = 1/4*(diff(u(x, y, t), x, x)+diff(u(x, y, t), y, y))-(1/10)*(diff(u(x, y, t), t))

iv[2] := (D[1](u))(0, y, t) = 0, (D[1](u))(1, y, t)+u(1, y, t) = 0, (D[2](u))(x, 0, t)-u(x, 0, t) = 0, (D[2](u))(x, 1, t) = 0, u(x, y, 0) = (1-(1/3)*x^2)*(1-(1/3)*(y-1)^2), (D[3](u))(x, y, 0) = 0

 

This problem is tricky ... There are three independent variables, therefore two eigenvalues (constants that appear separating variables by product) in the Sturm-Liouville problem. But after solving the separated system and also for the eigenvalues, the second eigenvalue is equal to the first one, and in addition cannot be expressed in terms of known functions, because the equation it solves cannot be inverted.

 

The solution to the problem is now computable in v.264

 

pdsolve([pde[2], iv[2]])

u(x, y, t) = `casesplit/ans`(Sum((1/6)*(exp((1/10)*t*(-200*lambda[n]^2+1)^(1/2))*(-200*lambda[n]^2+1)^(1/2)+exp((1/10)*t*(-200*lambda[n]^2+1)^(1/2))+(-200*lambda[n]^2+1)^(1/2)-1)*(cos(lambda[n]*y)*lambda[n]+sin(lambda[n]*y))*exp(-(1/20)*t*((-200*lambda[n]^2+1)^(1/2)+1))*(6*lambda[n]^2*cos(lambda[n])^2+cos(lambda[n])^2-5*lambda[n]^2-1)*cos(lambda[n]*x)*(lambda[n]^2*sin(lambda[n])-lambda[n]*cos(lambda[n])+sin(lambda[n]))/((-200*lambda[n]^2+1)^(1/2)*(-cos(lambda[n])^4+(sin(lambda[n])*lambda[n]-1)*cos(lambda[n])^3+(lambda[n]^2+sin(lambda[n])*lambda[n]+1)*cos(lambda[n])^2+(lambda[n]^2+2*sin(lambda[n])*lambda[n]+1)*cos(lambda[n])+lambda[n]*(lambda[n]+sin(lambda[n])))*(cos(lambda[n])-1)*lambda[n]^4), n = 0 .. infinity), {And(tan(lambda[n])*lambda[n]-1 = 0, 0 <= lambda[n])})

(2)

``


 

Download Update_v.264.mw

Edgardo S. Cheb-Terrab
Physics, Differential Equations and Mathematical Functions.

@Preben Alsholm 

You missed solve's option 'allsolutions' (one could argue it should be used by default, that is another story). With it, you have:

@Preben Alsholm 

I saw this comment by you only today. The distrubuted Maplesoft Physics Updates got updated with the _libraryversion number of Maple 2018.2.1 a week or so ago. Since then, the Warning message you saw is not displayed.

In regards to your comment, there is the "why" and the "what message".

The "why": as you know the Maplesoft Physics Updates pack not just updates to Physics but also to the Differential Equations and Mathematical Functions code, and bug fixes not in those three areas but necessary for them to work properly, some of these mentioned in previous Mapleprimes posts.

Now, Maple is a very interconnected set of software routines. If you use an adjustment in the code written for Maple 2018.2 with a different library, e.g. 2018.1, or the other way around, I cannot assure you that all the things will work as intended. That is the "why" of the message, that includes a sentence (not shown in your comment) that indicates from where you can download an update for the version of Maple you have.

I imagine you know the "why" above but didn't find the message appropriate anyway, So, what message would you find more appropriate for this situation?

Edgardo S. Cheb-Terrab
Physics, Differential Equations and Mathematical Functions, Maplesoft

@John Fredsted 

The current Maplesoft Physics Updates will not work well in previous Maple releases. Many of its features depend on things that only exist in Maple 2018, and several fixes only work as expected with the library of Maple 2018. I'd suggest you to install the latest Maplesoft Physics Updates for Maple 2017, but it does not include the change I did yesterday regarding this issue you mentioned.

Edgardo S. Cheb-Terrab
Physics, Differential Equations and Mathematical Functions, Maplesoft

@lucenalex 

I'm a bit overloaded today and tomorrow but will try to find time to show how to formulate the problem - hopefully before a couple of days from today.

Best

Edgardo S. Cheb-Terrab
Physics, Differential Equations and Mathematical Functions, Maplesoft

@Sofi 

The Maplesoft Physics Updates with code for solving the PDE & BC you posted only works for Maple 2018. Regarding the math (your other question), give a look at  this other Mapleprimes post where it is all explained.

Edgardo S. Cheb-Terrab
Physics, Differential Equations and Mathematical Functions, Maplesoft

@Mariusz Iwaniuk 

Note that the solution you posted, using new routines described in this other Mapleprimes post, actually shows two equations for lambda[n], one of which is sort of redundant (due to an intermediate call to simplify that split the single equation into two). This is adjusted in version 200 or higher of the Maplesoft Physics Updates, so that now the output has one single equation defining lambda[n].

Edgardo S. Cheb-Terrab
Physics, Differential Equations and Mathematical Functions, Maplesoft

@Rouben Rostamian

 

Just trying your worksheet today (version 195 of the Updates). I see that the other two issues you mentioned are also resolved.

pde := diff(u(x, y), x, x)+diff(u(x, y), y, y) = 0

bc1 := u(0, y) = 0, u(Pi, y) = 0, u(x, 0) = 0, u(x, Pi) = 0

 

Check the method that resolves the problem

infolevel[pdsolve] := 2

2

(1)

pdsolve({bc1, pde})

* trying method "_Fn" for 2nd order PDEs
   -> trying "linear_in_xt"
   -> trying "BC_equal_0"
* trying method "_Cn_cn" for 2nd order PDEs
   -> trying _Cn_cn
* trying method "Wave" for 2nd order PDEs
   -> trying "Cauchy"
   -> trying "SemiInfiniteDomain"
   -> trying "WithSourceTerm"
* trying method "Heat" for 2nd order PDEs
   -> trying "SemiInfiniteDomain"
   -> trying "WithSourceTerm"
* trying method "Series" for 2nd order PDEs
   -> trying "ThreeBCsincos"
   -> trying "FourBC"

   -> trying "ThreeBC"
   -> trying "ThreeBCPeriodic"
   -> trying "WithSourceTerm"
   -> trying "ThreeVars"
* trying method "Laplace" for 2nd order PDEs
   -> trying a Laplace transformation
* trying method "Fourier" for 2nd order PDEs
   -> trying a fourier transformation

   <- fourier transformation successful
<- method "Fourier" for 2nd order PDEs successful

 

u(x, y) = 0

(2)

By the way note (not documente yet) that since some versions of the Updates you can specify the method

pdsolve({bc1, pde}, method = Fourier)

* trying method "Fourier" for 2nd order PDEs
   -> trying a fourier transformation
   <- fourier transformation successful
<- method "Fourier" for 2nd order PDEs successful

 

u(x, y) = 0

(3)

On the methods:

indices(`pdsolve/BC/methods`)

[1], [2], [2, "Series"], [3], ["high_order"], [2, "_Fn"], [2, "Wave"], ["system"], [2, "Heat"]

(4)

These typically have submethods that in turn may more of that. For example

`pdsolve/BC/methods`[2]

["_Fn", "_Cn_cn", "Wave", "Heat", "Series", "Laplace", "Fourier", "Generic", "PolynomialSolutions", "Linear_diff_op"]

(5)

If the optional argument (method or methods) is a list, then only those methods and in the order you specify will be tried.

 

Some of the names you see today in these lists of methods will change to something more understandable.

 

The other problem you mentioned

bc2 := u(0, y) = 0, u(Pi, y) = sinh(Pi)*cos(y), u(x, 0) = 0, u(x, Pi) = 0

pdsolve({bc2, pde})

* trying method "_Fn" for 2nd order PDEs
   -> trying "linear_in_xt"
   -> trying "BC_equal_0"
* trying method "_Cn_cn" for 2nd order PDEs
* trying method "Wave" for 2nd order PDEs
   -> trying "Cauchy"
   -> trying "SemiInfiniteDomain"
   -> trying "WithSourceTerm"
* trying method "Heat" for 2nd order PDEs
   -> trying "SemiInfiniteDomain"
   -> trying "WithSourceTerm"
* trying method "Series" for 2nd order PDEs
   -> trying "ThreeBCsincos"
   -> trying "FourBC"

   <- trying "FourBC" successful
<- method "Series" for 2nd order PDEs successful

 

u(x, y) = Sum(2*sin(n*y)*exp(Pi*n)*n*sinh(Pi)*((-1)^n+1)*(exp(n*x)-exp(-n*x))/(Pi*(exp(2*Pi*n)-1)*(n^2-1)), n = 2 .. infinity)

(6)

The sum starts at n = 2.

NULL


Edgardo S. Cheb-Terrab
Physics, Differential Equations and Mathematical Functions, Maplesoft

Download The_other_two_issues.mw

 

 

@Mariusz Iwaniuk 

Thanks, version 194 includes also a fix to factor, but the advanced version of that fix got by accident into the code for 2018.2. It is corrected now in version 195 of the Maplesoft Physics Updates.

Edgardo S. Cheb-Terrab
Physics, Differential Equations and Mathematical Functions, Maplesoft

@Rouben Rostamian  

First of all, nice reply yours. It is exciting when the level of the discussion goes higher and at the same time stays respectful of other's (sometimes different) opinions.

On the contents: you are right in the distinction about classical and weak solutions; mainly: the solution could be discontinuous at the boundary conditions (e.g. in the example posted). So yes, I goofed here, prematurely discarding the problem because of the inconsistency of the boundary conditions posted.

On the other issues that resulted from this post: Maple's pdsolve does compute weak solutions, examples of that were posted, but the one of this post is not solved when it could - for instance, you show how - and the method (skipping details) is systematic. To cover this example and a family of related cases more work is necessary. We are working on that in this moment, among other things. I will post here again when it is ready.

@Carl Love 

Physics:-Fundiff implements functional differentiation, it computes a functional derivative, an operation at the realm of variational calculus. Because Physics is all about a principle of minimal action - a variational principle - it so happens that basically all formulations in classical and quantum physics require computing functional derivatives (of the action, and equating the result to 0). If one day you happen to work with variational calculus you will be surprised with the computational power behind Fundiff . Some of it is shown in several posts here in Mapleprimes about quantum mechanics and classical field theories.

If you are still not finding the computation familiar, check Wikipedia at https://en.wikipedia.org/wiki/Functional_derivative - to go to the point, search the page for Functional derivative o a function.

On the other hand, Physics:-diff implements, in a general/formal sense what in jet calculus is called jet space - frequently mentioned in the context of differential equations. Generally speaking, you consider independent variables, dependent variables, and their partial derivatives as a set of equivalent independent objects with respect to which you can differentiate. For example, all the Lie theory of symmetries of differential equations is based on this concept (to see the simplistic core of the implementation, give a look at ?PDEtools:-ToJet). The implementation of Physics:-diff in computer algebra is tricky due to the way Maple's diff (without Physics) got implemented day one, making necessary to implement together the concept of "computing the derivative with respect to x at fixed f(x, ...)" (see the help page), and in this way implement, for the first time in Maple, true partial derivatives, for real, as the ones you use when formulating thermodynamics for instance.

People sometimes think of these two commands as interrelated because, if you are functionally differentiating an integral, you can frequently express the result using a sum of terms all of which involve calls to Physics:-diff. But these two commands are actually different. Most of the things mentioned above are described in more details in the corresponding help pages.

By the way, it is remarkable that while these commands are at the very core of formulations in Physics, I actually coded them using Maple in 1994, they where my first programs, and it was the first time these operations got implemented in computer algebra - I even published a paper with them in Computer Physics Communications, still today there is no equivalent for them in Mathematica.

Edgardo S. Cheb-Terrab
Physics, Differential Equations and Mathematical Functions, Maplesoft.

 

There is clearly an issue here (BTW, good catch vv!), as usual splitting into two: "what is happening?" and "is that ok?". What is happening: Change is a wrapper for PDEtools:-dchange; debugging it you see that its input for 2J is not nearly the same as for J - mainly in the case of 2J only, the change of variables is performed first in the inner integral. The undesired output for J (not for 2J) appears while changing variables in J taking the iterated double integral as a multiple integral. "Is that OK?": generally speaking yes, however in cases like this changing variables one integral at a time takes you to a much simpler and readable result.

I'll give this a closer look later today or tomorrow and adjust the code to return the simpler form - then upload the adjustment within the next "Maplesoft Physics Updates for 2018.1" - will write here when it's ready.

Edgardo S. Cheb-Terrab
Physics, Differential Equations and Mathematical Functions, Maplesoft

Carl, very beautiful post, illustrative and insightful!

Best

Edgardo S. Cheb-Terrab
Physics, Differential Equations and Mathematical Functions, Maplesoft

 

 

@zhuxian 

Check the help page of PDEtools:-Library (the set of programming routines with which the PDEtools package is built), in particular, PDEtools:-Library:-Specialize_Cn and its optional arguments are what you want.

Edgardo S. Cheb-Terrab
Physics, Differential Equations and Mathematical Functions, Maplesoft

@digerdiga 
 

restart

"  with(Physics):"

 

Set your metric, and in doing so give a closer look to both the default signature and the coordinates for this metric:

g_[[15, 23, 2]]

" _______________________________________________________"

 

" `Systems of spacetime Coordinates are:` {X=(u,r,theta,phi)}"

 

" `Default differentiation variables for d_, D_ and dAlembertian are:` {X=(u,r,theta,phi)}"

 

"`The ``Eddington (1924)`` metric in coordinates ``Finkelstein(1958)`"

 

"`Parameters: `[m]"

 

"`Resetting the signature of spacetime from `(-` - - +`)` to `(-` + + +`)` in order to match the signature in the database of metrics`"

 

" _______________________________________________________"

 

Physics:-g_[mu, nu] = Matrix(%id = 18446744078371044758)

(1)

So: the signature and coordinates are, respectively, (- + + +) and (u, r, theta, phi).

 

Now: what I explained in the previous reply is not "how to change the metric" but "how to set the signature". Recalling, generally speaking, for a curved spacetime, the signature only tells the form of the Minkoski metric in an orthonormal frame (the tetrad system, that you can define locally, at any point of spacetime), and generally speaking it changes nothing in the form of the non-flat metric you choose ([15,23,2] in this example).

 

So you can change the signature and that will not change the metric. For example, this is the signature in this moment

Setup(signature)

[signature = `- + + +`]

(2)

Here I change it and the metric is still the same as (1)

Setup(signature = `+---`)

[signature = `+ - - -`]

(3)

"g_[]"

Physics:-g_[mu, nu] = Matrix(%id = 18446744078370950190)

(4)

You noticed this but interpreted it incorrectly: the above is correct and according to the design here, there is no reason to think that changing the signature should change the metric (as said it only changes the metric in the local, tetrad system of references).

 

On the other hand, it is the case that, computationally, we frequently want to change the signature and with that also redefine the metric, signs, lines and columns, even reorder the list of coordinates. The command for that is Redefine . As said in the previous reply: the value of computer algebra is proportional to your use of the help pages, and getting skilled in grabbing the information you need without reading the whole thing is part of the game...

 

Anyway, here is some of this valuable command. First reset the signature to what it was when you loaded the metric

"Setup(signature = `-+++`)"

[signature = `- + + +`]

(5)

The metric is the same as (1) and (4)

"g_[]"

Physics:-g_[mu, nu] = Matrix(%id = 18446744078370862494)

(6)

Now open the help page for Redefine and play with it a bit. Your desired signature does not change the order of the coordinates so the first thing I would explore is this:

Redefine(metric, tosignature = `+---`)

Matrix(%id = 18446744078370860798)

(7)

But the metric has not been set to (7), the above is only a preview of the redefinition you indicated, the metric it is still equal to (1):

"g_[]"

Physics:-g_[mu, nu] = Matrix(%id = 18446744078370862494)

(8)

As explained in the help page this permits you verifying that the redefinition you are exploring is indeed what you want.  Just to get confidence in the command, ask it to redefine everything, expecting that the coordinates will not change of course

"Redefine(all,tosignature = `+---`)"

[X], Matrix(%id = 18446744078483705358)

(9)

Indeed the coordinates are still X = (u, r, theta, phi). Assuming that the changed metric (7) , which is the same one shown in (9), is the one you want, then proceed according to what the help page of Redefine says and set this metric

"Redefine(setmetric,tosignature = `+---`)"

Matrix(%id = 18446744078483694646)

(10)

Check now the metric: it got set to what you wanted, it has an overall different sign with respect to (1)

"g_[]"

Physics:-g_[mu, nu] = Matrix(%id = 18446744078371043558)

(11)

Summary: the command for redefining the metric, changing signs and eventually lines and columns, by changing the signature, is Redefine.

 

On the other question, how to collect in the line element, note first that it already comes collected:

g_[line_element]

(r-2*m)*Physics:-d_(u)^2/r+2*Physics:-d_(u)*Physics:-d_(r)-r^2*Physics:-d_(theta)^2-r^2*sin(theta)^2*Physics:-d_(phi)^2

(12)

But suppose it were not collected, say as in

"expand((r-2*m)/r*Physics:-d_(u)^2+2*Physics:-d_(u)*Physics:-d_(r)-r^2*Physics:-d_(theta)^2-r^2*sin(theta)^2*Physics:-d_(phi)^2)"

Physics:-d_(u)^2-2*Physics:-d_(u)^2*m/r+2*Physics:-d_(u)*Physics:-d_(r)-r^2*Physics:-d_(theta)^2-r^2*sin(theta)^2*Physics:-d_(phi)^2

(13)

In the above, d(u)^2 is not collected. First you need to see "what" is to be collected, the simplest way is to use lprint

"lprint(Physics:-d_(u)^2-2/r*Physics:-d_(u)^2*m+2*Physics:-d_(u)*Physics:-d_(r)-r^2*Physics:-d_(theta)^2-r^2*sin(theta)^2*Physics:-d_(phi)^2)"

Physics:-d_(u)^2-2*Physics:-d_(u)^2*m/r+2*Physics:-d_(u)*Physics:-d_(r)-r^2*Physics:-d_(theta)^2-r^2*sin(theta)^2*Physics:-d_(phi)^2

 

 

You see it is the d_ command you need to collect. So this collects the way you wanted.

"collect(Physics:-d_(u)^2-2/r*Physics:-d_(u)^2*m+2*Physics:-d_(u)*Physics:-d_(r)-r^2*Physics:-d_(theta)^2-r^2*sin(theta)^2*Physics:-d_(phi)^2,d_)"

-r^2*sin(theta)^2*Physics:-d_(phi)^2+2*Physics:-d_(u)*Physics:-d_(r)-r^2*Physics:-d_(theta)^2+(1-2*m/r)*Physics:-d_(u)^2

(14)

However, as said, never use a command without giving a quick look at is help page. Check collect  and you will see the options you have.

``


 

Download Redefining_the_metric_by_changing_the_signature.mw

Edgardo S. Cheb-Terrab
Physics, Differential Equations and Mathematical Functions, Maplesoft

First 24 25 26 27 28 29 30 Last Page 26 of 60