Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
r""" Pseudo-Riemannian Metrics
The class :class:`PseudoRiemannianMetric` implements pseudo-Riemannian metrics on differentiable manifolds over `\RR`. The derived class :class:`PseudoRiemannianMetricParal` is devoted to metrics with values on a parallelizable manifold.
AUTHORS:
- Eric Gourgoulhon, Michal Bejger (2013-2015) : initial version - Pablo Angulo (2016): Schouten, Cotton and Cotton-York tensors
REFERENCES:
- [KN1963]_ - [Lee1997]_ - [ONe1983]_
""" #****************************************************************************** # Copyright (C) 2015 Eric Gourgoulhon <eric.gourgoulhon@obspm.fr> # Copyright (C) 2015 Michal Bejger <bejger@camk.edu.pl> # Copyright (C) 2016 Pablo Angulo <pang@cancamusa.net> # # Distributed under the terms of the GNU General Public License (GPL) # as published by the Free Software Foundation; either version 2 of # the License, or (at your option) any later version. # http://www.gnu.org/licenses/ #******************************************************************************
r""" Pseudo-Riemannian metric with values on an open subset of a differentiable manifold.
An instance of this class is a field of nondegenerate symmetric bilinear forms (metric field) along a differentiable manifold `U` with values on a differentiable manifold `M` over `\RR`, via a differentiable mapping `\Phi: U \rightarrow M`. The standard case of a metric field *on* a manifold corresponds to `U=M` and `\Phi = \mathrm{Id}_M`. Other common cases are `\Phi` being an immersion and `\Phi` being a curve in `M` (`U` is then an open interval of `\RR`).
A *metric* `g` is a field on `U`, such that at each point `p\in U`, `g(p)` is a bilinear map of the type:
.. MATH::
g(p):\ T_q M\times T_q M \longrightarrow \RR
where `T_q M` stands for the tangent space to the manifold `M` at the point `q=\Phi(p)`, such that `g(p)` is symmetric: `\forall (u,v)\in T_q M\times T_q M, \ g(p)(v,u) = g(p)(u,v)` and nondegenerate: `(\forall v\in T_q M,\ \ g(p)(u,v) = 0) \Longrightarrow u=0`.
.. NOTE::
If `M` is parallelizable, the class :class:`PseudoRiemannianMetricParal` should be used instead.
INPUT:
- ``vector_field_module`` -- module `\mathfrak{X}(U,\Phi)` of vector fields along `U` with values on `\Phi(U)\subset M` - ``name`` -- name given to the metric - ``signature`` -- (default: ``None``) signature `S` of the metric as a single integer: `S = n_+ - n_-`, where `n_+` (resp. `n_-`) is the number of positive terms (resp. number of negative terms) in any diagonal writing of the metric components; if ``signature`` is ``None``, `S` is set to the dimension of manifold `M` (Riemannian signature) - ``latex_name`` -- (default: ``None``) LaTeX symbol to denote the metric; if ``None``, it is formed from ``name``
EXAMPLES:
Standard metric on the sphere `S^2`::
sage: M = Manifold(2, 'S^2', start_index=1) sage: # The two open domains covered by stereographic coordinates (North and South): sage: U = M.open_subset('U') ; V = M.open_subset('V') sage: M.declare_union(U,V) # S^2 is the union of U and V sage: c_xy.<x,y> = U.chart() ; c_uv.<u,v> = V.chart() # stereographic coord sage: xy_to_uv = c_xy.transition_map(c_uv, (x/(x^2+y^2), y/(x^2+y^2)), ....: intersection_name='W', restrictions1= x^2+y^2!=0, ....: restrictions2= u^2+v^2!=0) sage: uv_to_xy = xy_to_uv.inverse() sage: W = U.intersection(V) # The complement of the two poles sage: eU = c_xy.frame() ; eV = c_uv.frame() sage: c_xyW = c_xy.restrict(W) ; c_uvW = c_uv.restrict(W) sage: eUW = c_xyW.frame() ; eVW = c_uvW.frame() sage: g = M.metric('g') ; g Riemannian metric g on the 2-dimensional differentiable manifold S^2
The metric is considered as a tensor field of type (0,2) on `S^2`::
sage: g.parent() Module T^(0,2)(S^2) of type-(0,2) tensors fields on the 2-dimensional differentiable manifold S^2
We define g by its components on domain U (factorizing them to have a nicer view)::
sage: g[eU,1,1], g[eU,2,2] = 4/(1+x^2+y^2)^2, 4/(1+x^2+y^2)^2 sage: g.display(eU) g = 4/(x^2 + y^2 + 1)^2 dx*dx + 4/(x^2 + y^2 + 1)^2 dy*dy
A matrix view of the components::
sage: g[eU,:] [4/(x^2 + y^2 + 1)^2 0] [ 0 4/(x^2 + y^2 + 1)^2]
The components of g on domain V expressed in terms of (u,v) coordinates are similar to those on domain U expressed in (x,y) coordinates, as we can check explicitly by asking for the component transformation on the common subdomain W::
sage: g.display(eVW, c_uvW) g = 4/(u^4 + v^4 + 2*(u^2 + 1)*v^2 + 2*u^2 + 1) du*du + 4/(u^4 + v^4 + 2*(u^2 + 1)*v^2 + 2*u^2 + 1) dv*dv
Therefore, we set::
sage: g[eV,1,1], g[eV,2,2] = 4/(1+u^2+v^2)^2, 4/(1+u^2+v^2)^2 sage: g[eV,1,1].factor() ; g[eV,2,2].factor() 4/(u^2 + v^2 + 1)^2 4/(u^2 + v^2 + 1)^2 sage: g.display(eV) g = 4/(u^2 + v^2 + 1)^2 du*du + 4/(u^2 + v^2 + 1)^2 dv*dv
At this stage, the metric is fully defined on the whole sphere. Its restriction to some subdomain is itself a metric (by default, it bears the same symbol)::
sage: g.restrict(U) Riemannian metric g on the Open subset U of the 2-dimensional differentiable manifold S^2 sage: g.restrict(U).parent() Free module T^(0,2)(U) of type-(0,2) tensors fields on the Open subset U of the 2-dimensional differentiable manifold S^2
The parent of `g|_U` is a free module because is `U` is a parallelizable domain, contrary to `S^2`. Actually, `g` and `g|_U` have different Python type::
sage: type(g) <class 'sage.manifolds.differentiable.metric.PseudoRiemannianMetric'> sage: type(g.restrict(U)) <class 'sage.manifolds.differentiable.metric.PseudoRiemannianMetricParal'>
As a field of bilinear forms, the metric acts on pairs of tensor fields, yielding a scalar field::
sage: a = M.vector_field('a') sage: a[eU,:] = [x, 2+y] sage: a.add_comp_by_continuation(eV, W, chart=c_uv) sage: b = M.vector_field('b') sage: b[eU,:] = [-y, x] sage: b.add_comp_by_continuation(eV, W, chart=c_uv) sage: s = g(a,b) ; s Scalar field g(a,b) on the 2-dimensional differentiable manifold S^2 sage: s.display() g(a,b): S^2 --> R on U: (x, y) |--> 8*x/(x^4 + y^4 + 2*(x^2 + 1)*y^2 + 2*x^2 + 1) on V: (u, v) |--> 8*(u^3 + u*v^2)/(u^4 + v^4 + 2*(u^2 + 1)*v^2 + 2*u^2 + 1)
The inverse metric is::
sage: ginv = g.inverse() ; ginv Tensor field inv_g of type (2,0) on the 2-dimensional differentiable manifold S^2 sage: ginv.parent() Module T^(2,0)(S^2) of type-(2,0) tensors fields on the 2-dimensional differentiable manifold S^2 sage: latex(ginv) g^{-1} sage: ginv.display(eU) # again the components are expanded inv_g = (1/4*x^4 + 1/4*y^4 + 1/2*(x^2 + 1)*y^2 + 1/2*x^2 + 1/4) d/dx*d/dx + (1/4*x^4 + 1/4*y^4 + 1/2*(x^2 + 1)*y^2 + 1/2*x^2 + 1/4) d/dy*d/dy sage: ginv.display(eV) inv_g = (1/4*u^4 + 1/4*v^4 + 1/2*(u^2 + 1)*v^2 + 1/2*u^2 + 1/4) d/du*d/du + (1/4*u^4 + 1/4*v^4 + 1/2*(u^2 + 1)*v^2 + 1/2*u^2 + 1/4) d/dv*d/dv
We have::
sage: ginv.restrict(U) is g.restrict(U).inverse() True sage: ginv.restrict(V) is g.restrict(V).inverse() True sage: ginv.restrict(W) is g.restrict(W).inverse() True
The volume form (Levi-Civita tensor) associated with `g`::
sage: eps = g.volume_form() ; eps 2-form eps_g on the 2-dimensional differentiable manifold S^2 sage: eps.display(eU) eps_g = 4/(x^4 + y^4 + 2*(x^2 + 1)*y^2 + 2*x^2 + 1) dx/\dy sage: eps.display(eV) eps_g = 4/(u^4 + v^4 + 2*(u^2 + 1)*v^2 + 2*u^2 + 1) du/\dv
The unique non-trivial component of the volume form is nothing but the square root of the determinant of g in the corresponding frame::
sage: eps[[eU,1,2]] == g.sqrt_abs_det(eU) True sage: eps[[eV,1,2]] == g.sqrt_abs_det(eV) True
The Levi-Civita connection associated with the metric `g`::
sage: nabla = g.connection() ; nabla Levi-Civita connection nabla_g associated with the Riemannian metric g on the 2-dimensional differentiable manifold S^2 sage: latex(nabla) \nabla_{g}
The Christoffel symbols `\Gamma^i_{\ \, jk}` associated with some coordinates::
sage: g.christoffel_symbols(c_xy) 3-indices components w.r.t. Coordinate frame (U, (d/dx,d/dy)), with symmetry on the index positions (1, 2) sage: g.christoffel_symbols(c_xy)[:] [[[-2*x/(x^2 + y^2 + 1), -2*y/(x^2 + y^2 + 1)], [-2*y/(x^2 + y^2 + 1), 2*x/(x^2 + y^2 + 1)]], [[2*y/(x^2 + y^2 + 1), -2*x/(x^2 + y^2 + 1)], [-2*x/(x^2 + y^2 + 1), -2*y/(x^2 + y^2 + 1)]]] sage: g.christoffel_symbols(c_uv)[:] [[[-2*u/(u^2 + v^2 + 1), -2*v/(u^2 + v^2 + 1)], [-2*v/(u^2 + v^2 + 1), 2*u/(u^2 + v^2 + 1)]], [[2*v/(u^2 + v^2 + 1), -2*u/(u^2 + v^2 + 1)], [-2*u/(u^2 + v^2 + 1), -2*v/(u^2 + v^2 + 1)]]]
The Christoffel symbols are nothing but the connection coefficients w.r.t. the coordinate frame::
sage: g.christoffel_symbols(c_xy) is nabla.coef(c_xy.frame()) True sage: g.christoffel_symbols(c_uv) is nabla.coef(c_uv.frame()) True
Test that `\nabla` is the connection compatible with `g`::
sage: t = nabla(g) ; t Tensor field nabla_g(g) of type (0,3) on the 2-dimensional differentiable manifold S^2 sage: t.display(eU) nabla_g(g) = 0 sage: t.display(eV) nabla_g(g) = 0 sage: t == 0 True
The Riemann curvature tensor of `g`::
sage: riem = g.riemann() ; riem Tensor field Riem(g) of type (1,3) on the 2-dimensional differentiable manifold S^2 sage: riem.display(eU) Riem(g) = 4/(x^4 + y^4 + 2*(x^2 + 1)*y^2 + 2*x^2 + 1) d/dx*dy*dx*dy - 4/(x^4 + y^4 + 2*(x^2 + 1)*y^2 + 2*x^2 + 1) d/dx*dy*dy*dx - 4/(x^4 + y^4 + 2*(x^2 + 1)*y^2 + 2*x^2 + 1) d/dy*dx*dx*dy + 4/(x^4 + y^4 + 2*(x^2 + 1)*y^2 + 2*x^2 + 1) d/dy*dx*dy*dx sage: riem.display(eV) Riem(g) = 4/(u^4 + v^4 + 2*(u^2 + 1)*v^2 + 2*u^2 + 1) d/du*dv*du*dv - 4/(u^4 + v^4 + 2*(u^2 + 1)*v^2 + 2*u^2 + 1) d/du*dv*dv*du - 4/(u^4 + v^4 + 2*(u^2 + 1)*v^2 + 2*u^2 + 1) d/dv*du*du*dv + 4/(u^4 + v^4 + 2*(u^2 + 1)*v^2 + 2*u^2 + 1) d/dv*du*dv*du
The Ricci tensor of `g`::
sage: ric = g.ricci() ; ric Field of symmetric bilinear forms Ric(g) on the 2-dimensional differentiable manifold S^2 sage: ric.display(eU) Ric(g) = 4/(x^4 + y^4 + 2*(x^2 + 1)*y^2 + 2*x^2 + 1) dx*dx + 4/(x^4 + y^4 + 2*(x^2 + 1)*y^2 + 2*x^2 + 1) dy*dy sage: ric.display(eV) Ric(g) = 4/(u^4 + v^4 + 2*(u^2 + 1)*v^2 + 2*u^2 + 1) du*du + 4/(u^4 + v^4 + 2*(u^2 + 1)*v^2 + 2*u^2 + 1) dv*dv sage: ric == g True
The Ricci scalar of `g`::
sage: r = g.ricci_scalar() ; r Scalar field r(g) on the 2-dimensional differentiable manifold S^2 sage: r.display() r(g): S^2 --> R on U: (x, y) |--> 2 on V: (u, v) |--> 2
In dimension 2, the Riemann tensor can be expressed entirely in terms of the Ricci scalar `r`:
.. MATH::
R^i_{\ \, jlk} = \frac{r}{2} \left( \delta^i_{\ \, k} g_{jl} - \delta^i_{\ \, l} g_{jk} \right)
This formula can be checked here, with the r.h.s. rewritten as `-r g_{j[k} \delta^i_{\ \, l]}`::
sage: delta = M.tangent_identity_field() sage: riem == - r*(g*delta).antisymmetrize(2,3) True
""" '_schouten', '_cotton', '_cotton_york')
latex_name=None): r""" Construct a metric.
TESTS::
sage: M = Manifold(2, 'M') sage: U = M.open_subset('U') ; V = M.open_subset('V') sage: M.declare_union(U,V) # M is the union of U and V sage: c_xy.<x,y> = U.chart() ; c_uv.<u,v> = V.chart() sage: xy_to_uv = c_xy.transition_map(c_uv, (x+y, x-y), ....: intersection_name='W', restrictions1= x>0, ....: restrictions2= u+v>0) sage: uv_to_xy = xy_to_uv.inverse() sage: W = U.intersection(V) sage: e_xy = c_xy.frame() ; e_uv = c_uv.frame() sage: XM = M.vector_field_module() sage: from sage.manifolds.differentiable.metric import \ ....: PseudoRiemannianMetric sage: g = PseudoRiemannianMetric(XM, 'g', signature=0); g Lorentzian metric g on the 2-dimensional differentiable manifold M sage: g[e_xy,0,0], g[e_xy,1,1] = -(1+x^2), 1+y^2 sage: g.add_comp_by_continuation(e_uv, W, c_uv) sage: TestSuite(g).run(skip=['_test_category', '_test_pickling'])
.. TODO::
- fix _test_pickling (in the superclass TensorField) - add a specific parent to the metrics, to fit with the category framework
""" name=name, latex_name=latex_name, sym=(0,1)) # signature: else: raise TypeError("the metric signature must be an integer") raise ValueError("metric signature out of range") if ndim%2 == 0: raise ValueError("the metric signature must be even") else: raise ValueError("the metric signature must be odd") # the pair (n_+, n_-): # Initialization of derived quantities:
r""" String representation of the object.
TESTS::
sage: M = Manifold(5, 'M') sage: g = M.metric('g') sage: g._repr_() 'Riemannian metric g on the 5-dimensional differentiable manifold M' sage: g = M.metric('g', signature=3) sage: g._repr_() 'Lorentzian metric g on the 5-dimensional differentiable manifold M' sage: g = M.metric('g', signature=1) sage: g._repr_() 'Pseudo-Riemannian metric g on the 5-dimensional differentiable manifold M'
""" else:
r""" Create an instance of the same class as ``self`` with the same signature.
TESTS::
sage: M = Manifold(5, 'M') sage: g = M.metric('g', signature=3) sage: g1 = g._new_instance(); g1 Lorentzian metric unnamed metric on the 5-dimensional differentiable manifold M sage: type(g1) == type(g) True sage: g1.parent() is g.parent() True sage: g1.signature() == g.signature() True
""" signature=self._signature, latex_name=r'\mbox{unnamed metric}')
r""" Initialize the derived quantities.
TESTS::
sage: M = Manifold(5, 'M') sage: g = M.metric('g') sage: g._init_derived()
""" # Initialization of quantities pertaining to the mother class: # inverse metric: latex_name=inv_latex_name, sym=(0,1))
r""" Delete the derived quantities.
TESTS::
sage: M = Manifold(5, 'M') sage: g = M.metric('g') sage: g._del_derived()
""" # First the derived quantities from the mother class are deleted: # The inverse metric is cleared: # The connection, Ricci scalar and Weyl tensor are reset to None: # The Schouten, Cotton and Cotton-York tensors are reset to None: # The dictionary of determinants over the various frames is cleared: # The volume form and the associated tensors is deleted:
r""" Delete the inverse metric.
TESTS::
sage: M = Manifold(5, 'M') sage: g = M.metric('g') sage: g._del_inverse()
"""
r""" Signature of the metric.
OUTPUT:
- signature `S` of the metric, defined as the integer `S = n_+ - n_-`, where `n_+` (resp. `n_-`) is the number of positive terms (resp. number of negative terms) in any diagonal writing of the metric components
EXAMPLES:
Signatures on a 2-dimensional manifold::
sage: M = Manifold(2, 'M') sage: g = M.metric('g') # if not specified, the signature is Riemannian sage: g.signature() 2 sage: h = M.metric('h', signature=0) sage: h.signature() 0
"""
r""" Return the restriction of the metric to some subdomain.
If the restriction has not been defined yet, it is constructed here.
INPUT:
- ``subdomain`` -- open subset `U` of the metric's domain (must be an instance of :class:`~sage.manifolds.differentiable.manifold.DifferentiableManifold`) - ``dest_map`` -- (default: ``None``) destination map `\Phi:\ U \rightarrow V`, where `V` is a subdomain of ``self._codomain`` (type: :class:`~sage.manifolds.differentiable.diff_map.DiffMap`) If None, the restriction of ``self._vmodule._dest_map`` to `U` is used.
OUTPUT:
- instance of :class:`PseudoRiemannianMetric` representing the restriction.
EXAMPLES::
sage: M = Manifold(5, 'M') sage: g = M.metric('g', signature=3) sage: U = M.open_subset('U') sage: g.restrict(U) Lorentzian metric g on the Open subset U of the 5-dimensional differentiable manifold M sage: g.restrict(U).signature() 3
See the top documentation of :class:`PseudoRiemannianMetric` for more examples.
""" # Construct the restriction at the tensor field level: # the type is correctly handled by TensorField.restrict, i.e. # resu is of type self.__class__, but the signature is not handled # by TensorField.restrict; we have to set it here: # Restrictions of derived quantities: resu.__setattr__(attr, derived.restrict(subdomain)) for eps in self._vol_forms: resu._vol_forms.append(eps.restrict(subdomain)) # NB: no initialization of resu._determinants nor # resu._sqrt_abs_dets # The restriction is ready:
r""" Defines the metric from a field of symmetric bilinear forms
INPUT:
- ``symbiform`` -- instance of :class:`~sage.manifolds.differentiable.tensorfield.TensorField` representing a field of symmetric bilinear forms
EXAMPLES:
Metric defined from a field of symmetric bilinear forms on a non-parallelizable 2-dimensional manifold::
sage: M = Manifold(2, 'M') sage: U = M.open_subset('U') ; V = M.open_subset('V') sage: M.declare_union(U,V) # M is the union of U and V sage: c_xy.<x,y> = U.chart() ; c_uv.<u,v> = V.chart() sage: xy_to_uv = c_xy.transition_map(c_uv, (x+y, x-y), intersection_name='W', ....: restrictions1= x>0, restrictions2= u+v>0) sage: uv_to_xy = xy_to_uv.inverse() sage: W = U.intersection(V) sage: eU = c_xy.frame() ; eV = c_uv.frame() sage: h = M.sym_bilin_form_field(name='h') sage: h[eU,0,0], h[eU,0,1], h[eU,1,1] = 1+x, x*y, 1-y sage: h.add_comp_by_continuation(eV, W, c_uv) sage: h.display(eU) h = (x + 1) dx*dx + x*y dx*dy + x*y dy*dx + (-y + 1) dy*dy sage: h.display(eV) h = (1/8*u^2 - 1/8*v^2 + 1/4*v + 1/2) du*du + 1/4*u du*dv + 1/4*u dv*du + (-1/8*u^2 + 1/8*v^2 + 1/4*v + 1/2) dv*dv sage: g = M.metric('g') sage: g.set(h) sage: g.display(eU) g = (x + 1) dx*dx + x*y dx*dy + x*y dy*dx + (-y + 1) dy*dy sage: g.display(eV) g = (1/8*u^2 - 1/8*v^2 + 1/4*v + 1/2) du*du + 1/4*u du*dv + 1/4*u dv*du + (-1/8*u^2 + 1/8*v^2 + 1/4*v + 1/2) dv*dv
""" raise TypeError("the argument must be a tensor field") raise TypeError("the argument must be of tensor type (0,2)") raise TypeError("the argument must be symmetric") raise TypeError("the symmetric bilinear form is not defined " + "on the metric domain") rst = self.restrict(symbiform._domain) rst.set(symbiform) else:
r""" Return the inverse metric.
OUTPUT:
- instance of :class:`~sage.manifolds.differentiable.tensorfield.TensorField` with tensor_type = (2,0) representing the inverse metric
EXAMPLES:
Inverse of the standard metric on the 2-sphere::
sage: M = Manifold(2, 'S^2', start_index=1) sage: U = M.open_subset('U') ; V = M.open_subset('V') sage: M.declare_union(U,V) # S^2 is the union of U and V sage: c_xy.<x,y> = U.chart() ; c_uv.<u,v> = V.chart() # stereographic coord. sage: xy_to_uv = c_xy.transition_map(c_uv, (x/(x^2+y^2), y/(x^2+y^2)), ....: intersection_name='W', restrictions1= x^2+y^2!=0, ....: restrictions2= u^2+v^2!=0) sage: uv_to_xy = xy_to_uv.inverse() sage: W = U.intersection(V) # the complement of the two poles sage: eU = c_xy.frame() ; eV = c_uv.frame() sage: g = M.metric('g') sage: g[eU,1,1], g[eU,2,2] = 4/(1+x^2+y^2)^2, 4/(1+x^2+y^2)^2 sage: g.add_comp_by_continuation(eV, W, c_uv) sage: ginv = g.inverse(); ginv Tensor field inv_g of type (2,0) on the 2-dimensional differentiable manifold S^2 sage: ginv.display(eU) inv_g = (1/4*x^4 + 1/4*y^4 + 1/2*(x^2 + 1)*y^2 + 1/2*x^2 + 1/4) d/dx*d/dx + (1/4*x^4 + 1/4*y^4 + 1/2*(x^2 + 1)*y^2 + 1/2*x^2 + 1/4) d/dy*d/dy sage: ginv.display(eV) inv_g = (1/4*u^4 + 1/4*v^4 + 1/2*(u^2 + 1)*v^2 + 1/2*u^2 + 1/4) d/du*d/du + (1/4*u^4 + 1/4*v^4 + 1/2*(u^2 + 1)*v^2 + 1/2*u^2 + 1/4) d/dv*d/dv
Let us check that ``ginv`` is indeed the inverse of ``g``::
sage: s = g.contract(ginv); s # contraction of last index of g with first index of ginv Tensor field of type (1,1) on the 2-dimensional differentiable manifold S^2 sage: s == M.tangent_identity_field() True
""" # Is the inverse metric up to date ? # update of the restriction
r""" Return the unique torsion-free affine connection compatible with ``self``.
This is the so-called Levi-Civita connection.
INPUT:
- ``name`` -- (default: ``None``) name given to the Levi-Civita connection; if ``None``, it is formed from the metric name - ``latex_name`` -- (default: ``None``) LaTeX symbol to denote the Levi-Civita connection; if ``None``, it is set to ``name``, or if the latter is None as well, it formed from the symbol `\nabla` and the metric symbol
OUTPUT:
- the Levi-Civita connection, as an instance of :class:`~sage.manifolds.differentiable.levi_civita_connection.LeviCivitaConnection`.
EXAMPLES:
Levi-Civita connection associated with the Euclidean metric on `\RR^3`::
sage: M = Manifold(3, 'R^3', start_index=1) sage: # Let us use spherical coordinates on R^3: sage: U = M.open_subset('U') # the complement of the half-plane (y=0, x>=0) sage: c_spher.<r,th,ph> = U.chart(r'r:(0,+oo) th:(0,pi):\theta ph:(0,2*pi):\phi') sage: g = U.metric('g') sage: g[1,1], g[2,2], g[3,3] = 1, r^2 , (r*sin(th))^2 # the Euclidean metric sage: g.connection() Levi-Civita connection nabla_g associated with the Riemannian metric g on the Open subset U of the 3-dimensional differentiable manifold R^3 sage: g.connection().display() # Nonzero connection coefficients Gam^r_th,th = -r Gam^r_ph,ph = -r*sin(th)^2 Gam^th_r,th = 1/r Gam^th_th,r = 1/r Gam^th_ph,ph = -cos(th)*sin(th) Gam^ph_r,ph = 1/r Gam^ph_th,ph = cos(th)/sin(th) Gam^ph_ph,r = 1/r Gam^ph_ph,th = cos(th)/sin(th)
Test of compatibility with the metric::
sage: Dg = g.connection()(g) ; Dg Tensor field nabla_g(g) of type (0,3) on the Open subset U of the 3-dimensional differentiable manifold R^3 sage: Dg == 0 True sage: Dig = g.connection()(g.inverse()) ; Dig Tensor field nabla_g(inv_g) of type (2,1) on the Open subset U of the 3-dimensional differentiable manifold R^3 sage: Dig == 0 True
""" LeviCivitaConnection else: latex_name = name latex_name=latex_name)
r""" Christoffel symbols of ``self`` with respect to a chart.
INPUT:
- ``chart`` -- (default: ``None``) chart with respect to which the Christoffel symbols are required; if none is provided, the default chart of the metric's domain is assumed.
OUTPUT:
- the set of Christoffel symbols in the given chart, as an instance of :class:`~sage.tensor.modules.comp.CompWithSym`
EXAMPLES:
Christoffel symbols of the flat metric on `\RR^3` with respect to spherical coordinates::
sage: M = Manifold(3, 'R3', r'\RR^3', start_index=1) sage: U = M.open_subset('U') # the complement of the half-plane (y=0, x>=0) sage: X.<r,th,ph> = U.chart(r'r:(0,+oo) th:(0,pi):\theta ph:(0,2*pi):\phi') sage: g = U.metric('g') sage: g[1,1], g[2,2], g[3,3] = 1, r^2, r^2*sin(th)^2 sage: g.display() # the standard flat metric expressed in spherical coordinates g = dr*dr + r^2 dth*dth + r^2*sin(th)^2 dph*dph sage: Gam = g.christoffel_symbols() ; Gam 3-indices components w.r.t. Coordinate frame (U, (d/dr,d/dth,d/dph)), with symmetry on the index positions (1, 2) sage: type(Gam) <class 'sage.tensor.modules.comp.CompWithSym'> sage: Gam[:] [[[0, 0, 0], [0, -r, 0], [0, 0, -r*sin(th)^2]], [[0, 1/r, 0], [1/r, 0, 0], [0, 0, -cos(th)*sin(th)]], [[0, 0, 1/r], [0, 0, cos(th)/sin(th)], [1/r, cos(th)/sin(th), 0]]] sage: Gam[1,2,2] -r sage: Gam[2,1,2] 1/r sage: Gam[3,1,3] 1/r sage: Gam[3,2,3] cos(th)/sin(th) sage: Gam[2,3,3] -cos(th)*sin(th)
Note that a better display of the Christoffel symbols is provided by the method :meth:`christoffel_symbols_display`::
sage: g.christoffel_symbols_display() Gam^r_th,th = -r Gam^r_ph,ph = -r*sin(th)^2 Gam^th_r,th = 1/r Gam^th_ph,ph = -cos(th)*sin(th) Gam^ph_r,ph = 1/r Gam^ph_th,ph = cos(th)/sin(th)
""" else:
latex_symbol=None, index_labels=None, index_latex_labels=None, coordinate_labels=True, only_nonzero=True, only_nonredundant=True): r""" Display the Christoffel symbols w.r.t. to a given chart, one per line.
The output is either text-formatted (console mode) or LaTeX-formatted (notebook mode).
INPUT:
- ``chart`` -- (default: ``None``) chart with respect to which the Christoffel symbols are defined; if none is provided, the default chart of the metric's domain is assumed. - ``symbol`` -- (default: ``None``) string specifying the symbol of the connection coefficients; if ``None``, 'Gam' is used - ``latex_symbol`` -- (default: ``None``) string specifying the LaTeX symbol for the components; if ``None``, '\\Gamma' is used - ``index_labels`` -- (default: ``None``) list of strings representing the labels of each index; if ``None``, coordinate symbols are used except if ``coordinate_symbols`` is set to ``False``, in which case integer labels are used - ``index_latex_labels`` -- (default: ``None``) list of strings representing the LaTeX labels of each index; if ``None``, coordinate LaTeX symbols are used, except if ``coordinate_symbols`` is set to ``False``, in which case integer labels are used - ``coordinate_labels`` -- (default: ``True``) boolean; if ``True``, coordinate symbols are used by default (instead of integers) - ``only_nonzero`` -- (default: ``True``) boolean; if ``True``, only nonzero connection coefficients are displayed - ``only_nonredundant`` -- (default: ``True``) boolean; if ``True``, only nonredundant (w.r.t. the symmetry of the last two indices) connection coefficients are displayed
EXAMPLES:
Christoffel symbols of the flat metric on `\RR^3` with respect to spherical coordinates::
sage: M = Manifold(3, 'R3', r'\RR^3', start_index=1) sage: U = M.open_subset('U') # the complement of the half-plane (y=0, x>=0) sage: X.<r,th,ph> = U.chart(r'r:(0,+oo) th:(0,pi):\theta ph:(0,2*pi):\phi') sage: g = U.metric('g') sage: g[1,1], g[2,2], g[3,3] = 1, r^2, r^2*sin(th)^2 sage: g.display() # the standard flat metric expressed in spherical coordinates g = dr*dr + r^2 dth*dth + r^2*sin(th)^2 dph*dph sage: g.christoffel_symbols_display() Gam^r_th,th = -r Gam^r_ph,ph = -r*sin(th)^2 Gam^th_r,th = 1/r Gam^th_ph,ph = -cos(th)*sin(th) Gam^ph_r,ph = 1/r Gam^ph_th,ph = cos(th)/sin(th)
To list all nonzero Christoffel symbols, including those that can be deduced by symmetry, use ``only_nonredundant=False``::
sage: g.christoffel_symbols_display(only_nonredundant=False) Gam^r_th,th = -r Gam^r_ph,ph = -r*sin(th)^2 Gam^th_r,th = 1/r Gam^th_th,r = 1/r Gam^th_ph,ph = -cos(th)*sin(th) Gam^ph_r,ph = 1/r Gam^ph_th,ph = cos(th)/sin(th) Gam^ph_ph,r = 1/r Gam^ph_ph,th = cos(th)/sin(th)
Listing all Christoffel symbols (except those that can be deduced by symmetry), including the vanishing one::
sage: g.christoffel_symbols_display(only_nonzero=False) Gam^r_r,r = 0 Gam^r_r,th = 0 Gam^r_r,ph = 0 Gam^r_th,th = -r Gam^r_th,ph = 0 Gam^r_ph,ph = -r*sin(th)^2 Gam^th_r,r = 0 Gam^th_r,th = 1/r Gam^th_r,ph = 0 Gam^th_th,th = 0 Gam^th_th,ph = 0 Gam^th_ph,ph = -cos(th)*sin(th) Gam^ph_r,r = 0 Gam^ph_r,th = 0 Gam^ph_r,ph = 1/r Gam^ph_th,th = 0 Gam^ph_th,ph = cos(th)/sin(th) Gam^ph_ph,ph = 0
Using integer labels::
sage: g.christoffel_symbols_display(coordinate_labels=False) Gam^1_22 = -r Gam^1_33 = -r*sin(th)^2 Gam^2_12 = 1/r Gam^2_33 = -cos(th)*sin(th) Gam^3_13 = 1/r Gam^3_23 = cos(th)/sin(th)
""" symbol=symbol, latex_symbol=latex_symbol, index_labels=index_labels, index_latex_labels=index_latex_labels, coordinate_labels=coordinate_labels, only_nonzero=only_nonzero, only_nonredundant=only_nonredundant)
r""" Return the Riemann curvature tensor associated with the metric.
This method is actually a shortcut for ``self.connection().riemann()``
The Riemann curvature tensor is the tensor field `R` of type (1,3) defined by
.. MATH::
R(\omega, u, v, w) = \left\langle \omega, \nabla_u \nabla_v w - \nabla_v \nabla_u w - \nabla_{[u, v]} w \right\rangle
for any 1-form `\omega` and any vector fields `u`, `v` and `w`.
INPUT:
- ``name`` -- (default: ``None``) name given to the Riemann tensor; if none, it is set to "Riem(g)", where "g" is the metric's name - ``latex_name`` -- (default: ``None``) LaTeX symbol to denote the Riemann tensor; if none, it is set to "\\mathrm{Riem}(g)", where "g" is the metric's name
OUTPUT:
- the Riemann curvature tensor `R`, as an instance of :class:`~sage.manifolds.differentiable.tensorfield.TensorField`
EXAMPLES:
Riemann tensor of the standard metric on the 2-sphere::
sage: M = Manifold(2, 'S^2', start_index=1) sage: U = M.open_subset('U') # the complement of a meridian (domain of spherical coordinates) sage: c_spher.<th,ph> = U.chart(r'th:(0,pi):\theta ph:(0,2*pi):\phi') sage: a = var('a') # the sphere radius sage: g = U.metric('g') sage: g[1,1], g[2,2] = a^2, a^2*sin(th)^2 sage: g.display() # standard metric on the 2-sphere of radius a: g = a^2 dth*dth + a^2*sin(th)^2 dph*dph sage: g.riemann() Tensor field Riem(g) of type (1,3) on the Open subset U of the 2-dimensional differentiable manifold S^2 sage: g.riemann()[:] [[[[0, 0], [0, 0]], [[0, sin(th)^2], [-sin(th)^2, 0]]], [[[0, (cos(th)^2 - 1)/sin(th)^2], [1, 0]], [[0, 0], [0, 0]]]]
In dimension 2, the Riemann tensor can be expressed entirely in terms of the Ricci scalar `r`:
.. MATH::
R^i_{\ \, jlk} = \frac{r}{2} \left( \delta^i_{\ \, k} g_{jl} - \delta^i_{\ \, l} g_{jk} \right)
This formula can be checked here, with the r.h.s. rewritten as `-r g_{j[k} \delta^i_{\ \, l]}`::
sage: g.riemann() == \ ....: -g.ricci_scalar()*(g*U.tangent_identity_field()).antisymmetrize(2,3) True
Using SymPy as symbolic engine::
sage: M.set_calculus_method('sympy') sage: g = U.metric('g') sage: g[1,1], g[2,2] = a**2, a**2*sin(th)**2 sage: g.riemann()[:] [[[[0, 0], [0, 0]], [[0, sin(2*th)/(2*tan(th)) - cos(2*th)], [-sin(2*th)/(2*tan(th)) + cos(2*th), 0]]], [[[0, -1], [1, 0]], [[0, 0], [0, 0]]]]
"""
r""" Return the Ricci tensor associated with the metric.
This method is actually a shortcut for ``self.connection().ricci()``
The Ricci tensor is the tensor field `Ric` of type (0,2) defined from the Riemann curvature tensor `R` by
.. MATH::
Ric(u, v) = R(e^i, u, e_i, v)
for any vector fields `u` and `v`, `(e_i)` being any vector frame and `(e^i)` the dual coframe.
INPUT:
- ``name`` -- (default: ``None``) name given to the Ricci tensor; if none, it is set to "Ric(g)", where "g" is the metric's name - ``latex_name`` -- (default: ``None``) LaTeX symbol to denote the Ricci tensor; if none, it is set to "\\mathrm{Ric}(g)", where "g" is the metric's name
OUTPUT:
- the Ricci tensor `Ric`, as an instance of :class:`~sage.manifolds.differentiable.tensorfield.TensorField` of tensor type (0,2) and symmetric
EXAMPLES:
Ricci tensor of the standard metric on the 2-sphere::
sage: M = Manifold(2, 'S^2', start_index=1) sage: U = M.open_subset('U') # the complement of a meridian (domain of spherical coordinates) sage: c_spher.<th,ph> = U.chart(r'th:(0,pi):\theta ph:(0,2*pi):\phi') sage: a = var('a') # the sphere radius sage: g = U.metric('g') sage: g[1,1], g[2,2] = a^2, a^2*sin(th)^2 sage: g.display() # standard metric on the 2-sphere of radius a: g = a^2 dth*dth + a^2*sin(th)^2 dph*dph sage: g.ricci() Field of symmetric bilinear forms Ric(g) on the Open subset U of the 2-dimensional differentiable manifold S^2 sage: g.ricci()[:] [ 1 0] [ 0 sin(th)^2] sage: g.ricci() == a^(-2) * g True
"""
r""" Return the Ricci scalar associated with the metric.
The Ricci scalar is the scalar field `r` defined from the Ricci tensor `Ric` and the metric tensor `g` by
.. MATH::
r = g^{ij} Ric_{ij}
INPUT:
- ``name`` -- (default: ``None``) name given to the Ricci scalar; if none, it is set to "r(g)", where "g" is the metric's name - ``latex_name`` -- (default: ``None``) LaTeX symbol to denote the Ricci scalar; if none, it is set to "\\mathrm{r}(g)", where "g" is the metric's name
OUTPUT:
- the Ricci scalar `r`, as an instance of :class:`~sage.manifolds.differentiable.scalarfield.DiffScalarField`
EXAMPLES:
Ricci scalar of the standard metric on the 2-sphere::
sage: M = Manifold(2, 'S^2', start_index=1) sage: U = M.open_subset('U') # the complement of a meridian (domain of spherical coordinates) sage: c_spher.<th,ph> = U.chart(r'th:(0,pi):\theta ph:(0,2*pi):\phi') sage: a = var('a') # the sphere radius sage: g = U.metric('g') sage: g[1,1], g[2,2] = a^2, a^2*sin(th)^2 sage: g.display() # standard metric on the 2-sphere of radius a: g = a^2 dth*dth + a^2*sin(th)^2 dph*dph sage: g.ricci_scalar() Scalar field r(g) on the Open subset U of the 2-dimensional differentiable manifold S^2 sage: g.ricci_scalar().display() # The Ricci scalar is constant: r(g): U --> R (th, ph) |--> 2/a^2
""" r"\right)"
r""" Return the Weyl conformal tensor associated with the metric.
The Weyl conformal tensor is the tensor field `C` of type (1,3) defined as the trace-free part of the Riemann curvature tensor `R`
INPUT:
- ``name`` -- (default: ``None``) name given to the Weyl conformal tensor; if ``None``, it is set to "C(g)", where "g" is the metric's name - ``latex_name`` -- (default: ``None``) LaTeX symbol to denote the Weyl conformal tensor; if ``None``, it is set to "\\mathrm{C}(g)", where "g" is the metric's name
OUTPUT:
- the Weyl conformal tensor `C`, as an instance of :class:`~sage.manifolds.differentiable.tensorfield.TensorField`
EXAMPLES:
Checking that the Weyl tensor identically vanishes on a 3-dimensional manifold, for instance the hyperbolic space `H^3`::
sage: M = Manifold(3, 'H^3', start_index=1) sage: U = M.open_subset('U') # the complement of the half-plane (y=0, x>=0) sage: X.<rh,th,ph> = U.chart(r'rh:(0,+oo):\rho th:(0,pi):\theta ph:(0,2*pi):\phi') sage: g = U.metric('g') sage: b = var('b') sage: g[1,1], g[2,2], g[3,3] = b^2, (b*sinh(rh))^2, (b*sinh(rh)*sin(th))^2 sage: g.display() # standard metric on H^3: g = b^2 drh*drh + b^2*sinh(rh)^2 dth*dth + b^2*sin(th)^2*sinh(rh)^2 dph*dph sage: C = g.weyl() ; C Tensor field C(g) of type (1,3) on the Open subset U of the 3-dimensional differentiable manifold H^3 sage: C == 0 True
""" raise ValueError("the Weyl tensor is not defined for a " + "manifold of dimension n <= 2") self._vmodule._dest_map) # First index of the Ricci tensor raised with the metric
r""" Return the Schouten tensor associated with the metric.
The Schouten tensor is the tensor field `Sc` of type (0,2) defined from the Ricci curvature tensor `Ric` (see :meth:`ricci`) and the scalar curvature `r` (see :meth:`ricci_scalar`) and the metric `g` by
.. MATH::
Sc(u, v) = \frac{1}{n-2}\left(Ric(u, v) + \frac{r}{2(n-1)}g(u,v) \right)
for any vector fields `u` and `v`.
INPUT:
- ``name`` -- (default: ``None``) name given to the Schouten tensor; if none, it is set to "Schouten(g)", where "g" is the metric's name - ``latex_name`` -- (default: ``None``) LaTeX symbol to denote the Schouten tensor; if none, it is set to "\\mathrm{Schouten}(g)", where "g" is the metric's name
OUTPUT:
- the Schouten tensor `Sc`, as an instance of :class:`~sage.manifolds.differentiable.tensorfield.TensorField` of tensor type (0,2) and symmetric
EXAMPLES:
Schouten tensor of the left invariant metric of Heisenberg's Nil group::
sage: M = Manifold(3, 'Nil', start_index=1) sage: X.<x,y,z> = M.chart() sage: g = M.riemannian_metric('g') sage: g[1,1], g[2,2], g[2,3], g[3,3] = 1, 1+x^2, -x, 1 sage: g.display() g = dx*dx + (x^2 + 1) dy*dy - x dy*dz - x dz*dy + dz*dz sage: g.schouten() Field of symmetric bilinear forms Schouten(g) on the 3-dimensional differentiable manifold Nil sage: g.schouten().display() Schouten(g) = -3/8 dx*dx + (5/8*x^2 - 3/8) dy*dy - 5/8*x dy*dz - 5/8*x dz*dy + 5/8 dz*dz
""" raise ValueError("the Schouten tensor is only defined for a " + "manifold of dimension >= 3")
r""" Return the Cotton conformal tensor associated with the metric. The tensor has type (0,3) and is defined in terms of the Schouten tensor `S` (see :meth:`schouten`):
.. MATH::
C_{ijk} = (n-2) \left(\nabla_k S_{ij} - \nabla_j S_{ik}\right)
INPUT:
- ``name`` -- (default: ``None``) name given to the Cotton conformal tensor; if ``None``, it is set to "Cot(g)", where "g" is the metric's name - ``latex_name`` -- (default: ``None``) LaTeX symbol to denote the Cotton conformal tensor; if ``None``, it is set to "\\mathrm{Cot}(g)", where "g" is the metric's name
OUTPUT:
- the Cotton conformal tensor `Cot`, as an instance of :class:`~sage.manifolds.differentiable.tensorfield.TensorField`
EXAMPLES:
Checking that the Cotton tensor identically vanishes on a conformally flat 3-dimensional manifold, for instance the hyperbolic space `H^3`::
sage: M = Manifold(3, 'H^3', start_index=1) sage: U = M.open_subset('U') # the complement of the half-plane (y=0, x>=0) sage: X.<rh,th,ph> = U.chart(r'rh:(0,+oo):\rho th:(0,pi):\theta ph:(0,2*pi):\phi') sage: g = U.metric('g') sage: b = var('b') sage: g[1,1], g[2,2], g[3,3] = b^2, (b*sinh(rh))^2, (b*sinh(rh)*sin(th))^2 sage: g.display() # standard metric on H^3: g = b^2 drh*drh + b^2*sinh(rh)^2 dth*dth + b^2*sin(th)^2*sinh(rh)^2 dph*dph sage: Cot = g.cotton() ; Cot # long time Tensor field Cot(g) of type (0,3) on the Open subset U of the 3-dimensional differentiable manifold H^3 sage: Cot == 0 # long time True
""" n = self._ambient_domain.dimension() if n < 3: raise ValueError("the Cotton tensor is only defined for a " + "manifold of dimension >= 3") if self._cotton is None: nabla = self.connection() s = self.schouten() cot = 2*(n-2)*nabla(s).antisymmetrize(1,2) name = name or 'Cot(' + self._name + ')' latex_name = latex_name or r'\mathrm{Cot}(' + self._latex_name + ')' cot.set_name(name=name, latex_name=latex_name) self._cotton = cot return self._cotton
r""" Return the Cotton-York conformal tensor associated with the metric. The tensor has type (0,2) and is only defined for manifolds of dimension 3. It is defined in terms of the Cotton tensor `C` (see :meth:`cotton`) or the Schouten tensor `S` (see :meth:`schouten`):
.. MATH::
CY_{ij} = \frac{1}{2} \epsilon^{kl}_{\ \ \, i} C_{jlk} = \epsilon^{kl}_{\ \ \, i} \nabla_k S_{lj}
INPUT:
- ``name`` -- (default: ``None``) name given to the Cotton-York tensor; if ``None``, it is set to "CY(g)", where "g" is the metric's name - ``latex_name`` -- (default: ``None``) LaTeX symbol to denote the Cotton-York tensor; if ``None``, it is set to "\\mathrm{CY}(g)", where "g" is the metric's name
OUTPUT:
- the Cotton-York conformal tensor `CY`, as an instance of :class:`~sage.manifolds.differentiable.tensorfield.TensorField`
EXAMPLES:
Compute the determinant of the Cotton-York tensor for the Heisenberg group with the left invariant metric::
sage: M = Manifold(3, 'Nil', start_index=1) sage: X.<x,y,z> = M.chart() sage: g = M.riemannian_metric('g') sage: g[1,1], g[2,2], g[2,3], g[3,3] = 1, 1+x^2, -x, 1 sage: g.display() g = dx*dx + (x^2 + 1) dy*dy - x dy*dz - x dz*dy + dz*dz sage: CY = g.cotton_york() ; CY # long time Tensor field CY(g) of type (0,2) on the 3-dimensional differentiable manifold Nil sage: CY.display() # long time CY(g) = 1/2 dx*dx + (-x^2 + 1/2) dy*dy + x dy*dz + x dz*dy - dz*dz sage: det(CY[:]) # long time -1/4
""" n = self._ambient_domain.dimension() if n != 3: raise ValueError("the Cotton-York tensor is only defined for a " + "manifold of dimension 3") if self._cotton_york is None: cot = self.cotton() eps = self.volume_form(2) cy = eps.contract(0, 1, cot, 2, 1)/2 name = name or 'CY(' + self._name + ')' latex_name = latex_name or r'\mathrm{CY}(' + self._latex_name + ')' cy.set_name(name=name, latex_name=latex_name) self._cotton_york = cy return self._cotton_york
r""" Determinant of the metric components in the specified frame.
INPUT:
- ``frame`` -- (default: ``None``) vector frame with respect to which the components `g_{ij}` of the metric are defined; if ``None``, the default frame of the metric's domain is used. If a chart is provided instead of a frame, the associated coordinate frame is used
OUTPUT:
- the determinant `\det (g_{ij})`, as an instance of :class:`~sage.manifolds.differentiable.scalarfield.DiffScalarField`
EXAMPLES:
Metric determinant on a 2-dimensional manifold::
sage: M = Manifold(2, 'M', start_index=1) sage: X.<x,y> = M.chart() sage: g = M.metric('g') sage: g[1,1], g[1, 2], g[2, 2] = 1+x, x*y , 1-y sage: g[:] [ x + 1 x*y] [ x*y -y + 1] sage: s = g.determinant() # determinant in M's default frame sage: s.expr() -x^2*y^2 - (x + 1)*y + x + 1
A shortcut is ``det()``::
sage: g.det() == g.determinant() True
The notation ``det(g)`` can be used::
sage: det(g) == g.determinant() True
Determinant in a frame different from the default's one::
sage: Y.<u,v> = M.chart() sage: ch_X_Y = X.transition_map(Y, [x+y, x-y]) sage: ch_X_Y.inverse() Change of coordinates from Chart (M, (u, v)) to Chart (M, (x, y)) sage: g.comp(Y.frame())[:, Y] [ 1/8*u^2 - 1/8*v^2 + 1/4*v + 1/2 1/4*u] [ 1/4*u -1/8*u^2 + 1/8*v^2 + 1/4*v + 1/2] sage: g.determinant(Y.frame()).expr() -1/4*x^2*y^2 - 1/4*(x + 1)*y + 1/4*x + 1/4 sage: g.determinant(Y.frame()).expr(Y) -1/64*u^4 - 1/64*v^4 + 1/32*(u^2 + 2)*v^2 - 1/16*u^2 + 1/4*v + 1/4
A chart can be passed instead of a frame::
sage: g.determinant(X) is g.determinant(X.frame()) True sage: g.determinant(Y) is g.determinant(Y.frame()) True
The metric determinant depends on the frame::
sage: g.determinant(X.frame()) == g.determinant(Y.frame()) False
Using SymPy as symbolic engine::
sage: M.set_calculus_method('sympy') sage: g = M.metric('g') sage: g[1,1], g[1, 2], g[2, 2] = 1+x, x*y , 1-y sage: s = g.determinant() # determinant in M's default frame sage: s.expr() -x**2*y**2 + x - y*(x + 1) + 1
""" # frame is actually a chart and is changed to the associated # coordinate frame: # a new computation is necessary # TODO: do the computation without the 'SR' enforcement for j in manif.irange()] for i in manif.irange()] )
r""" Square root of the absolute value of the determinant of the metric components in the specified frame.
INPUT:
- ``frame`` -- (default: ``None``) vector frame with respect to which the components `g_{ij}` of ``self`` are defined; if ``None``, the domain's default frame is used. If a chart is provided, the associated coordinate frame is used
OUTPUT:
- `\sqrt{|\det (g_{ij})|}`, as an instance of :class:`~sage.manifolds.differentiable.scalarfield.DiffScalarField`
EXAMPLES:
Standard metric in the Euclidean space `\RR^3` with spherical coordinates::
sage: M = Manifold(3, 'M', start_index=1) sage: U = M.open_subset('U') # the complement of the half-plane (y=0, x>=0) sage: c_spher.<r,th,ph> = U.chart(r'r:(0,+oo) th:(0,pi):\theta ph:(0,2*pi):\phi') sage: g = U.metric('g') sage: g[1,1], g[2,2], g[3,3] = 1, r^2, (r*sin(th))^2 sage: g.display() g = dr*dr + r^2 dth*dth + r^2*sin(th)^2 dph*dph sage: g.sqrt_abs_det().expr() r^2*sin(th)
Metric determinant on a 2-dimensional manifold::
sage: M = Manifold(2, 'M', start_index=1) sage: X.<x,y> = M.chart() sage: g = M.metric('g') sage: g[1,1], g[1, 2], g[2, 2] = 1+x, x*y , 1-y sage: g[:] [ x + 1 x*y] [ x*y -y + 1] sage: s = g.sqrt_abs_det() ; s Scalar field on the 2-dimensional differentiable manifold M sage: s.expr() sqrt(-x^2*y^2 - (x + 1)*y + x + 1)
Determinant in a frame different from the default's one::
sage: Y.<u,v> = M.chart() sage: ch_X_Y = X.transition_map(Y, [x+y, x-y]) sage: ch_X_Y.inverse() Change of coordinates from Chart (M, (u, v)) to Chart (M, (x, y)) sage: g[Y.frame(),:,Y] [ 1/8*u^2 - 1/8*v^2 + 1/4*v + 1/2 1/4*u] [ 1/4*u -1/8*u^2 + 1/8*v^2 + 1/4*v + 1/2] sage: g.sqrt_abs_det(Y.frame()).expr() 1/2*sqrt(-x^2*y^2 - (x + 1)*y + x + 1) sage: g.sqrt_abs_det(Y.frame()).expr(Y) 1/8*sqrt(-u^4 - v^4 + 2*(u^2 + 2)*v^2 - 4*u^2 + 16*v + 16)
A chart can be passed instead of a frame::
sage: g.sqrt_abs_det(Y) is g.sqrt_abs_det(Y.frame()) True
The metric determinant depends on the frame::
sage: g.sqrt_abs_det(X.frame()) == g.sqrt_abs_det(Y.frame()) False
Using SymPy as symbolic engine::
sage: M.set_calculus_method('sympy') sage: g = M.metric('g') sage: g[1,1], g[1, 2], g[2, 2] = 1+x, x*y , 1-y sage: g.sqrt_abs_det().expr() sqrt(-x**2*y**2 - x*y + x - y + 1) sage: g.sqrt_abs_det(Y.frame()).expr() sqrt(-x**2*y**2 - x*y + x - y + 1)/2 sage: g.sqrt_abs_det(Y.frame()).expr(Y) sqrt(-u**4 + 2*u**2*v**2 - 4*u**2 - v**4 + 4*v**2 + 16*v + 16)/8
""" # frame is actually a chart and is changed to the associated # coordinate frame: # a new computation is necessary
r""" Volume form (Levi-Civita tensor) `\epsilon` associated with the metric.
This assumes that the manifold is orientable.
The volume form `\epsilon` is a `n`-form (`n` being the manifold's dimension) such that for any vector basis `(e_i)` that is orthonormal with respect to the metric,
.. MATH::
\epsilon(e_1,\ldots,e_n) = \pm 1
There are only two such `n`-forms, which are opposite of each other. The volume form `\epsilon` is selected such that the domain's default frame is right-handed with respect to it.
INPUT:
- ``contra`` -- (default: 0) number of contravariant indices of the returned tensor
OUTPUT:
- if ``contra = 0`` (default value): the volume `n`-form `\epsilon`, as an instance of :class:`~sage.manifolds.differentiable.diff_form.DiffForm` - if ``contra = k``, with `1\leq k \leq n`, the tensor field of type (k,n-k) formed from `\epsilon` by raising the first k indices with the metric (see method :meth:`~sage.manifolds.differentiable.tensorfield.TensorField.up`); the output is then an instance of :class:`~sage.manifolds.differentiable.tensorfield.TensorField`, with the appropriate antisymmetries, or of the subclass :class:`~sage.manifolds.differentiable.multivectorfield.MultivectorField` if `k=n`
EXAMPLES:
Volume form on `\RR^3` with spherical coordinates::
sage: M = Manifold(3, 'M', start_index=1) sage: U = M.open_subset('U') # the complement of the half-plane (y=0, x>=0) sage: c_spher.<r,th,ph> = U.chart(r'r:(0,+oo) th:(0,pi):\theta ph:(0,2*pi):\phi') sage: g = U.metric('g') sage: g[1,1], g[2,2], g[3,3] = 1, r^2, (r*sin(th))^2 sage: g.display() g = dr*dr + r^2 dth*dth + r^2*sin(th)^2 dph*dph sage: eps = g.volume_form() ; eps 3-form eps_g on the Open subset U of the 3-dimensional differentiable manifold M sage: eps.display() eps_g = r^2*sin(th) dr/\dth/\dph sage: eps[[1,2,3]] == g.sqrt_abs_det() True sage: latex(eps) \epsilon_{g}
The tensor field of components `\epsilon^i_{\ \, jk}` (``contra=1``)::
sage: eps1 = g.volume_form(1) ; eps1 Tensor field of type (1,2) on the Open subset U of the 3-dimensional differentiable manifold M sage: eps1.symmetries() no symmetry; antisymmetry: (1, 2) sage: eps1[:] [[[0, 0, 0], [0, 0, r^2*sin(th)], [0, -r^2*sin(th), 0]], [[0, 0, -sin(th)], [0, 0, 0], [sin(th), 0, 0]], [[0, 1/sin(th), 0], [-1/sin(th), 0, 0], [0, 0, 0]]]
The tensor field of components `\epsilon^{ij}_{\ \ k}` (``contra=2``)::
sage: eps2 = g.volume_form(2) ; eps2 Tensor field of type (2,1) on the Open subset U of the 3-dimensional differentiable manifold M sage: eps2.symmetries() no symmetry; antisymmetry: (0, 1) sage: eps2[:] [[[0, 0, 0], [0, 0, sin(th)], [0, -1/sin(th), 0]], [[0, 0, -sin(th)], [0, 0, 0], [1/(r^2*sin(th)), 0, 0]], [[0, 1/sin(th), 0], [-1/(r^2*sin(th)), 0, 0], [0, 0, 0]]]
The tensor field of components `\epsilon^{ijk}` (``contra=3``)::
sage: eps3 = g.volume_form(3) ; eps3 3-vector field on the Open subset U of the 3-dimensional differentiable manifold M sage: eps3.tensor_type() (3, 0) sage: eps3.symmetries() no symmetry; antisymmetry: (0, 1, 2) sage: eps3[:] [[[0, 0, 0], [0, 0, 1/(r^2*sin(th))], [0, -1/(r^2*sin(th)), 0]], [[0, 0, -1/(r^2*sin(th))], [0, 0, 0], [1/(r^2*sin(th)), 0, 0]], [[0, 1/(r^2*sin(th)), 0], [-1/(r^2*sin(th)), 0, 0], [0, 0, 0]]] sage: eps3[1,2,3] 1/(r^2*sin(th)) sage: eps3[[1,2,3]] * g.sqrt_abs_det() == 1 True
""" # a new computation is necessary # The result is constructed on the vector field module, # so that dest_map is taken automatically into account: latex_name=r'\epsilon_{'+self._latex_name+r'}') # Tensors related to the Levi-Civita one by index rising: # restoring the antisymmetry after the up operation:
r""" Compute the Hodge dual of a differential form with respect to the metric.
If the differential form is a `p`-form `A`, its *Hodge dual* with respect to the metric `g` is the `(n-p)`-form `*A` defined by
.. MATH::
*A_{i_1\ldots i_{n-p}} = \frac{1}{p!} A_{k_1\ldots k_p} \epsilon^{k_1\ldots k_p}_{\qquad\ i_1\ldots i_{n-p}}
where `n` is the manifold's dimension, `\epsilon` is the volume `n`-form associated with `g` (see :meth:`volume_form`) and the indices `k_1,\ldots, k_p` are raised with `g`.
INPUT:
- ``pform``: a `p`-form `A`; must be an instance of :class:`~sage.manifolds.differentiable.scalarfield.DiffScalarField` for `p=0` and of :class:`~sage.manifolds.differentiable.diff_form.DiffForm` or :class:`~sage.manifolds.differentiable.diff_form.DiffFormParal` for `p\geq 1`.
OUTPUT:
- the `(n-p)`-form `*A`
EXAMPLES:
Hodge dual of a 1-form in the Euclidean space `R^3`::
sage: M = Manifold(3, 'M', start_index=1) sage: X.<x,y,z> = M.chart() sage: g = M.metric('g') sage: g[1,1], g[2,2], g[3,3] = 1, 1, 1 sage: a = M.one_form('A') sage: var('Ax Ay Az') (Ax, Ay, Az) sage: a[:] = (Ax, Ay, Az) sage: sa = g.hodge_star(a) ; sa 2-form *A on the 3-dimensional differentiable manifold M sage: sa.display() *A = Az dx/\dy - Ay dx/\dz + Ax dy/\dz sage: ssa = g.hodge_star(sa) ; ssa 1-form **A on the 3-dimensional differentiable manifold M sage: ssa.display() **A = Ax dx + Ay dy + Az dz sage: ssa == a # must hold for a Riemannian metric in dimension 3 True
Hodge dual of a 0-form (scalar field) in `R^3`::
sage: f = M.scalar_field(function('F')(x,y,z), name='f') sage: sf = g.hodge_star(f) ; sf 3-form *f on the 3-dimensional differentiable manifold M sage: sf.display() *f = F(x, y, z) dx/\dy/\dz sage: ssf = g.hodge_star(sf) ; ssf Scalar field **f on the 3-dimensional differentiable manifold M sage: ssf.display() **f: M --> R (x, y, z) |--> F(x, y, z) sage: ssf == f # must hold for a Riemannian metric True
Hodge dual of a 0-form in Minkowksi spacetime::
sage: M = Manifold(4, 'M') sage: X.<t,x,y,z> = M.chart() sage: g = M.lorentzian_metric('g') sage: g[0,0], g[1,1], g[2,2], g[3,3] = -1, 1, 1, 1 sage: g.display() # Minkowski metric g = -dt*dt + dx*dx + dy*dy + dz*dz sage: var('f0') f0 sage: f = M.scalar_field(f0, name='f') sage: sf = g.hodge_star(f) ; sf 4-form *f on the 4-dimensional differentiable manifold M sage: sf.display() *f = f0 dt/\dx/\dy/\dz sage: ssf = g.hodge_star(sf) ; ssf Scalar field **f on the 4-dimensional differentiable manifold M sage: ssf.display() **f: M --> R (t, x, y, z) |--> -f0 sage: ssf == -f # must hold for a Lorentzian metric True
Hodge dual of a 1-form in Minkowksi spacetime::
sage: a = M.one_form('A') sage: var('At Ax Ay Az') (At, Ax, Ay, Az) sage: a[:] = (At, Ax, Ay, Az) sage: a.display() A = At dt + Ax dx + Ay dy + Az dz sage: sa = g.hodge_star(a) ; sa 3-form *A on the 4-dimensional differentiable manifold M sage: sa.display() *A = -Az dt/\dx/\dy + Ay dt/\dx/\dz - Ax dt/\dy/\dz - At dx/\dy/\dz sage: ssa = g.hodge_star(sa) ; ssa 1-form **A on the 4-dimensional differentiable manifold M sage: ssa.display() **A = At dt + Ax dx + Ay dy + Az dz sage: ssa == a # must hold for a Lorentzian metric in dimension 4 True
Hodge dual of a 2-form in Minkowksi spacetime::
sage: F = M.diff_form(2, 'F') sage: var('Ex Ey Ez Bx By Bz') (Ex, Ey, Ez, Bx, By, Bz) sage: F[0,1], F[0,2], F[0,3] = -Ex, -Ey, -Ez sage: F[1,2], F[1,3], F[2,3] = Bz, -By, Bx sage: F[:] [ 0 -Ex -Ey -Ez] [ Ex 0 Bz -By] [ Ey -Bz 0 Bx] [ Ez By -Bx 0] sage: sF = g.hodge_star(F) ; sF 2-form *F on the 4-dimensional differentiable manifold M sage: sF[:] [ 0 Bx By Bz] [-Bx 0 Ez -Ey] [-By -Ez 0 Ex] [-Bz Ey -Ex 0] sage: ssF = g.hodge_star(sF) ; ssF 2-form **F on the 4-dimensional differentiable manifold M sage: ssF[:] [ 0 Ex Ey Ez] [-Ex 0 -Bz By] [-Ey Bz 0 -Bx] [-Ez -By Bx 0] sage: ssF.display() **F = Ex dt/\dx + Ey dt/\dy + Ez dt/\dz - Bz dx/\dy + By dx/\dz - Bx dy/\dz sage: F.display() F = -Ex dt/\dx - Ey dt/\dy - Ez dt/\dz + Bz dx/\dy - By dx/\dz + Bx dy/\dz sage: ssF == -F # must hold for a Lorentzian metric in dimension 4 True
Test of the standard identity
.. MATH::
*(A\wedge B) = \epsilon(A^\sharp, B^\sharp, ., .)
where `A` and `B` are any 1-forms and `A^\sharp` and `B^\sharp` the vectors associated to them by the metric `g` (index raising)::
sage: b = M.one_form('B') sage: var('Bt Bx By Bz') (Bt, Bx, By, Bz) sage: b[:] = (Bt, Bx, By, Bz) ; b.display() B = Bt dt + Bx dx + By dy + Bz dz sage: epsilon = g.volume_form() sage: g.hodge_star(a.wedge(b)) == epsilon.contract(0,a.up(g)).contract(0,b.up(g)) True
""" format_unop_latex else: latex_name=format_unop_latex(r'\star ', pform._latex_name))
#******************************************************************************
r""" Pseudo-Riemannian metric with values on a parallelizable manifold.
An instance of this class is a field of nondegenerate symmetric bilinear forms (metric field) along a differentiable manifold `U` with values in a parallelizable manifold `M` over `\RR`, via a differentiable mapping `\Phi: U \rightarrow M`. The standard case of a metric field *on* a manifold corresponds to `U=M` and `\Phi = \mathrm{Id}_M`. Other common cases are `\Phi` being an immersion and `\Phi` being a curve in `M` (`U` is then an open interval of `\RR`).
A *metric* `g` is a field on `U`, such that at each point `p\in U`, `g(p)` is a bilinear map of the type:
.. MATH::
g(p):\ T_q M\times T_q M \longrightarrow \RR
where `T_q M` stands for the tangent space to manifold `M` at the point `q=\Phi(p)`, such that `g(p)` is symmetric: `\forall (u,v)\in T_q M\times T_q M, \ g(p)(v,u) = g(p)(u,v)` and nondegenerate: `(\forall v\in T_q M,\ \ g(p)(u,v) = 0) \Longrightarrow u=0`.
.. NOTE::
If `M` is not parallelizable, the class :class:`PseudoRiemannianMetric` should be used instead.
INPUT:
- ``vector_field_module`` -- free module `\mathfrak{X}(U,\Phi)` of vector fields along `U` with values on `\Phi(U)\subset M` - ``name`` -- name given to the metric - ``signature`` -- (default: ``None``) signature `S` of the metric as a single integer: `S = n_+ - n_-`, where `n_+` (resp. `n_-`) is the number of positive terms (resp. number of negative terms) in any diagonal writing of the metric components; if ``signature`` is ``None``, `S` is set to the dimension of manifold `M` (Riemannian signature) - ``latex_name`` -- (default: ``None``) LaTeX symbol to denote the metric; if ``None``, it is formed from ``name``
EXAMPLES:
Metric on a 2-dimensional manifold::
sage: M = Manifold(2, 'M', start_index=1) sage: c_xy.<x,y> = M.chart() sage: g = M.metric('g') ; g Riemannian metric g on the 2-dimensional differentiable manifold M sage: latex(g) g
A metric is a special kind of tensor field and therefore inheritates all the properties from class :class:`~sage.manifolds.differentiable.tensorfield.TensorField`::
sage: g.parent() Free module T^(0,2)(M) of type-(0,2) tensors fields on the 2-dimensional differentiable manifold M sage: g.tensor_type() (0, 2) sage: g.symmetries() # g is symmetric: symmetry: (0, 1); no antisymmetry
Setting the metric components in the manifold's default frame::
sage: g[1,1], g[1,2], g[2,2] = 1+x, x*y, 1-x sage: g[:] [ x + 1 x*y] [ x*y -x + 1] sage: g.display() g = (x + 1) dx*dx + x*y dx*dy + x*y dy*dx + (-x + 1) dy*dy
Metric components in a frame different from the manifold's default one::
sage: c_uv.<u,v> = M.chart() # new chart on M sage: xy_to_uv = c_xy.transition_map(c_uv, [x+y, x-y]) ; xy_to_uv Change of coordinates from Chart (M, (x, y)) to Chart (M, (u, v)) sage: uv_to_xy = xy_to_uv.inverse() ; uv_to_xy Change of coordinates from Chart (M, (u, v)) to Chart (M, (x, y)) sage: M.atlas() [Chart (M, (x, y)), Chart (M, (u, v))] sage: M.frames() [Coordinate frame (M, (d/dx,d/dy)), Coordinate frame (M, (d/du,d/dv))] sage: g[c_uv.frame(),:] # metric components in frame c_uv.frame() expressed in M's default chart (x,y) [ 1/2*x*y + 1/2 1/2*x] [ 1/2*x -1/2*x*y + 1/2] sage: g.display(c_uv.frame()) g = (1/2*x*y + 1/2) du*du + 1/2*x du*dv + 1/2*x dv*du + (-1/2*x*y + 1/2) dv*dv sage: g[c_uv.frame(),:,c_uv] # metric components in frame c_uv.frame() expressed in chart (u,v) [ 1/8*u^2 - 1/8*v^2 + 1/2 1/4*u + 1/4*v] [ 1/4*u + 1/4*v -1/8*u^2 + 1/8*v^2 + 1/2] sage: g.display(c_uv.frame(), c_uv) g = (1/8*u^2 - 1/8*v^2 + 1/2) du*du + (1/4*u + 1/4*v) du*dv + (1/4*u + 1/4*v) dv*du + (-1/8*u^2 + 1/8*v^2 + 1/2) dv*dv
The inverse metric is obtained via :meth:`inverse`::
sage: ig = g.inverse() ; ig Tensor field inv_g of type (2,0) on the 2-dimensional differentiable manifold M sage: ig[:] [ (x - 1)/(x^2*y^2 + x^2 - 1) x*y/(x^2*y^2 + x^2 - 1)] [ x*y/(x^2*y^2 + x^2 - 1) -(x + 1)/(x^2*y^2 + x^2 - 1)] sage: ig.display() inv_g = (x - 1)/(x^2*y^2 + x^2 - 1) d/dx*d/dx + x*y/(x^2*y^2 + x^2 - 1) d/dx*d/dy + x*y/(x^2*y^2 + x^2 - 1) d/dy*d/dx - (x + 1)/(x^2*y^2 + x^2 - 1) d/dy*d/dy
""" latex_name=None): r""" Construct a metric on a parallelizable manifold.
TESTS::
sage: M = Manifold(2, 'M') sage: X.<x,y> = M.chart() # makes M parallelizable sage: XM = M.vector_field_module() sage: from sage.manifolds.differentiable.metric import \ ....: PseudoRiemannianMetricParal sage: g = PseudoRiemannianMetricParal(XM, 'g', signature=0); g Lorentzian metric g on the 2-dimensional differentiable manifold M sage: g[0,0], g[1,1] = -(1+x^2), 1+y^2 sage: TestSuite(g).run(skip='_test_category')
.. TODO::
- add a specific parent to the metrics, to fit with the category framework
""" name=name, latex_name=latex_name, sym=(0,1)) # signature: else: raise TypeError("the metric signature must be an integer") raise ValueError("metric signature out of range") if ndim%2 == 0: raise ValueError("the metric signature must be even") else: raise ValueError("the metric signature must be odd") # the pair (n_+, n_-): # Initialization of derived quantities:
r""" Initialize the derived quantities.
TESTS::
sage: M = Manifold(3, 'M') sage: X.<x,y,z> = M.chart() # makes M parallelizable sage: g = M.metric('g') sage: g._init_derived()
""" # Initialization of quantities pertaining to the mother classes:
r""" Delete the derived quantities.
INPUT:
- ``del_restrictions`` -- (default: True) determines whether the restrictions of ``self`` to subdomains are deleted.
TESTS::
sage: M = Manifold(3, 'M') sage: X.<x,y,z> = M.chart() # makes M parallelizable sage: g = M.metric('g') sage: g._del_derived(del_restrictions=False) sage: g._del_derived()
""" # The derived quantities from the mother classes are deleted:
r""" Delete the inverse metric.
TESTS::
sage: M = Manifold(3, 'M') sage: X.<x,y,z> = M.chart() # makes M parallelizable sage: g = M.metric('g') sage: g._del_inverse()
"""
r""" Return the restriction of the metric to some subdomain.
If the restriction has not been defined yet, it is constructed here.
INPUT:
- ``subdomain`` -- open subset `U` of ``self._domain`` (must be an instance of :class:`~sage.manifolds.differentiable.manifold.DifferentiableManifold`) - ``dest_map`` -- (default: ``None``) destination map `\Phi:\ U \rightarrow V`, where `V` is a subdomain of ``self._codomain`` (type: :class:`~sage.manifolds.differentiable.diff_map.DiffMap`) If None, the restriction of ``self._vmodule._dest_map`` to `U` is used.
OUTPUT:
- instance of :class:`PseudoRiemannianMetricParal` representing the restriction.
EXAMPLES:
Restriction of a Lorentzian metric on `\RR^2` to the upper half plane::
sage: M = Manifold(2, 'M') sage: X.<x,y> = M.chart() sage: g = M.lorentzian_metric('g') sage: g[0,0], g[1,1] = -1, 1 sage: U = M.open_subset('U', coord_def={X: y>0}) sage: gU = g.restrict(U); gU Lorentzian metric g on the Open subset U of the 2-dimensional differentiable manifold M sage: gU.signature() 0 sage: gU.display() g = -dx*dx + dy*dy
""" # Construct the restriction at the tensor field level: # the type is correctly handled by TensorFieldParal.restrict, i.e. # resu is of type self.__class__, but the signature is not handled # by TensorFieldParal.restrict; we have to set it here: # Restrictions of derived quantities: resu._ricci_scalar = self._ricci_scalar.restrict(subdomain) resu._weyl = self._weyl.restrict(subdomain) for eps in self._vol_forms: resu._vol_forms.append(eps.restrict(subdomain)) # NB: no initialization of resu._determinants nor # resu._sqrt_abs_dets # The restriction is ready:
r""" Define the metric from a field of symmetric bilinear forms.
INPUT:
- ``symbiform`` -- instance of :class:`~sage.manifolds.differentiable.tensorfield_paral.TensorFieldParal` representing a field of symmetric bilinear forms
EXAMPLES::
sage: M = Manifold(2, 'M') sage: X.<x,y> = M.chart() sage: s = M.sym_bilin_form_field(name='s') sage: s[0,0], s[0,1], s[1,1] = 1+x^2, x*y, 1+y^2 sage: g = M.metric('g') sage: g.set(s) sage: g.display() g = (x^2 + 1) dx*dx + x*y dx*dy + x*y dy*dx + (y^2 + 1) dy*dy
""" raise TypeError("the argument must be a tensor field with " + "values on a parallelizable domain") raise TypeError("the argument must be of tensor type (0,2)") raise TypeError("the argument must be symmetric") raise TypeError("the symmetric bilinear form and the metric are " + "not defined on the same vector field module")
r""" Return the inverse metric.
OUTPUT:
- instance of :class:`~sage.manifolds.differentiable.tensorfield_paral.TensorFieldParal` with tensor_type = (2,0) representing the inverse metric
EXAMPLES:
Inverse metric on a 2-dimensional manifold::
sage: M = Manifold(2, 'M', start_index=1) sage: c_xy.<x,y> = M.chart() sage: g = M.metric('g') sage: g[1,1], g[1,2], g[2,2] = 1+x, x*y, 1-x sage: g[:] # components in the manifold's default frame [ x + 1 x*y] [ x*y -x + 1] sage: ig = g.inverse() ; ig Tensor field inv_g of type (2,0) on the 2-dimensional differentiable manifold M sage: ig[:] [ (x - 1)/(x^2*y^2 + x^2 - 1) x*y/(x^2*y^2 + x^2 - 1)] [ x*y/(x^2*y^2 + x^2 - 1) -(x + 1)/(x^2*y^2 + x^2 - 1)]
If the metric is modified, the inverse metric is automatically updated::
sage: g[1,2] = 0 ; g[:] [ x + 1 0] [ 0 -x + 1] sage: g.inverse()[:] [ 1/(x + 1) 0] [ 0 -1/(x - 1)]
Using SymPy as symbolic engine::
sage: M.set_calculus_method('sympy') sage: g[1,1], g[1,2], g[2,2] = 1+x, x*y, 1-x sage: g[:] # components in the manifold's default frame [ x + 1 x*y] [ x*y -x + 1] sage: g.inverse()[:] [ (x - 1)/(x**2*y**2 + x**2 - 1) x*y/(x**2*y**2 + x**2 - 1)] [ x*y/(x**2*y**2 + x**2 - 1) -(x + 1)/(x**2*y**2 + x**2 - 1)]
""" # Is the inverse metric up to date ? # the computation is necessary output_formatter=fmodule._output_formatter) # of the inverse (keys: comp. indices) # TODO: do the computation without the 'SR' enforcement [[self.comp(frame)[i, j, chart].expr(method='SR') for j in range(si, nsi)] for i in range(si, nsi)]) except (KeyError, ValueError): continue
r""" Return the metric's Ricci scalar.
The Ricci scalar is the scalar field `r` defined from the Ricci tensor `Ric` and the metric tensor `g` by
.. MATH::
r = g^{ij} Ric_{ij}
INPUT:
- ``name`` -- (default: ``None``) name given to the Ricci scalar; if none, it is set to "r(g)", where "g" is the metric's name - ``latex_name`` -- (default: ``None``) LaTeX symbol to denote the Ricci scalar; if none, it is set to "\\mathrm{r}(g)", where "g" is the metric's name
OUTPUT:
- the Ricci scalar `r`, as an instance of :class:`~sage.manifolds.differentiable.scalarfield.DiffScalarField`
EXAMPLES:
Ricci scalar of the standard metric on the 2-sphere::
sage: M = Manifold(2, 'S^2', start_index=1) sage: U = M.open_subset('U') # the complement of a meridian (domain of spherical coordinates) sage: c_spher.<th,ph> = U.chart(r'th:(0,pi):\theta ph:(0,2*pi):\phi') sage: a = var('a') # the sphere radius sage: g = U.metric('g') sage: g[1,1], g[2,2] = a^2, a^2*sin(th)^2 sage: g.display() # standard metric on the 2-sphere of radius a: g = a^2 dth*dth + a^2*sin(th)^2 dph*dph sage: g.ricci_scalar() Scalar field r(g) on the Open subset U of the 2-dimensional differentiable manifold S^2 sage: g.ricci_scalar().display() # The Ricci scalar is constant: r(g): U --> R (th, ph) |--> 2/a^2
""" else: self._ricci_scalar._name = name self._latex_name + r"\right)" else: self._ricci_scalar._latex_name = latex_name |