ecterrab

14540 Reputation

24 Badges

20 years, 22 days

MaplePrimes Activity


These are answers submitted by ecterrab

Hi Christopher,

Note first that the Maplesoft R&D Physics webpage has the Physics Updates for the previous releases Maple 2017, Maple 2016 and Maple 2015 - only the latest one for each of them. This webpage was cleaned-up removing updates for some older releases.

Not exactly what you asked but, regarding maple 2018, the first release that has updates distributed through the Maple cloud, if you 'view' the package you see the instructions to install any previous version of the Physics Updates for Maple 2018. Being that we are at version 192, that means you have an archive of 191 previous versions. This repository may remain available as such when Maple 2019 becomes the current release next year.

On the versions for Maple 2018, additionally, not mentioned in that 'view' of the package's update, a few updates ago I introduced the following: Physics:-Version, that tells you the version on the cloud and in your computer so that you can see whether you have the latest version, now also installs any version. For example, if you have the latest update, input Physics:-Version(123) and you will have version 123 installed (and remember that in version 123 Physics:-Version does not have this capability of installing the package).

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

Hi,

First of all, assuming you are using Maple 2018, be sure you have installed the free upgrade to Maple 2018.2 and also have the latest Maplesoft Physics Updates (version 192 or higher) installed in order to reproduce what is shown below.

 

with(Physics)

 

The following is what you entered first, and indicates that A and C act on the same space, then B and C act on the same space, but A and B do not act on the same space - hence A and B commute, but not A and C or B and C. All bases are discrete with dimensions Na, Nb and Na Nb

 

Setup(hilbertspaces = {{A, C}, {B, C}}, quantumbasisdimension = {A = 1 .. Na, B = 1 .. Nb, C = 1 .. Na*Nb}, quantumdiscretebasis = {A, B, C})

[disjointedspaces = {{A, C}, {B, C}}, quantumbasisdimension = {A = 1 .. Na, B = 1 .. Nb, C = 1 .. Na*Nb}, quantumdiscretebasis = {A, B, C}]

 

Regarding your first question: indicating the dimensions of A and B does not determine the dimension of C, so indeed you need to tell the system what is C's dimension, if you are going to use that (it mainly serves the purpose of computing projectors and their use), not otherwise.

 

About your intention of creating a bracket rule such that `⊗`(`⊗`(Bra(A, i), Bra(B, j)), Ket(C, t)) = X[i, j](t) you cannot define a bracket rule of this sort (two Bras times one Ket, that is not really a Bracket)  but you can define single Kets that belong to C and that have two quantum numbers, ie representing a tensor product, and you can assign dimensions to each of the quantum numbers as you wish. Quoting from the help page for Setup,

 

"It is valid to state different dimensions for different quantum numbers of Kets of a single basis. For that purpose, if the label of the basis is A, you can set the dimension of each quantum number entering quantumbasisdimension = {A[1] = N, A[2] = M}, where N and M are positive integers, algebraic expressions representing them, or ranges of the form a .. b. This dimension of a quantum basis is used to construct projectors when you call the Projector  command."

 

For example, add another quantum space label, say F, and set the dimensions accordingly, noting that your dependence on t is a continuous dependence (in your post I see it as a function, not as an index)

 

Setup(redo, hilbertspaces = {{A, C, F}, {B, C, F}}, quantumbasisdimension = {A = 1 .. Na, B = 1 .. Nb, C[1] = Na, C[2] = Nb}, quantumcontinuousbasis = F)

[disjointedspaces = {{A, C, F}, {B, C, F}}, quantumbasisdimension = {A = 1 .. Na, B = 1 .. Nb, C[1] = Na, C[2] = Nb}, quantumcontinuousbasis = {F}]

So your bases look like this

`𝕀__A` := Projector(Ket(A, i)); `𝕀__B` := Projector(Ket(B, i)); `𝕀__C` := Projector(Ket(C, i, j)); `𝕀__F` := Projector(Ket(F, t))

Sum(Physics:-`*`(Physics:-Ket(A, i), Physics:-Bra(A, i)), i = 1 .. Na)

Sum(Physics:-`*`(Physics:-Ket(B, i), Physics:-Bra(B, i)), i = 1 .. Nb)

Sum(Sum(Physics:-`*`(Physics:-Ket(C, i, j), Physics:-Bra(C, i, j)), i = 0 .. Na-1), j = 0 .. Nb-1)

Int(Physics:-`*`(Physics:-Ket(F, t), Physics:-Bra(F, t)), t = -infinity .. infinity)

For instance,

Bra(A, n).`𝕀__A`

Physics:-Bra(A, n)

Bra(F, tau).`𝕀__F`

Physics:-Bra(F, tau)

From you post I understand that your bracket rule would then be

%Bracket(Bra(C, i, j), Ket(F, t)) = X[i, j](t)

%Bracket(Physics:-Bra(C, i, j), Physics:-Ket(F, t)) = X[i, j](t)

Setup(%Bracket(Physics[Bra](C, i, j), Physics[Ket](F, t)) = X[i, j](t))

[bracketrules = {%Bracket(%Bra(C, i, j), %Ket(F, t)) = X[i, j](t)}]

Now you have

Bra(C, a, b).Ket(F, tau)

X[a, b](tau)

And with that,

Bra(C, a, b).`𝕀__F`

Int(X[a, b](t)*Physics:-Bra(F, t), t = -infinity .. infinity)

 

This answers the question on how to do what I understood you were trying to do. Note next that you can set some rules that have to do with C states being tensor products of A and B states, For example, if

 

Ket(C, j, k) = Ket(A, j)*Ket(B, k)

Physics:-Ket(C, j, k) = Physics:-`*`(Physics:-Ket(A, j), Physics:-Ket(B, k))

Then the bracket rules would be

Bra(A, i).%, Bra(B, i).%

Physics:-Bracket(Physics:-Bra(A, i), Physics:-Ket(C, j, k)) = Physics:-KroneckerDelta[i, j]*Physics:-Ket(B, k), Physics:-Bracket(Physics:-Bra(B, i), Physics:-Ket(C, j, k)) = Physics:-KroneckerDelta[i, k]*Physics:-Ket(A, j)

Set these rules

Setup(%)

[bracketrules = {%Bracket(%Bra(A, i), %Ket(C, j, k)) = Physics:-KroneckerDelta[i, j]*Physics:-Ket(B, k), %Bracket(%Bra(B, i), %Ket(C, j, k)) = Physics:-KroneckerDelta[i, k]*Physics:-Ket(A, j), %Bracket(%Bra(C, i, j), %Ket(F, t)) = X[i, j](t)}]

Now you have

Bra(A, m)*Bra(B, n).Ket(C, r, s)

Physics:-KroneckerDelta[n, s]*Physics:-KroneckerDelta[m, r]

And hence

Bra(A, n).`𝕀__C`

Sum(Physics:-`*`(Physics:-Ket(B, j), Physics:-Bra(C, n, j)), j = 0 .. Nb-1)

Bra(B, n).`𝕀__C`

Sum(Physics:-`*`(Physics:-Ket(A, i), Physics:-Bra(C, i, n)), i = 0 .. Na-1)

If you tell how would A and B Bras contract with F Kets that can also be implemented as a bracket rule, after which you can use them and Simplify will take them into account when simplifying expressions. Without these rules the system doesn't know how to compute these other ones

Bra(A, n).`𝕀__F`

Int(Physics:-Bracket(Physics:-Bra(A, n), Physics:-Ket(F, t))*Physics:-Bra(F, t), t = -infinity .. infinity)

Bra(B, n).`𝕀__F`

Int(Physics:-Bracket(Physics:-Bra(B, n), Physics:-Ket(F, t))*Physics:-Bra(F, t), t = -infinity .. infinity)

 

I'm not sure I interpreted your question correctly. If the above does not answer your question I'd need more details on the problem you are trying to formulate in order to suggest you something else.

NULL


 

Download algebra_rules_(reviewed).mw

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


 

Hi

This problem, that can be considered as one with discontinuous boundary conditions (see my previous reply pointing at their classical incompatibility) is now solved.

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

BC := u(0, y) = 0, u(Pi, y) = sinh(Pi)*cos(y), u(x, 0) = sin(x), u(x, Pi) = -sinh(x)

pdsolve([PDE, BC])

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

(1)

The improvement is available to everybody within the Maplesoft Physics Updates version 190 or higher. (Note that since its version 188 the Physics Updates requires Maple 2018.2, the free upgrade for those having Maple 2018.)

 

Actually, there is a large number of improvements in exact solutions to PDE & BC problems. I mentioned one in a previous post on Sturm-Liouville problems with RootOf eigenvalues. I will try to summarize next week the new types of problems that can now be solved in maple 2018.2 and beyond.


 

Download weak_solution_-_problem_with_discontinuous_boundary_conditions.mw

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

Hi; John is correct: besides the most relevant tensors of General Relativity, the idea in Physics is to provide flexible ways for you to define, using Physics :-Define (it is worth looking at its help page), any tensor. That is, in any dimension, even a tensor in spacetime, in only space, or in the tetrad system of references (for that purpose pass the corresponding indices in the definition sent to Define). That, of course, allows you for defining not just the Cotton tensor but any tensor that makes sense in your formulation or intermediate computations. 

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

Hi

This problem was one in MultiSet, not in DifferentialGeometry. Anyway, it is fixed now, and the fix uploaded for everybody, available within the Maplesoft Physics Updates version 172 or higher. With the Updates installed, you get


 

restart

with(DifferentialGeometry); with(Tensor); with(LieAlgebras); with(GroupActions); with(Library)

DGsetup([x, y], M)

g := evalDG((`&t`(dx, dx)+`&t`(dy, dy))/y^2)

_DG([["tensor", M, [["cov_bas", "cov_bas"], []]], [[[1, 1], 1/y^2], [[2, 2], 1/y^2]]])

(1)

kv := KillingVectors(g)

[_DG([["vector", M, []], [[[1], (1/2)*x^2-(1/2)*y^2], [[2], y*x]]]), _DG([["vector", M, []], [[[1], x], [[2], y]]]), _DG([["vector", M, []], [[[1], 1]]])]

(2)

L := LieAlgebraData(kv, Poincare)

_DG([["LieAlgebra", Poincare, [3]], [[[1, 2, 1], -1], [[1, 3, 2], -1], [[2, 3, 3], -1]]])

(3)

DGsetup(L); MultiplicationTable("LieTable")

Query(Poincare, "Semisimple")

true

(4)

CSA := CartanSubalgebra(Poincare)

[_DG([["vector", Poincare, []], [[[2], 1]]])]

(5)

RSD := RootSpaceDecomposition(CSA)

table( [( [-1] ) = _DG([["vector", Poincare, []], [[[3], 1]]]), ( [1] ) = _DG([["vector", Poincare, []], [[[1], 1]]]) ] )

(6)

NULL


 

Download error_(reviewed).mw

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

Hi,

Your coordinates are [x, y, z, t], and the line element of your metric is 4*dx*dt-dy^2-dz^2+4*u*(dt^2)). Look closer ... the line elements does not contain the coordinates. Therefore all the derivatives of the metric with respect to the coordinates are equal to zero, so all the Christoffel symbols (proportional to those derivatives) too, and with it, all of the Ricci, Riemman and Weyl tensors are equal to 0.

In the UPDATE I read that u is a function of the variables (not specified ... I imagine you meant to say coordinates x,y,z,t?). If that is the case, definitely this is the problem. To answer that question, give it a try, enter change u by u(x,y,z,t) in your input, press enter, and you will see that all the general relativity tensors, starting with Christoffel, now have nonzero components.

And that resolves the issue. Tip: use CompactDisplay for u(X), and you can type u(X), no need to type u(x,y,z,t).

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


 

Hi

There exists no solution for this problem. The boundary conditions are incompatible. Your PDE and BC are

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

BC := u(0, y) = 0, u(Pi, y) = sinh(Pi)*cos(y), u(x, 0) = sin(x), u(x, Pi) = -sinh(x)


Always check the compatibility of the BC before proceeding, evaluating them at the different evaluation points you see in them, in this case x = 0, x = Pi, y = 0, y = Pi, which does not exhaust the compatibility test but frequently helps detecting an issue. So evaluate BC, for instance, at y = 0, x = Pi

eval([BC], [y = 0, x = Pi])

[u(0, 0) = 0, u(Pi, 0) = sinh(Pi), u(Pi, 0) = 0, u(Pi, Pi) = -sinh(Pi)]

(1)

These two

[u(0, 0) = 0, u(Pi, 0) = sinh(Pi), u(Pi, 0) = 0, u(Pi, Pi) = -sinh(Pi)][2 .. 3]

[u(Pi, 0) = sinh(Pi), u(Pi, 0) = 0]

(2)

cannot be both true:

evalf([u(Pi, 0) = sinh(Pi), u(Pi, 0) = 0])

[u(Pi, 0) = 11.54873936, u(Pi, 0) = 0.]

(3)

So I'd say the answer someone posted before given by Mathematica (now I see it deleted?) can never be correct - in that there exists no solution for this problem, and the refusal of Maple's pdsolve/numeric to return something is correct, and Maple's pdsolve exact solutions returning NULL is also correct.

 


 

Download incompatible_boundary_conditions.mw

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


 

with(Physics); Setup(mathematicalnotation = true)

In your input, where you wrote

metric = f(dx1^2+dx2^2+dx3^2+dx4^2)

 

you are saying that the metric is an arbitrary function of the sum of the squares of the differentials of the coordinates. That is not a diagonal metric as I understand you intended to say. It is also not a line element quadratic in the differential of the coordinates (because it an arbitrary function f of that), and it is also not a tensorial expression with two free indices. That is what the error message is telling you, but for one typo at the end: it should say 'none', instead of repeating [%d_(x1), %d_(x2), %d_(x3), %d_(x4)].

 

So I'd suggest first for you to update your library with the Maplesoft Physics Updates where this typo in the error message is already corrected, and for your input you will see:

Setup(Coordinatesystem = (X = [x1, x2, x3, x4]), metric = f(dx1^2+dx2^2+dx3^2+dx4^2))

`* Partial match of  '`*Coordinatesystem*`' against keyword '`*coordinatesystems*`' `

 

`Default differentiation variables for d_, D_ and dAlembertian are:`*{X = (x1, x2, x3, x4)}

 

`Systems of spacetime Coordinates are:`*{X = (x1, x2, x3, x4)}

 

Error, (in Physics:-NormalizeMetricInfo) expected definition of a metric as an expression quadratic in [%d_(x1), %d_(x2), %d_(x3), %d_(x4)] or a tensorial algebraic expression with two free indices; received one that is not quadratic in [%d_(x1), %d_(x2), %d_(x3), %d_(x4)] and with free indices: none

 

Now the message makes sense. A comment to the side: the default way of setting coordinates automatically uses x1, x2, x3, x4 when the label is X ;so, you can but you do not need to indicate each coordinate as you did, Setup(coordinates = X, metric = ds^2) is sufficient.

 

On to your question: how do you enter a diagonal metric where the elements in the diagonal are arbitrary?

 

Just write the line element. As we know, it is an expression quadratic in the differentials of the coordinates, which in the context of Setup can be entered as the letter d followed by the coordinate (for example x1) with no spaces in between, as you do with paper and pencil

Parse:-ConvertTo1D, "first argument to _Inert_ASSIGN must be assignable"

f[1](X)*dx1^2+f[2](X)*dx2^2+f[3](X)*dx3^2+f[4](X)*dx4^2

(1)

Use CompactDisplay to avoid redundant display of functionality (not necessary, but for me it is always desired)

CompactDisplay(f[1](X)*dx1^2+f[2](X)*dx2^2+f[3](X)*dx3^2+f[4](X)*dx4^2)

` f`(X)*`will now be displayed as`*f

(2)

Now set the metric

Setup(metric = ds^2)

[metric = {(1, 1) = f[1](X), (2, 2) = f[2](X), (3, 3) = f[3](X), (4, 4) = f[4](X)}]

(3)

That answers your question.

 

 

What follows are some comments about syntax and functionality, things that I find useul to know or take into account.

 

When setting the metric, check it out, to be sure you got what you intended to get; these are the covariant and contravariant components

g_[]

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

(4)

"g_[~]"

Physics:-g_[`~mu`, `~nu`] = Matrix(%id = 18446744078251369406)

(5)

The functionality of the f__j functions is not displayed, but it is there. To see it for the previous expression use show

show

Physics:-g_[`~mu`, `~nu`] = Matrix(%id = 18446744078251369406)

(6)

That is the usefulness of CompactDisplay. You can also set the compact display ON and OFF around the worksheet as you see it convenient by just entering ON and OFF; see the help page CompactDisplay .

 

Next, in (1) I entered the line element writing it explicitly. Here is a version of it writing less that works the same way

Parse:-ConvertTo1D, "first argument to _Inert_ASSIGN must be assignable"

f[1](X)*Physics:-d_(x1)^2+f[2](X)*Physics:-d_(x2)^2+f[3](X)*Physics:-d_(x3)^2+f[4](X)*Physics:-d_(x4)^2

(7)

In the above the d(X[j]) is an full computational representation of the differential (see the help page for d_ )

 

You  can also indicate the metric as a symmetric matrix. In your example, that could be

piecewise(mu = nu, h[mu](X), 0)

piecewise(mu = nu, h[mu](X), 0)

(8)

I am using a different letter, h only to make evident whether the new definition is kicking in

g_[mu, nu] = Matrix(4, proc (mu, nu) options operator, arrow; piecewise(mu = nu, h[mu](X), 0) end proc)

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

(9)

"Setup(?)"

[metric = {(1, 1) = h[1](X), (2, 2) = h[2](X), (3, 3) = h[3](X), (4, 4) = h[4](X)}]

(10)

g_[]

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

(11)

Finally, the coordinates x1, x2, x3, x4 although ok, look a bit cumbersome to me. I prefer to use the literal subscripted

Coordinates(X = [x__1, x__2, x__3, x__4])

`Default differentiation variables for d_, D_ and dAlembertian are:`*{X = (x__1, x__2, x__3, x__4)}

 

`Systems of spacetime Coordinates are:`*{X = (x__1, x__2, x__3, x__4)}

 

{X}

(12)

or even just Cartesian coordinates

Coordinates(X = cartesian)

`Default differentiation variables for d_, D_ and dAlembertian are:`*{X = (x, y, z, t)}

 

`Systems of spacetime Coordinates are:`*{X = (x, y, z, t)}

 

{X}

(13)

and as you see you can define the coordinates using Setup or using Coordinates. Likewise, you can define the metric using Setup or using Define.

 

If you redefine the coordinates in the middle of the computation as I did above, you need to redefine the metric, that it is still a function of the previous coordinates

g_[]

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

(14)

To redefine the metric just re-enter the lines (8), (9), (10).

NULL


 

Download Diagonal_metric.mw

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

Joe already answered regarding the problem in your DifferentialGeometry input. Although a great package, it suffers of this issue where the input syntax is not entirely natural, or with the natural meaning in context, eg. in the context of evalDG, dr dr could be dr &t dr, after you have said that [r, theta] are coordinates in DGsetup - this frequently bites people.

An alternative to these subtleties is to use the Physics package, more modern and free of those issues. This is the same computation with it:
 

with(Physics)

Setup(dimension = 2, coordinates = (X = [r, theta]), metric = r^2*`dθ`^2+dr^2, quiet)

[coordinatesystems = {X}, dimension = 2, metric = {(1, 1) = 1, (2, 2) = r^2}]

(1)

Christoffel[nonzero]

Physics:-Christoffel[alpha, mu, nu] = {(1, 2, 2) = -r, (2, 1, 2) = r, (2, 2, 1) = r}

(2)

seq(Christoffel[j, mu, nu, matrix], j = [1, 2])

Physics:-Christoffel[1, mu, nu] = Matrix(%id = 18446744078349940366), Physics:-Christoffel[2, mu, nu] = Matrix(%id = 18446744078349941446)

(3)

 

The main feature in the above being that everything you write is in natural physics notation (basically, as you would write it with paper and pencil).

Download 2D_spacetime_in_polar_coordinates.mw

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

You are correct vv. This is already fixed and the fix uploaded for everybody within the Maplesoft Physics Updates (version 146 or higher).

(On the way I implemented copy & paste of intats, something pending that never worked. And being at that, also copy and paste of ODESolStruc, PDESolStruc and `casesplit/ans`, none of these worked before, they are used in the output of dsolve, pdsolve, casesplit and for exact solutions to PDE & Boundary Conditions when the eigenvalues have RootOf form - see recent post.)

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

Hi, check the help page for PDEtools:-PolynomialSolutions, it is the command you are looking for. There are also some related commands mentioned in that page. Another command with impressive functionality, that I believe could be of use in your case is PDEtools:-InvariantSolutions. That command allows you to specify that you want to use only polynomial infinitesimals - the simpest case would be linear - and also to specify their degree, all in all resulting in a computational process that is not expensive and leading to all kinds of concrete and simple solutions (when they exist), like polynomial ones.

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

@Preben Alsholm 

This problem you mention, dsolve's solution

ode:=diff(x(t),t)=piecewise(t=1,0,x(t));
sol:=dsolve({ode,x(0)=1});

was a problem in convert/Heaviside, not actually in dsolve. By the way, I implemented, in dsolve, not convert, a workaround, so that the class of problem (in convert/Heaviside) behind this example doesn't hit dsolve's solutions. That change/workaround is included for everybody in the Maplesoft Physics Updates for 2018.1 since some versions ago (current version is 119).

NOTE after posting: oops! I see you mentioned this is a problem in convert/Heaviside too. Anyway this problem is not hitting dsolve anymore.


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


 

This issue, the unexpected presence of undefined in the solutions, is fixed, and the fix is available to everybody within the Physics Updates (version 99) at the Maplesoft R&D Physics webpage. With that Update in place all the PDE & BC problems posse are solved exactly.  Below I am including a verification of the solutions using pdetest on the PDE and the boundary conditions.

 

This is the input you posted:

 

restart

with(PDEtools); umain := 1-exp(-y)+sum(A^i*u[i](y, z), i = 1 .. 5); vmain := -1+sum(A^i*v[i](y, z), i = 1 .. 5); wmain := sum(A^i*w[i](y, z), i = 1 .. 5); pde[main] := diff(vmain, y)+diff(wmain, z); pde[main[2]] := vmain*(diff(wmain, y))+wmain*(diff(wmain, z))-(diff(wmain, `$`(y, 2)))

1-exp(-y)+A*u[1](y, z)+A^2*u[2](y, z)+A^3*u[3](y, z)+A^4*u[4](y, z)+A^5*u[5](y, z)

 

-1+A*v[1](y, z)+A^2*v[2](y, z)+A^3*v[3](y, z)+A^4*v[4](y, z)+A^5*v[5](y, z)

 

A*w[1](y, z)+A^2*w[2](y, z)+A^3*w[3](y, z)+A^4*w[4](y, z)+A^5*w[5](y, z)

 

A*(diff(v[1](y, z), y))+A^2*(diff(v[2](y, z), y))+A^3*(diff(v[3](y, z), y))+A^4*(diff(v[4](y, z), y))+A^5*(diff(v[5](y, z), y))+A*(diff(w[1](y, z), z))+A^2*(diff(w[2](y, z), z))+A^3*(diff(w[3](y, z), z))+A^4*(diff(w[4](y, z), z))+A^5*(diff(w[5](y, z), z))

 

(-1+A*v[1](y, z)+A^2*v[2](y, z)+A^3*v[3](y, z)+A^4*v[4](y, z)+A^5*v[5](y, z))*(A*(diff(w[1](y, z), y))+A^2*(diff(w[2](y, z), y))+A^3*(diff(w[3](y, z), y))+A^4*(diff(w[4](y, z), y))+A^5*(diff(w[5](y, z), y)))+(A*w[1](y, z)+A^2*w[2](y, z)+A^3*w[3](y, z)+A^4*w[4](y, z)+A^5*w[5](y, z))*(A*(diff(w[1](y, z), z))+A^2*(diff(w[2](y, z), z))+A^3*(diff(w[3](y, z), z))+A^4*(diff(w[4](y, z), z))+A^5*(diff(w[5](y, z), z)))-A*(diff(diff(w[1](y, z), y), y))-A^2*(diff(diff(w[2](y, z), y), y))-A^3*(diff(diff(w[3](y, z), y), y))-A^4*(diff(diff(w[4](y, z), y), y))-A^5*(diff(diff(w[5](y, z), y), y))

(1)

sys[1] := [A*coeff(pde[main[2]], A) = 0, w[1](0, z) = sin(z), w[1](infinity, z) = 0]

[A*(-(diff(w[1](y, z), y))-(diff(diff(w[1](y, z), y), y))) = 0, w[1](0, z) = sin(z), w[1](infinity, z) = 0]

(2)

First solution by pdsolve is OK

pdsolve(sys[1])

w[1](y, z) = sin(z)*exp(-y)

(3)

pdetest(w[1](y, z) = sin(z)*exp(-y), sys[1])

[0, 0, 0]

(4)

w[1] := proc (y, z) options operator, arrow, function_assign; sin(z)*exp(-y) end proc; pdex[1] := diff(w[1](y, z), z)+diff(v[1](y, z), y) = 0; pdexbound := v[1](0, z) = 0; sysx[1] := [pdex[1], pdexbound]

proc (y, z) options operator, arrow, function_assign; sin(z)*exp(-y) end proc

 

cos(z)*exp(-y)+diff(v[1](y, z), y) = 0

 

v[1](0, z) = 0

 

[cos(z)*exp(-y)+diff(v[1](y, z), y) = 0, v[1](0, z) = 0]

(5)

This solution is also correct

pdsolve(sysx[1])

v[1](y, z) = cos(z)*(exp(-y)-1)

(6)

pdetest(v[1](y, z) = cos(z)*(exp(-y)-1), sysx[1])

[0, 0]

(7)

v[1] := proc (y, z) options operator, arrow, function_assign; cos(z)*(exp(-y)-1) end proc

proc (y, z) options operator, arrow, function_assign; cos(z)*(exp(-y)-1) end proc

(8)

I prefer to present the following loop using 1D input to be more transparent: I am adding a verification using pdetest for each solution computed with pdsolve. There is no undefined around, and every solution is correct, matching the PDE and the initial conditions. So all the PDE & BC problems formulated can and are solved exactly.

for i from 2 to 3 do
    `_____________________________________________________________________________`;
    pde[i] := coeff(pde[main[2]], A^i) = 0;
    sys[i] := [pde[i], w[i](0, z) = 0, w[i](infinity, z) = 0];
    sol[i] := pdsolve(sys[i]);
    verification = pdetest(sol[i], sys[i]);
    assign(op(0, lhs(pdsolve(sys[i]))) = unapply(rhs(pdsolve(sys[i])), op(lhs(pdsolve(sys[i])))));

    `_____________________________________________________________________________`;
    sysx[i] := [diff(w[i](y, z), z)+diff(v[i](y, z), y) = 0, v[i](0, z) = 0];
    solx[i] := pdsolve(sysx[i]);
    verification = pdetest(solx[i], sysx[i]);

    assign(op(0, lhs(pdsolve(sysx[i]))) = unapply(rhs(pdsolve(sysx[i])), op(lhs(pdsolve(sysx[i])))))
end do;

_____________________________________________________________________________

 

-(diff(w[2](y, z), y))-cos(z)*(exp(-y)-1)*sin(z)*exp(-y)+sin(z)*(exp(-y))^2*cos(z)-(diff(diff(w[2](y, z), y), y)) = 0

 

[-(diff(w[2](y, z), y))-cos(z)*(exp(-y)-1)*sin(z)*exp(-y)+sin(z)*(exp(-y))^2*cos(z)-(diff(diff(w[2](y, z), y), y)) = 0, w[2](0, z) = 0, w[2](infinity, z) = 0]

 

w[2](y, z) = -(1/2)*exp(-y)*sin(2*z)*y

 

verification = [0, 0, 0]

 

_____________________________________________________________________________

 

[-exp(-y)*cos(2*z)*y+diff(v[2](y, z), y) = 0, v[2](0, z) = 0]

 

v[2](y, z) = -cos(2*z)*(exp(-y)*y+exp(-y)-1)

 

verification = [0, 0]

 

_____________________________________________________________________________

 

-(diff(w[3](y, z), y))+cos(z)*(exp(-y)-1)*((1/2)*exp(-y)*sin(2*z)*y-(1/2)*exp(-y)*sin(2*z))+cos(2*z)*(exp(-y)*y+exp(-y)-1)*sin(z)*exp(-y)-sin(z)*(exp(-y))^2*cos(2*z)*y-(1/2)*(exp(-y))^2*sin(2*z)*y*cos(z)-(diff(diff(w[3](y, z), y), y)) = 0

 

[-(diff(w[3](y, z), y))+cos(z)*(exp(-y)-1)*((1/2)*exp(-y)*sin(2*z)*y-(1/2)*exp(-y)*sin(2*z))+cos(2*z)*(exp(-y)*y+exp(-y)-1)*sin(z)*exp(-y)-sin(z)*(exp(-y))^2*cos(2*z)*y-(1/2)*(exp(-y))^2*sin(2*z)*y*cos(z)-(diff(diff(w[3](y, z), y), y)) = 0, w[3](0, z) = 0, w[3](infinity, z) = 0]

 

w[3](y, z) = (1/8)*((y^2+4*y-1)*sin(3*z)+sin(z)*(y-1)*(y-3))*exp(-y)-(3/8)*exp(-2*y)*(sin(z)-(1/3)*sin(3*z))

 

verification = [0, 0, 0]

 

_____________________________________________________________________________

 

[(1/8)*(3*(y^2+4*y-1)*cos(3*z)+cos(z)*(y-1)*(y-3))*exp(-y)-(3/8)*exp(-2*y)*(cos(z)-cos(3*z))+diff(v[3](y, z), y) = 0, v[3](0, z) = 0]

 

v[3](y, z) = (1/16)*((6*y^2+36*y+30)*cos(3*z)+2*cos(z)*(y-1)^2)*exp(-y)+(1/16)*(3*exp(-2*y)-33)*cos(3*z)-(3/16)*exp(-2*y)*cos(z)+(1/16)*cos(z)

 

verification = [0, 0]

(9)

``


 

Download NoUndefinedAndExactSolutionsAreCorrect.mw

 

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



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

Hi, your input contains d*(t^2), it should be dt^2, and that is all. To the side: note that, in Maple, gamma is a constant, so using gamma as a function is prone to a problem at some point, so I'd consider replacing that letter by another letter.

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

First 29 30 31 32 33 34 35 Last Page 31 of 59