ecterrab

13431 Reputation

24 Badges

19 years, 362 days

MaplePrimes Activity


These are answers submitted by ecterrab

Hi Shillis

This kind of operation is simpler to perform using the Physics package. Suppose the metric is given by some matrix M that is 4 x 4, or by the square of a line element ds2. Then enter Setup(metric = M) or Setup(metric = ds2) to set the metric. You could also set the coordinates you are using together, for instance as in Setup(coordinates = [x,y,z,t], metric = M).

Then just use Christoffel[1,1,1] to get the all covariant 1,1,1 component you want to pick out, or Christoffel[~1,1,1] to get the component where the 1st index is contravariant (so prefixed with ~). You see that the somehow artificial distinction between 1st and 2nd kind dissapears, you just indicate whether the 1st index is covariant or contravariant.  Check the help page ?Christoffel for examples, it is pretty straightforward. To see for instance the matrix corresponding to Christoffel[3, mu, nu] enter Christoffel[3, mu,nu, matrix], to see the definition (only in Maple 2015) enter Christoffel[definition], or if you have a previous version of Maple check the definition directly in the help page for Physics:-Christoffel. For example:

with(Physics):

M := Matrix(4, 4, {(1, 1) = r/(-r+2*m), (2, 2) = -r^2, (3, 3) = -r^2*sin(theta)^2, (4, 4) = (r-2*m)/r}, shape = symmetric)

M := Matrix(4, 4, {(1, 1) = r/(-r+2*m), (1, 2) = 0, (1, 3) = 0, (1, 4) = 0, (2, 1) = 0, (2, 2) = -r^2, (2, 3) = 0, (2, 4) = 0, (3, 1) = 0, (3, 2) = 0, (3, 3) = -r^2*sin(theta)^2, (3, 4) = 0, (4, 1) = 0, (4, 2) = 0, (4, 3) = 0, (4, 4) = (r-2*m)/r})

(1)

Physics:-Setup(coordinates = spherical, metric = M)

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

 

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

 

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

 

[coordinatesystems = {X}, metric = {(1, 1) = r/(-r+2*m), (2, 2) = -r^2, (3, 3) = -r^2*sin(theta)^2, (4, 4) = (r-2*m)/r}]

(2)

Christoffel[1, 1, 1]

m/(-r+2*m)^2

(3)

Christoffel[`~1`, 1, 1]

m/(r*(-r+2*m))

(4)

"seq(Christoffel[j, alpha, beta, matrix], j=[~1, ~2, ~3, ~4]);"

Physics:-Christoffel[`~1`, alpha, beta] = Matrix(%id = 18446744078475670222), Physics:-Christoffel[`~2`, alpha, beta] = Matrix(%id = 18446744078607044726), Physics:-Christoffel[`~3`, alpha, beta] = Matrix(%id = 18446744078607046766), Physics:-Christoffel[`~4`, alpha, beta] = Matrix(%id = 18446744078607036414)

(5)

 

 

Download PickOutChristoffelComponents.mw

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

Hi zia

I gave a look at your worksheet, where you use RootFinding:-NextZero. There is one limitation in that command: if, as the second argument, you send to it something containing 0.0*I, it will refuse to work, while if you simplify that second argument, this spurious imaginary part dissapears. So, first thing: instead of 

R[k] := RootFinding:-NextZero(K, R[k-1], maxdistance = 10000, guardDigits = 3)

use

R[k] := RootFinding:-NextZero(K, simplify(R[k-1]), maxdistance = 10000, guardDigits = 3)

The second limitation is that RootFinding:-NextZero can return FAIL, and then your conditional 

for k to 20 while R[k-1]<10000 do 

won't work, because the system cannot determine wheter FAIL < 10000. Use instead 

for k to 20 while evalb(R[k-1]<10000) = true do 

Putting altogether, use the following and you do get the results you are looking for:

for k to 20 while evalb(R[k-1] < 10000) = true do
    R[k] := RootFinding:-NextZero(K, simplify(R[k-1]), maxdistance = 10000, guardDigits = 3) e
od

=============================================

Edit: originally evalb(R[k-1] < 10000), corrected to evalb(R[k-1] < 10000) = true

=============================================

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

Hi Inacio

Note that in d_[mu](A[nu]) you are applying an operator, i.e. in computer language A[nu] is the argument of a 'function' called d_[mu]. On the other hand, in A[nu] d_[mu] you are multiplying two objects, i.e. in computer language this is the`*` (product) operator being used. These two operations are different, of course, and the Commutator, actually computes the `.` (dot operator `.`) product, that depending on the case (arguments passed to Commutator) is transformed into the `*` or Vectors:-`.` (scalar vector product) operator.

To have Commutator behaving as you were expecting we would need to decide, on the fly, when to apply or multiply two objectswhich may make Commutator a bit unpredictable, although in the example you are posted (one of the arguments is a recognized differential operator) this can be implemented. I will think about this and post again here in this thread.

Meantime, the answer from all the above is: no, you cannot use Commutator with the purpose you are mentioning.

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

@Michael_Watson 

Please note that I cannot compare an image (that is what you posted) with a Maple expression (as in a worksheet, that is what I asked), unless I read your image and try to retype everything, which is time consumming and prone to errors. In any case attached is the worksheet you posted originally, without changing anything, just executed via menus (Edit -> Execute -> Worksheet). In this reviewed worksheet you see that (19) is not 0 = 0 as you posted originally. From this reviewed worksheet you can also copy the result (19) and paste it into the worksheet you used to obtain the image you presented in your reply, then try simplifying one expression minus the other one to see if if you get 0 or what. Hope it helps.

Question_algsubs_3.27.15_(reviewed).mw

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

Hi,

In the eqdN you show there are parenthesis that are open and not closed (e.g. in eqd3 and eqd4), or you open with `(` but close with `}`. You need to correct that first, I suggest you entering one equation at a time and only entering the next one after you are sure the previous one is correctly entered.

Then, do not enclose each equation within {}; In Maple, {} is used represent sets. Copying from the error message you posted (it only contains three equations, not four), correct input for your system of equations is:

> sys := {diff(v(x), x) = V(x), 1/250*x*(diff(U(x), x))+(x-8*x*V(x)+3/500)*U(x)+(3/2-12*V(x)-8*x*(diff(V(x), x)))*u(x) = 0, 1/125*x*V(x)^2+2*(3/5+x)*v(x)*V(x)+v(x)*u(x)+v(x)^2 = 0}
and for this system you get:

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


Hi mpre,

You cannot force to find the result you want, but the result you get is basically the same:

with(Physics); -1; t := Physics:-Intc(Physics:-`*`(Physics:-`*`(Physics:-`*`(Dirac(k1+k2+k3), phi(k1)), phi(k2)), phi(k3)), k1, k2, k3)

Int(Int(Int(Dirac(k1+k2+k3)*phi(k1)*phi(k2)*phi(k3), k1 = -infinity .. infinity), k2 = -infinity .. infinity), k3 = -infinity .. infinity)

(1)

Physics:-Fundiff(t, phi(-k));

2*(Int(phi(k3)*phi(-k3+k), k3 = -infinity .. infinity))+Int(phi(-k2+k)*phi(k2), k2 = -infinity .. infinity)

(2)

Here combine  and Simplify  accomplish the same

Physics:-Simplify(2*(Int(phi(k3)*phi(-k3+k), k3 = -infinity .. infinity))+Int(phi(-k2+k)*phi(k2), k2 = -infinity .. infinity))

Int(3*phi(k3)*phi(-k3+k), k3 = -infinity .. infinity)

(3)

The result above is is the same you want to find but for one integration:

Physics:-`*`(3, Physics:-Intc(Physics:-`*`(Physics:-`*`(Dirac(k2+k3-k), phi(k2)), phi(k3)), k2, k3));

3*(Int(Int(Dirac(k2+k3-k)*phi(k2)*phi(k3), k2 = -infinity .. infinity), k3 = -infinity .. infinity))

(4)

value(3*(Int(Int(Dirac(k2+k3-k)*phi(k2)*phi(k3), k2 = -infinity .. infinity), k3 = -infinity .. infinity)))

3*(int(phi(k3)*phi(-k3+k), k3 = -infinity .. infinity))

(5)

You see that (5) is equal to (3). The difference between (3) and (4) is then that in (3) all the Dirac functions got integrated, that is the two Dirac functions that appear through functional differentiation multiplied by the original Dirac(k1+k2+k3) already found in the integrand of t. In the result you want to find ((4)) you seem to want one of the three Dirac functions to remain not integrated. There is no way to tell the code to integrate two of them but not the third one.

 

Now, just out of curiosity, why do want to integrate two Dirac functions but not the third one?


Download Fundiff_Intc.mw

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


Hi jschulzb,

 

The input lines that follow are taken from your worksheet:

restart:
with(Physics):
Setup(mathematicalnotation = true, dimension = 2):
Coordinates(X):

`The dimension and signature of the tensor space are set to: [2, - +] `

 

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

 

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

(1)

# Method 1: Using Define and various differential operators
Define(z):
z :=sqrt(R^2-X[mu]*X[mu]);
# d_[mu](z(X));
# d_[1](z(X))
diff(z, x1);
# diff(z, X[mu]);

`Defined objects with tensor properties`

 

(R^2-Physics:-SpaceTimeVector[mu](X)^2)^(1/2)

 

-Physics:-g_[1, mu]*Physics:-SpaceTimeVector[`~mu`](X)/(R^2-Physics:-SpaceTimeVector[mu](X)^2)^(1/2)

(2)

From these four differentiations you input, the first two, in comments, are not correct: enter z(X) and you see it is not what you meant:

z(X)

(R(X)^2-(Physics:-SpaceTimeVector[mu](X))(X)^2)^(1/2)

(3)

Then the third input, diff(z, x1), not in comments, returns a correct result (2). The fourth input, diff(z, X[mu]), in comments, is correct input, but the output is off by a factor 2 as you say, and the same happens with your Method #2, again quoting from your worksheet:

# Method #2: Using functions

z2 := mu -> sqrt(R^2-X[mu]*X[mu]);

proc (mu) options operator, arrow; sqrt(Physics:-`^`(R, 2)-Physics:-`*`((X)[mu], (X)[mu])) end proc

(4)

Note that your z in your Method 1 is equal to your z2(mu)in your Method 2, so we have only one issue, not two.

z = z2(mu)

(R^2-Physics:-SpaceTimeVector[mu](X)^2)^(1/2) = (R^2-Physics:-SpaceTimeVector[mu](X)^2)^(1/2)

(5)

Now, the wrong factor 2 is due to a combination of things. First, in diff(sqrt(R^2-X[mu]*X[mu]), X[mu]) you repeat the index mu twice (Einstein's rule requires that contracted indices are not repeated more than once). Had you entered diff(sqrt(R^2-X[nu]*X[nu]), X[mu]) that represents the same object, the wrong factor 2 would not be there. Second, the code scans for these incorrect repetition of indices before proceeding, but the first repetition of mu is within a function (sqrt) and for that reason the code missed the occurrence (it shouldn't).

 

This issue is fixed. To have the fix installed in your computer you need to update your Physics Library with the one available for download on the Maplesoft R&D Physics webpage. With the update, you have, from your Method 1:

diff(z, X[mu])

-Physics:-SpaceTimeVector[`~mu`](X)/(R^2-Physics:-SpaceTimeVector[nu](X)^2)^(1/2)

(6)

 


Download PhysicsDiffBug_(reviewed).mw

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


Suppose your function

f(x) = 1/(x^2+1)^2

f(x) = 1/(x^2+1)^2

(1)

There is actually no computer algebra function but a rational expression on the right-hand side (this is a subtlety of language, but relevant in what follows), so first transform this expression into an actual function: for all rational expressions, a conversion to a hypergeometric function is possible as in:

convert(f(x) = 1/(x^2+1)^2, hypergeom, include = powers)

f(x) = hypergeom([2], [], -hypergeom([-2], [], 1-x))

(2)

At this point, you can transform this into a differential equation, via

PDEtools:-dpolyform(f(x) = hypergeom([2], [], -hypergeom([-2], [], 1-x)), no_Fn)

`casesplit/ans`([diff(diff(f(x), x), x) = (3/2)*(diff(f(x), x))^2/f(x)+(diff(f(x), x))/x], [diff(f(x), x) <> 0, (diff(f(x), x))*x+4*f(x) <> 0])

(3)

In this way you have expressed "the derivative of a function or expression in terms of the function (and its derivatives) itself". By the way, verify that this result (3) is correct, that is: that this result relates the function and its deriatives

odetest(f(x) = 1/(x^2+1)^2, `casesplit/ans`([diff(diff(f(x), x), x) = (3/2)*(diff(f(x), x))^2/f(x)+(diff(f(x), x))/x], [diff(f(x), x) <> 0, (diff(f(x), x))*x+4*f(x) <> 0]))

[0]

(4)

You can do this basically for all the functions of the mathematical language (but for the GAMMA function and a few related ones), or for any mathematical expression that involes mathematical functions, powers (fractional or not) and their compositions, etc. For example:

f(x) = exp(exp(x^(1/2)))/x^2

f(x) = exp(exp(x^(1/2)))/x^2

(5)

PDEtools:-dpolyform(f(x) = exp(exp(x^(1/2)))/x^2, no_Fn)

`casesplit/ans`([(diff(diff(f(x), x), x))^2 = (2*(diff(f(x), x))^2/f(x)-(diff(f(x), x))/x+2*f(x)/x^2)*(diff(diff(f(x), x), x))-(diff(f(x), x))^4/f(x)^2+(diff(f(x), x))^3/(x*f(x))+(1/4)*(x-9)*(diff(f(x), x))^2/x^2+(x+1)*f(x)*(diff(f(x), x))/x^3+(x-1)*f(x)^2/x^4], [2*(diff(diff(f(x), x), x))*f(x)*x^2-2*(diff(f(x), x))^2*x^2+(diff(f(x), x))*f(x)*x-2*f(x)^2 <> 0])

(6)

odetest(f(x) = exp(exp(x^(1/2)))/x^2, `casesplit/ans`([(diff(diff(f(x), x), x))^2 = (2*(diff(f(x), x))^2/f(x)-(diff(f(x), x))/x+2*f(x)/x^2)*(diff(diff(f(x), x), x))-(diff(f(x), x))^4/f(x)^2+(diff(f(x), x))^3/(x*f(x))+(1/4)*(x-9)*(diff(f(x), x))^2/x^2+(x+1)*f(x)*(diff(f(x), x))/x^3+(x-1)*f(x)^2/x^4], [2*(diff(diff(f(x), x), x))*f(x)*x^2-2*(diff(f(x), x))^2*x^2+(diff(f(x), x))*f(x)*x-2*f(x)^2 <> 0]))

[0]

(7)

You can also do this with functions of many variables: the result expresses the partial derivatives in terms of the function and the other partial derivatives typically in the form of a PDE system. For example:

f(x, y) = tan(y^2+x)

f(x, y) = tan(y^2+x)

(8)

PDEtools:-dpolyform(f(x, y) = tan(y^2+x), no_Fn)

`casesplit/ans`([diff(f(x, y), x) = f(x, y)^2+1, diff(f(x, y), y) = 2*f(x, y)^2*y+2*y], [f(x, y)^2+1 <> 0, f(x, y) <> 0])

(9)

pdetest(f(x, y) = tan(y^2+x), `casesplit/ans`([diff(f(x, y), x) = f(x, y)^2+1, diff(f(x, y), y) = 2*f(x, y)^2*y+2*y], [f(x, y)^2+1 <> 0, f(x, y) <> 0]))

[0, 0]

(10)

``


Download connection_between_functions_and_their_derivatives.mw

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

Enter your expression, highlight the input, right click and: 2D-math -> convert to -> atomic identifier, and that produces what you want, that is an arbitrary expression as a singlie variable name. This image illustrates the final result:

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

Hi Trace

There is a problem with the input syntax, you are: missing a closing parenthesis, have extra spaces after d in 'd theta' and 'd phi' (should be dtheta and dphi), and you missed a multiplication sign after r in the last term of the square of the line element. Please give a look at the revised worksheet attached, where all this is adjusted, and the nonzero components of the Ricci tensor shown.

non_zero_components_of_ricci_tensor_(reviewed).mw

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

Hi Michael
Unfortunately, there seems to be a GUI (Graphical User Interface) issue around, related to using the default document mode. Concretely, if instead of expand((13)) you open a prompt and execute > expand((13)); the result is not 0 = 0. See the attached worksheet.

I will give a further look tomorrow to see if there is a programatic (library) way to workaround this GUI issue until it is fixed.

Question_R12_(reviewed).mw

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

Hi Codge

You know, the sections plotted are projections up to 3 dimension while the phase space has more than 3 dimensions; there is no programatic way to restrict the value of one of the variables (dimensions) out of sight - in your example, p1. I didn't investigate your problem further but it springs to my mind: why don't you analyze the problem trying to determine the range of t such that p1 > 0? Then you can specify that range instead of using t = -5000..5000 as you are using. For the purpose of this investigation, you can write Hamilton equations (you can use DEtools[hamilton_eqs]), then a numerical or exact solution of the system of equations will let you visually or algebraically determine the range of t for which p1 > 0.

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

Hi

In the post on "Computer Algebra for Theoretical Physics" you see how to define a tensor whose components are the partial derivatives of something else, or a matrix, or just defined by an algebraic tensorial expression. Check please the "Examples" in the sections for "Classical Field Theory" and "General Relativity". In another previous post, "Physics, one year of developments", in the section entitled "Tensors", within the examples, equation (48) shows what you ask, the definition of a tensor in terms of the partial derivatives of another tensor. The first section on Simplification, what concerns tensors, also contains material similar to what you ask. Finally in the help page ?Physics,Examples that comes with the Maple system, there is a section on "Tensors in Special and General Relativity" that contains answers to your question and more details on how to work with tensors in Maple in general.

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

There is an issue with Dirac and these equations that I'd need to give a look at the internals, not for today. Meantime, try as shown in this image and it resolves the problem; i.e. let dsolve tackle it without knowing who is 'delta'; use the option explicit to get the solution built (see ?dsolve,details where this is expllained):

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


Hi MathyMaple

This is again a combination of premature evaluation problems. The first one is in the power operation. In Maple, debatable or not, we have

0^0

1

(1)

But then,

0^m

0

(2)

Naturally, if in a second step m assume the value0, "it is too late". But this is resolved when you load Physics, that among other things also redefines the power operation to handle noncommutative operands. So,

with(Physics)

0^m

Physics:-`^`(0, m)

(3)

:) So that issue is not there anymore. The second premature evaluation problem happens in sum, but that is also addressed within the Physics package, that allows for redefining sum: just in order for it to be free of these evaluation issues.

Setup(redefinesum = true)

[redefinesum = true]

(4)

Try now your example:

X := [0, 1, 2, 3]

[0, 1, 2, 3]

(5)

n := 2

2

(6)

sum(sum(a[k]*X[i]^k, k = 0 .. n), i = 1 .. nops(X))

4*a[0]+6*a[1]+14*a[2]

(7)

sum(sum(a[k]*X[i]^k, i = 1 .. nops(X)), k = 0 .. n)

4*a[0]+6*a[1]+14*a[2]

(8)

To have the above working you need to update your Physics library with the one distributed in the Maplesoft R&D webpage.


Download DoubleSumAndPowerEvaluationProblem.mw

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

First 38 39 40 41 42 43 44 Last Page 40 of 55