Advanced Search

ecterrab

13431 Reputation

24 Badges

19 years, 361 days
Contact ecterrab

MaplePrimes Activity


These are answers submitted by ecterrab

Physics:-Define...

ecterrab 13431
November 02 2018
3 0

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

Fixed in the Physics Updates version 172...

ecterrab 13431
November 01 2018
4 1

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

look closer at your metric ......

ecterrab 13431
October 23 2018
3 1

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

The boundary conditions are incompatible...

ecterrab 13431
October 22 2018
2 3


 

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

Diagonal metric...

ecterrab 13431
October 21 2018
3 0


 

> 

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

Physics...

ecterrab 13431
October 19 2018
6 0

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

Fixed....

ecterrab 13431
October 17 2018
3 1

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

PDEtools:-PolynomialSolutions...

ecterrab 13431
September 19 2018
3 3

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 is a problem i...

ecterrab 13431
September 18 2018
2 0

@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

Fixed, no undefined around. All of the P...

ecterrab 13431
September 05 2018
3 0


 

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

Physics:-Vectors package...

ecterrab 13431
September 01 2018
2 0



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

dt^2, not d*(t^2)...

ecterrab 13431
August 25 2018
3 0

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

Fixed....

ecterrab 13431
August 13 2018
5 5


 

Back in town, this one is fixed, the fix available for download on the MapleCloud, within the Maplesoft Physics Updates version 83. The problem is solvable, the method is explained on page 10 of Polyanin's book "Linear Partial Differential Equations for Engineers and Scientists", although this solution is not straightforward to use in that it involves a symbolic composition of a differential operator.

 

This is the PDE and its solution, first convert the floating point numbers to exact numbers (this is about exact solutions)

> 

PDE := convert(diff(u(t, x, y), t)+.5*(diff(u(t, x, y), x, x))+diff(u(t, x, y), x, y)+.5*(diff(u(t, x, y), y, y))+diff(u(t, x, y), x)+diff(u(t, x, y), y)+(1+(x+y)^2)*u(t, x, y) = 0, rational)

diff(u(t, x, y), t)+(1/2)*(diff(diff(u(t, x, y), x), x))+diff(diff(u(t, x, y), x), y)+(1/2)*(diff(diff(u(t, x, y), y), y))+diff(u(t, x, y), x)+diff(u(t, x, y), y)+(1+(x+y)^2)*u(t, x, y) = 0

 (1)
> 

bc := u(0, x, y) = 1

u(0, x, y) = 1

 (2)
> 

pdsolve([PDE, bc])

u(t, x, y) = 1+Sum(t^n*((proc (U) options operator, arrow; -(1/2)*(diff(diff(U, x), x))-(diff(diff(U, x), y))-(1/2)*(diff(diff(U, y), y))-(diff(U, x))-(diff(U, y))-(1+(x+y)^2)*U end proc)@@n)(1)/factorial(n), n = 1 .. infinity)

 (3)

By eye this solution satisfies the condition at t = 0, but testing it against the PDE itself is not straightforward and in fact beyond what pdetest can do today:

> 

pdetest(u(t, x, y) = 1+Sum(t^n*((proc (U) options operator, arrow; -(1/2)*(diff(diff(U, x), x))-(diff(diff(U, x), y))-(1/2)*(diff(diff(U, y), y))-(diff(U, x))-(diff(U, y))-(1+(x+y)^2)*U end proc)@@n)(1)/factorial(n), n = 1 .. infinity), [PDE, bc])

[sum(t^(n-1)*((proc (U) options operator, arrow; -(1/2)*(diff(diff(U, x), x))-(diff(diff(U, x), y))-(1/2)*(diff(diff(U, y), y))-(diff(U, x))-(diff(U, y))-(1+(x+y)^2)*U end proc)@@n)(1)*(t*x^2+2*t*x*y+t*y^2+n+t)/factorial(n), n = 1 .. infinity)+x^2+2*x*y+y^2+1, 0]

 (4)

The idea behind the method, however is simple: take diff(u(t, x, y), t) as the left-hand side of PDE

> 

pde := isolate(PDE, diff(u(t, x, y), t))

diff(u(t, x, y), t) = -(1/2)*(diff(diff(u(t, x, y), x), x))-(diff(diff(u(t, x, y), x), y))-(1/2)*(diff(diff(u(t, x, y), y), y))-(diff(u(t, x, y), x))-(diff(u(t, x, y), y))-(1+(x+y)^2)*u(t, x, y)

 (5)

You see that the right-hand side does not involve t (but through u(t, x, y)). Rewrite the right-hand side as a differential operator M that depends only on x and y (not  t)

> 

diff(u(t, x, y), t) = M(u(t, x, y))

diff(u(t, x, y), t) = M(u(t, x, y))

 (6)

where

> 

M := unapply(subs(u(t, x, y) = U, rhs(diff(u(t, x, y), t) = -(1/2)*(diff(diff(u(t, x, y), x), x))-(diff(diff(u(t, x, y), x), y))-(1/2)*(diff(diff(u(t, x, y), y), y))-(diff(u(t, x, y), x))-(diff(u(t, x, y), y))-(1+(x+y)^2)*u(t, x, y))), U)

proc (U) options operator, arrow; -(1/2)*(diff(diff(U, x), x))-(diff(diff(U, x), y))-(1/2)*(diff(diff(U, y), y))-(diff(U, x))-(diff(U, y))-(1+(x+y)^2)*U end proc

 (7)

This PDE problem, without any boundary condition, admits a particular solution in terms of an arbitrary function _F1(x, y) (see Polyanin's book)

> 

sol := u(t, x, y) = _F1(x, y)+Sum(t^n*(M@@n)(_F1(x, y))/factorial(n), n = 1 .. infinity)

u(t, x, y) = _F1(x, y)+Sum(t^n*(M@@n)(_F1(x, y))/factorial(n), n = 1 .. infinity)

 (8)

By eye this solution also satisfies u(0, x, y) = _F1(x, y), from where (3) is a solution to the system formed by PDE and bc.

> 

``


 

Download PDE_solution_(composed_differential_operator).mw

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

Interesting suggestion ......

ecterrab 13431
July 30 2018
1 3

Hi Adam,

It is not possible, but it is an interesting suggestion. But depending on the scope of this - e.g. is it for you or everybody - writing a friendly interface for this is not just a snap of fingers. A relevant feature/design of the FunctionAdvisor is that it gathers information from several places in the library (e.g. simplify, diff, expand, ...), and it is not just a database; for instance the DE is calculated not just for a function but for arbitrary compositions of them, using PDEtools:-dpolyform, etc.

Anyway, it is a valuable suggestion - in the end, we are adding functions almost at every release, and a simple interface would also simplify development. I will think about.

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

Physics...

ecterrab 13431
July 26 2018
4 0


 

> 

with(Physics)

Set the coordinates

> 

Coordinates(X = [u, v, y, z])

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

  

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

  

{X}

 (1)

Enter the line element

> 

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

du*dv+F(y, z)*du^2+dy^2+dz^2

 (2)

Set the metric using this line element

> 

Setup(metric = du*dv+F(y, z)*du^2+dy^2+dz^2)

[metric = {(1, 1) = F(y, z), (1, 2) = 1/2, (3, 3) = 1, (4, 4) = 1}]

 (3)

That is all, these are the nonzero Christoffel symbol

> 

Christoffel[nonzero]

Physics:-Christoffel[alpha, mu, nu] = {(1, 1, 3) = (1/2)*(diff(F(y, z), y)), (1, 1, 4) = (1/2)*(diff(F(y, z), z)), (1, 3, 1) = (1/2)*(diff(F(y, z), y)), (1, 4, 1) = (1/2)*(diff(F(y, z), z)), (3, 1, 1) = -(1/2)*(diff(F(y, z), y)), (4, 1, 1) = -(1/2)*(diff(F(y, z), z))}

 (4)
> 

Christoffel[`~alpha`, mu, nu, nonzero]

Physics:-Christoffel[`~alpha`, mu, nu] = {(2, 1, 3) = diff(F(y, z), y), (2, 1, 4) = diff(F(y, z), z), (2, 3, 1) = diff(F(y, z), y), (2, 4, 1) = diff(F(y, z), z), (3, 1, 1) = -(1/2)*(diff(F(y, z), y)), (4, 1, 1) = -(1/2)*(diff(F(y, z), z))}

 (5)

See the help page for Physics:-Christoffel

 

To see the nonzero components of the Ricci tensor, it is the same

> 

Ricci[nonzero]

Physics:-Ricci[mu, nu] = {(1, 1) = -(1/2)*(diff(diff(F(y, z), y), y))-(1/2)*(diff(diff(F(y, z), z), z))}

 (6)

Since this is a 2-tensor you can directly see all the components in matricial form

> 

Ricci[]

Physics:-Ricci[mu, nu] = Matrix(%id = 18446744078352828286)

 (7)
> 

"Ricci[~]"

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

 (8)

For the Ricci scalar

> 

Ricci[scalar]

0

 (9)

For the five Ricci scalars in the context of the Newman-Penrose formalism it is the same. to have a more compact display, use

> 

CompactDisplay(F(y, z))

` F`(y, z)*`will now be displayed as`*F

 (10)

So now the functionality of F is there but not displayed, and derivatives are displayed indexed, for example

> 

"Ricci[~]"

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

 (11)

Try now the scalars

> 

Ricci[scalars]

Phi__00 = -(diff(diff(F(y, z), y), y)+diff(diff(F(y, z), z), z))*(F(y, z)+(F(y, z)^2+1)^(1/2))/(8*F(y, z)*(F(y, z)^2+1)^(1/2)+8*F(y, z)^2+8), Phi__01 = (1/8)*(diff(diff(F(y, z), y), y)+diff(diff(F(y, z), z), z))*(F(y, z)+(F(y, z)^2+1)^(1/2))^(1/2)*(-F(y, z)+(F(y, z)^2+1)^(1/2))^(1/2)/((F(y, z)*(F(y, z)^2+1)^(1/2)+F(y, z)^2+1)^(1/2)*(F(y, z)*(F(y, z)^2+1)^(1/2)-F(y, z)^2-1)^(1/2)), Phi__02 = -(diff(diff(F(y, z), y), y)+diff(diff(F(y, z), z), z))*(-F(y, z)+(F(y, z)^2+1)^(1/2))/(8*F(y, z)*(F(y, z)^2+1)^(1/2)-8*F(y, z)^2-8), Phi__11 = (1/8)*(diff(diff(F(y, z), y), y)+diff(diff(F(y, z), z), z))*(F(y, z)^2+1)^(1/2)/((F(y, z)*(F(y, z)^2+1)^(1/2)+F(y, z)^2+1)*(F(y, z)*(F(y, z)^2+1)^(1/2)-F(y, z)^2-1)), Phi__12 = (1/8)*(diff(diff(F(y, z), y), y)+diff(diff(F(y, z), z), z))*(F(y, z)+(F(y, z)^2+1)^(1/2))^(1/2)*(-F(y, z)+(F(y, z)^2+1)^(1/2))^(1/2)/((F(y, z)*(F(y, z)^2+1)^(1/2)+F(y, z)^2+1)^(1/2)*(F(y, z)*(F(y, z)^2+1)^(1/2)-F(y, z)^2-1)^(1/2)), Phi__22 = -(diff(diff(F(y, z), y), y)+diff(diff(F(y, z), z), z))*(F(y, z)+(F(y, z)^2+1)^(1/2))/(8*F(y, z)*(F(y, z)^2+1)^(1/2)+8*F(y, z)^2+8), Lambda = 0

 (12)

And the same is for any of the tensors used in general relativity or one that you can define; see the Physics:-Define  help page

> 

``


 

Download Christoffel_symbols_and_Ricci_tensor.mw

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

First 25 26 27 28 29 30 31 Last Page 27 of 55