Ask Your Question
0

Misusing sage.manifolds?

asked 2018-03-07 20:42:50 +0200

Richard_L gravatar image

updated 2018-03-08 18:57:24 +0200

Two misbehaviors are apparent in the following: 1) objects not displaying in declared default frame, and 2) unexpected zero after applying Laplace-Beltrami operator to a scalar function in one of the two charts/frames. Q: Am I setting up the charts and frames incorrectly? Q: Missing something else? (I have heavily elided the sagemath output.) Thank you.

Preliminaries:

sage: version()
'SageMath version 8.2.beta7, Release Date: 2018-03-03'
sage: Parallelism().set(nproc=2)
sage: Manifold.options.omit_function_arguments=True
sage: from sage.manifolds.utilities import ExpressionNice
sage: var('rho12,rho13,rho23', domain='real')
sage: assume(rho12>0, rho13>0, rho23>0)
sage: var('r12,r13,r23', domain='real')
sage: assume(r12>0, r13>0, r23>0, r23<r12+r13, r13<r12+r23, r12<r13+r23)
sage: var('m1 m2 m3', domain='real')
sage: var('mu12,mu13,mu23', domain='real')
sage: assume(m1>0, m2>0, m3>0)
sage: m1=1; m2=1; m3=1
sage: mu12 = (m1*m2)/(m1+m2); mu23 = (m2*m3)/(m2+m3); mu13 = (m1*m3)/(m1+m3)
sage: rho12 = r12^2; rho13 = r13^2; rho23 = r23^2

sage: Ginv = matrix([[2*rho12/mu12, 1/m1*(rho12+rho13-rho23), 1/m2*(rho12+rho23-rho13)],\
...           [1/m1*(rho12+rho13-rho23), 2*rho13/mu13, 1/m3*(rho13+rho23-rho12)],\
...           [1/m2*(rho12+rho23-rho13), 1/m3*(rho13+rho23-rho12), 2*rho23/mu23]])
...
sage: Ginv
[               4*r12^2  r12^2 + r13^2 - r23^2  r12^2 - r13^2 + r23^2]
[ r12^2 + r13^2 - r23^2                4*r13^2 -r12^2 + r13^2 + r23^2]
[ r12^2 - r13^2 + r23^2 -r12^2 + r13^2 + r23^2                4*r23^2]

sage: G = Ginv.inverse(); #G.simplify_full()

Establish manifold, Riemann metric, charts and frames...

sage: M = Manifold(3,'R^3',field='real',start_index=1)
sage: U = M.open_subset('U')
sage: Rho.<r12,r13,r23> = U.chart("r12:(0,+oo) r13:(0,+oo) r23:(0,+oo)")
sage: Rho.add_restrictions([r23<r12+r13, r13<r12+r23, r12<r13+r23])
sage: g = M.riemannian_metric('g');
sage: g[:] = G[:].simplify_full()
sage: ginv = g.inverse(); ginv.display()
inv_g = 4*r12^2 d/dr12*d/dr12 + (r12^2 + r13^2 - r23^2) d/dr12*d/dr13 + (r12^2 - r13^2 + r23^2) d/dr12*d/dr23 + (r12^2 + r13^2 - r23^2) d/dr13*d/dr12 + 4*r13^2 d/dr13*d/dr13 + (-r12^2 + r13^2 + r23^2) d/dr13*d/dr23 + (r12^2 - r13^2 + r23^2) d/dr23*d/dr12 + (-r12^2 + r13^2 + r23^2) d/dr23*d/dr13 + 4*r23^2 d/dr23*d/dr23

sage: Rho.domain(); Rho.parent()
Open subset U of the 3-dimensional differentiable manifold R^3
<class 'sage.manifolds.differentiable.chart.RealDiffChart'>

sage: # Pekeris coordinates
sage: Tau.<u1,u2,u3> = U.chart("u1:(0,+oo) u2:(0,+oo) u3:(0,+oo)")
sage: Rho_Tau = Rho.transition_map(Tau, ((r12+r13-r23)/2, (r12-r13+r23)/2, (-r12+r13+r23)/2))
sage: Tau_Rho = Rho_Tau.inverse()
sage: U.set_default_chart(Tau)
sage: U.set_default_frame(Tau.frame())
sage: Rho_Tau.display(); Tau_Rho.display();
u1 = 1/2*r12 + 1/2*r13 - 1/2*r23
u2 = 1/2*r12 - 1/2*r13 + 1/2*r23
u3 = -1/2*r12 + 1/2*r13 + 1/2*r23
r12 = u1 + u2
r13 = u1 + u3
r23 = u2 + u3

sage: U.atlas(); Rho.frame(); Tau.frame()
[Chart (U, (r12, r13, r23)), Chart (U, (u1, u2, u3))]
Coordinate frame (U, (d/dr12,d/dr13,d/dr23))
Coordinate frame (U, (d/du1,d/du2,d/du3))

sage: g.display_comp(frame=Tau.frame())
g_u1,u1 = 1/48*(3*u2^4 + 12*u2*u3^3 + 3*u3^4 + 2*(8*u1*u2 + 9*u2^2)*u3^2 + 4*(4*u1^2*u2 + 4*u1*u2^2 + 3*u2^3)*u3)/(u1*u2*u3^4 + 2*(u1^2*u2 + u1*u2^2)*u3^3 + (2*u1^3*u2 + 3*u1^2*u2^2 + 2*u1*u2^3)*u3^2 + (u1^4*u2 + 2*u1^3*u2^2 + 2*u1^2*u2^3 + u1*u2^4)*u3) 
... [elided]
g_u3,u3 = 1/48*(3*u1^4 + 12*u1^3*u2 + 18*u1^2*u2^2 + 12*u1*u2^3 + 3*u2^4 + 16*u1*u2*u3^2 + 16*(u1^2*u2 + u1*u2^2)*u3)/(u1*u2*u3^4 + 2*(u1^2*u2 + u1*u2^2)*u3^3 + (2*u1^3*u2 + 3*u1^2*u2^2 + 2*u1*u2^3)*u3^2 + (u1^4*u2 + 2*u1^3*u2^2 + 2*u1^2*u2^3 + u1*u2^4)*u3)

First problem:

sage: #Why does this cell not compute and display in the Tau frame (declared default above)?
sage: ginv = g.inverse(); ginv.display_comp()
inv_g^r12,r12 = 4*u1^2 + 8*u1*u2 + 4*u2^2 
inv_g^r12,r13 = 2*u1^2 + 2*u1*u2 + 2*(u1 - u2)*u3 
inv_g^r12,r23 = 2*u1*u2 + 2*u2^2 - 2*(u1 - u2)*u3 
inv_g^r13,r12 = 2*u1^2 + 2*u1*u2 + 2*(u1 - u2)*u3 
inv_g^r13,r13 = 4*u1^2 + 8*u1*u3 + 4*u3^2 
inv_g^r13,r23 = -2*u1*u2 + 2*(u1 + u2)*u3 + 2*u3^2 
inv_g^r23,r12 = 2*u1*u2 + 2*u2^2 - 2*(u1 - u2)*u3 
inv_g^r23,r13 = -2*u1*u2 + 2*(u1 + u2)*u3 + 2*u3^2 
inv_g^r23,r23 = 4*u2^2 + 8*u2*u3 + 4*u3^2 

sage: #Display in Tau frame can be forced.
sage: ginv.display_comp(frame=Tau.frame())
inv_g^u1,u1 = 3*u1^2 + 3*u1*u2 + u2^2 + (3*u1 - u2)*u3 + u3^2 
inv_g^u1,u2 = u1*u2 - (u1 + u2)*u3 - u3^2 
inv_g^u1,u3 = -u1*u2 - u2^2 + (u1 - u2)*u3 
inv_g^u2,u1 = u1*u2 - (u1 + u2)*u3 - u3^2 
inv_g^u2,u2 = u1^2 + 3*u1*u2 + 3*u2^2 - (u1 - 3*u2)*u3 + u3^2 
inv_g^u2,u3 = -u1^2 - u1*u2 - (u1 - u2)*u3 
inv_g^u3,u1 = -u1*u2 - u2^2 + (u1 - u2)*u3 
inv_g^u3,u2 = -u1^2 - u1*u2 - (u1 - u2)*u3 
inv_g^u3,u3 = u1^2 - u1*u2 + u2^2 + 3*(u1 + u2)*u3 + 3*u3^2 

sage: nabla = g.connection()
sage: # nabla also not displayed in Tau frame!
sage: nabla.display()
Gam^r12_r12,r12 = -1/6*(8*r12^11 - r12^10*r13 - 8*r12^9*r13^2 + 3*r12^8*r13^3 - 24*r12^7*r13^4 - 2*r12^6*r13^5 + 40*r12^5*r13^6 - 2*r12^4*r13^7 - 16*r12^3*r13^8 + 3*r12^2*r13^9 - r13^11 + 13*r13*r23^10 - r23^11 + (3*r12^2 + 29*r13^2)*r23^9 - (16*r12^3 + 29*r12^2*r13 + 96*r12*r13^2 + 17*r13^3)*r23^8 - 2*(r12^4 - 4*r12^2*r13^2 + 23*r13^4)*r23^7 + 2*(20*r12^5 + 5*r12^4*r13 + 32*r12^3*r13^2 + 12*r12^2*r13^3 + 48*r12*r13^4 + 11*r13^5)*r23^6 - 2*(r12^6 + 7*r12^4*r13^2 + 3*r12^2*r13^4 - 11*r13^6)*r23^5 - 2*(12*r12^7 - 7*r12^6*r13 - 28*r12^5*r13^2 - 3*r12^4*r13^3 - 144*r12^3*r13^4 + 3*r12^2*r13^5 - 48*r12*r13^6 + 23*r13^7)*r23^4 + (3*r12^8 - 16*r12^6*r13^2 + 6*r12^4*r13^4 + 24*r12^2*r13^6 - 17*r13^8)*r23^3 - (8*r12^9 + 7*r12^8*r13 + 16*r12^7*r13^2 + 16*r12^6*r13^3 - 56*r12^5*r13^4 + 14*r12^4*r13^5 - 64*r12^3*r13^6 - 8*r12^2*r13^7 + 96*r12*r13^8 - 29*r13^9)*r23^2 - (r12^10 + 7*r12^8*r13^2 - 14*r12^6*r13^4 - 10*r12^4*r13^6 + 29*r12^2*r13^8 - 13*r13^10)*r23)/(r12^12 - 2*r12^10*r13^2 - r12^8*r13^4 + 4*r12^6*r13^6 - r12^4*r13^8 - 2*r12^2*r13^10 + r13^12 + r23^12 - 2*(r12^2 + r13^2)*r23^10 - (r12^4 + 10*r12^2*r13^2 + r13^4)*r23^8 + 4*(r12^6 + 3*r12^4*r13^2 + 3*r12^2*r13^4 + r13^6)*r23^6 - (r12^8 - 12*r12^6*r13^2 - 42*r12^4*r13^4 - 12*r12^2*r13^6 + r13^8)*r23^4 - 2*(r12^10 + 5*r12^8*r13^2 - 6*r12^6*r13^4 - 6*r12^4*r13^6 + 5*r12^2*r13^8 + r13^10)*r23^2) 
... [elided]

Second problem: [Appears to me to be incomplete conversion to the Tau frame.]

sage: psi = U.scalar_field({Rho: function('Psi')(r12,r13,r23)}, name='psi',latex_name='\Psi'); psi
Scalar field psi on the Open subset U of the 3-dimensional differentiable manifold R^3

sage: LBpsi = nabla(nabla(psi).up(g)).trace(); LBpsi
Scalar field on the Open subset U of the 3-dimensional differentiable manifold R^3

sage: DelPsi = LBpsi.expr().factor(); ExpressionNice(DelPsi)
1/2*(8*u1^6*u2*u3*d^2(Psi)/d(u1 + u2)^2 + 32*u1^5*u2^2*u3*d^2(Psi)/d(u1 + u2)^2 + 56*u1^4*u2^3*u3*d^2(Psi)/d(u1 + u2)^2 + 56*u1^3*u2^4*u3*d^2(Psi)/d(u1 + u2)^2 + 32*u1^2*u2^5*u3*d^2(Psi)/d(u1 + u2)^2 + 8*u1*u2^6*u3*d^2(Psi)/d(u1 + u2)^2 + 16*u1^5*u2*u3^2*d^2(Psi)/d(u1 + u2)^2 + 56*u1^4*u2^2*u3^2*d^2(Psi)/d(u1 + u2)^2 + 
... [elided]
u3) + 4*u2^2*u3^5*d(Psi)/d(u2 + u3) - u1*u3^6*d(Psi)/d(u2 + u3) + u2*u3^6*d(Psi)/d(u2 + u3))/((u1^2 + u1*u2 + u2^2 + u1*u3 + u2*u3 + u3^2)*(u1 + u2 + u3)*u1*u2*u3)

Leading to the unexpected zero:

sage: psi = U.scalar_field({Rho: function('Psi')(u1,u2,u3)}, name='psi',latex_name='\Psi'); psi
Scalar field psi on the Open subset U of the 3-dimensional differentiable manifold R^3

sage: LBpsi = nabla(nabla(psi).up(g)).trace(); LBpsi
Scalar field zero on the Open subset U of the 3-dimensional differentiable manifold R^3
edit retag flag offensive close merge delete

1 Answer

Sort by ยป oldest newest most voted
1

answered 2018-03-08 00:23:39 +0200

eric_g gravatar image

The first problem arises because the domain of ginv is M and not U (since g is defined as g = M.riemannian_metric('g')) and the default frame of M is Rho.frame(). So, you should set:

sage: M.set_default_frame(Tau.frame())

The second problem lies in the definition of psi:

sage: psi = U.scalar_field({Rho: function('Psi')(u1,u2,u3)}, name='psi',latex_name='\Psi')

You are declaring the coordinate expression of psi in the chart Rho, but (u1,u2,u3) are coordinates of the chart Tau, so that function('Psi')(u1,u2,u3) appears as a constant function. For instance, one has diff(function('Psi')(u1,u2,u3), r12) = 0. Hence the Laplace-Beltrami operator applied to psi yields zero.

edit flag offensive delete link more

Comments

Thank you, Eric! Silly oversights on my part. For simplicity, I have decided to define everything on U instead of M. Further, I now remove the set_default_..., making frame and chart explicit wherever .display() or .expr() allows.

Richard_L gravatar imageRichard_L ( 2018-03-08 18:55:06 +0200 )edit

I have also changed this post's title.

Richard_L gravatar imageRichard_L ( 2018-03-08 18:57:53 +0200 )edit

You're welcome.

eric_g gravatar imageeric_g ( 2018-03-08 22:20:18 +0200 )edit

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower

Stats

Asked: 2018-03-07 20:42:50 +0200

Seen: 262 times

Last updated: Mar 08 '18