Via the command type, we get the type of the Python object corresponding to M (here the Python class DifferentiableManifold_with_category):
# In[7]: type(M) # We can also ask for the category of M and see that it is the category of smooth manifolds over $\mathbb{R}$: # In[8]: category(M) #The indices on the manifold are generated by the method irange(), to be used in loops:
# In[9]: for i in M.irange(): print(i) #If the parameter start_index had not been specified, the default range of the indices would have been $\{0,1,2\}$ instead:
# In[10]: M0 = Manifold(3, 'M', r'\mathcal{M}') for i in M0.irange(): print(i) #Let us assume that the manifold $\mathcal{M}$ can be covered by a single chart (other cases are discussed below); the chart is declared as follows:
# In[11]: X.The chart is displayed as a pair formed by the open set covered by it (here the whole manifold) and the coordinates:
# In[12]: print(X) # In[13]: X #The coordinates can be accessed individually, by means of their indices, following the convention defined by start_index=1 in the manifold's definition:
# In[14]: X[1] # In[15]: X[2] # In[16]: X[3] #The full set of coordinates is obtained by means of the operator [:]:
# In[17]: X[:] # Thanks to the operator `Real-valued functions of the chart coordinates (mathematically speaking, functions defined on the chart codomain) are generated via the method function() acting on the chart:
# In[20]: f = X.function(x+y^2+z^3) ; f # In[21]: f.display() # In[22]: f(1,2,3) # They belong to SageMath class `ChartFunction`: # In[23]: type(f) #and differ from SageMath standard symbolic functions by automatic simplifications in all operations. For instance, adding the two symbolic functions
# In[24]: f0(x,y,z) = cos(x)^2 ; g0(x,y,z) = sin(x)^2 #results in
# In[25]: f0 + g0 # while the sum of the corresponding functions in the class `ChartFunction` is automatically simplified: # In[26]: f1 = X.function(cos(x)^2) ; g1 = X.function(sin(x)^2) f1 + g1 #To get the same output with symbolic functions, one has to invoke the method simplify_trig():
# In[27]: (f0 + g0).simplify_trig() #Another difference regards the display; if we ask for the symbolic function f0, we get:
# In[28]: f0 #while if we ask for the chart function f1, we get only the coordinate expression:
# In[29]: f1 #To get an output similar to that of f0, one should call the method display():
# In[30]: f1.display() # Note that the method `expr()` returns the underlying symbolic expression: # In[31]: f1.expr() # In the present case, it is a SymPy object: # In[32]: type(f1.expr()) # ### Introducing a second chart on the manifold # # Let us first consider an open subset of $\mathcal{M}$, for instance the complement $U$ of the region defined by $\{y=0, x\geq 0\}$ (note that `(y!=0, x<0)` stands for $y\not=0$ OR $x<0$; the condition $y\not=0$ AND $x<0$ would have been written `[y!=0, x<0]` instead): # In[33]: U = M.open_subset('U', coord_def={X: (y!=0, x<0)}) # Let us call `X_U` the restriction of the chart `X` to the open subset $U$: # In[34]: X_U = X.restrict(U) ; X_U #We introduce another chart on $U$, with spherical-type coordinates $(r,\theta,\phi)$:
# In[35]: Y.The function chart() has now some argument; it is a string, which contains specific LaTeX symbols, hence the prefix 'r' to it (for raw string). It also contains the coordinate ranges, since they are different from the default value, which is $(-\infty, +\infty)$. For a given coordinate, the various fields are separated by the character ':' and a space character separates the coordinates. Note that for the coordinate $r$, there are only two fields, since the LaTeX symbol has not to be specified. The LaTeX symbols are used for the outputs:
# In[36]: th, ph # In[37]: Y[2], Y[3] #The declared coordinate ranges are now known to Sage, as we may check by means of the command assumptions():
# In[38]: assumptions() #They are used in simplifications:
# In[39]: simplify(abs(r)) # In[40]: simplify(abs(x)) # no simplification occurs since x can take any value in R #After having been declared, the chart Y can be fully specified by its relation to the chart X_U, via a transition map:
# In[41]: transit_Y_to_X = Y.transition_map(X_U, [r*sin(th)*cos(ph), r*sin(th)*sin(ph), r*cos(th)]) # In[42]: transit_Y_to_X # In[43]: transit_Y_to_X.display() #The inverse of the transition map can be specified by means of the method set_inverse():
# In[44]: transit_Y_to_X.set_inverse(sqrt(x^2+y^2+z^2), atan2(sqrt(x^2+y^2),z), atan2(y, x)) transit_Y_to_X.inverse().display() #At this stage, the manifold's atlas (the "user atlas", not the maximal atlas!) contains three charts:
# In[45]: M.atlas() #The first chart defined on the manifold is considered as the manifold's default chart (it can be changed by the method set_default_chart()):
# In[46]: M.default_chart() #Each open subset has its own atlas (since an open subset of a manifold is a manifold by itself):
# In[47]: U.atlas() # In[48]: U.default_chart() # We can draw the chart $Y$ in terms of the chart $X$. # Let us first define a viewer for 3D plots (use `'threejs'` or `'jmol'` for interactive 3D graphics): # In[49]: viewer3D = 'threejs' # must be 'threejs', 'jmol', 'tachyon' or None (default) # The plot shows lines of constant coordinates from the $Y$ chart in a "Cartesian frame" based on the $X$ coordinates: # In[50]: graph = Y.plot(X) show(graph, viewer=viewer3D, online=True) # The command plot() allows for many options, to control the number of coordinate lines to be drawn, their style and color, as well as the coordinate ranges (cf. the [list of all options](http://doc.sagemath.org/html/en/reference/manifolds/sage/manifolds/chart.html#sage.manifolds.chart.RealChart.plot)): # In[51]: graph = Y.plot(X, ranges={r:(1,2), th:(0,pi/2)}, number_values=4, color={r:'blue', th:'green', ph:'red'}) show(graph, aspect_ratio=1, viewer=viewer3D, online=True) #Conversly, the chart $X|_{U}$ can be plotted in terms of the chart $Y$ (this is not possible for the whole chart $X$ since its domain is larger than that of chart $Y$):
# In[52]: #graph = X_U.plot(Y) #show(graph, viewer=viewer3D, online=True, axes_labels=['r','theta','phi']) #A point on $\mathcal{M}$ is defined by its coordinates in a given chart:
# In[53]: p = M.point((1,2,-1), chart=X, name='p') ; print(p) ; p #Since $X=(\mathcal{M}, (x,y,z))$ is the manifold's default chart, its name can be omitted:
# In[54]: p = M.point((1,2,-1), name='p') ; print(p) ; p #Of course, $p$ belongs to $\mathcal{M}$:
# In[55]: p in M #It is also in $U$:
# In[56]: p in U #Indeed the coordinates of $p$ have $y\not=0$:
# In[57]: p.coord(X) #Note in passing that since $X$ is the default chart on $\mathcal{M}$, its name can be omitted in the arguments of coord():
# In[58]: p.coord() #The coordinates of $p$ can also be obtained by letting the chart acting of the point (from the very definition of a chart!):
# In[59]: X(p) #Let $q$ be a point with $y = 0$ and $x \geq 0$:
# In[60]: q = M.point((1,0,2), name='q') #This time, the point does not belong to $U$:
# In[61]: q in U #Accordingly, we cannot ask for the coordinates of $q$ in the chart $Y=(U, (r,\theta,\phi))$:
# In[62]: try: q.coord(Y) except ValueError as exc: print("Error: " + str(exc)) #but we can for point $p$:
# In[63]: p.coord(Y) # In[64]: Y(p) #Points can be compared:
# In[65]: q == p # In[66]: p1 = U.point((sqrt(6), pi-atan(sqrt(5)), atan(2)), Y) p1 == p #In SageMath's terminology, points are elements, whose parents are the manifold on which they have been defined:
# In[67]: p.parent() # In[68]: q.parent() # In[69]: p1.parent() #A scalar field is a differentiable mapping $U \longrightarrow \mathbb{R}$, where $U$ is an open subset of $\mathcal{M}$.
#The scalar field is defined by its expressions in terms of charts covering its domain (in general more than one chart is necessary to cover all the domain):
# In[70]: f = U.scalar_field({X_U: x+y^2+z^3}, name='f') ; print(f) #The coordinate expressions of the scalar field are passed as a Python dictionary, with the charts as keys, hence the writing {X_U: x+y^2+z^3}.
#Since in the present case, there is only one chart in the dictionary, an alternative writing is
# In[71]: f = U.scalar_field(x+y^2+z^3, chart=X_U, name='f') ; print(f) #Since X_U is the domain's default chart, it can be omitted in the above declaration:
# In[72]: f = U.scalar_field(x+y^2+z^3, name='f') ; print(f) # As a mapping $U\subset\mathcal{M}\longrightarrow\mathbb{R}$, a scalar field acts on points, not on coordinates: # In[73]: f(p) # The method `display()` provides the expression of the scalar field in terms of a given chart: # In[74]: f.display(X_U) # If no argument is provided, the method `display()` shows the coordinate expression of the scalar field in all the charts defined on the domain (except for *subcharts*, i.e. the restrictions of some chart to a subdomain): # In[75]: f.display() #Note that the expression of $f$ in terms of the coordinates $(r,\theta,\phi)$ has not been provided by the user but has been automatically computed by means of the change-of-coordinate formula declared above in the transition map.
# In[76]: f.display(Y) # In each chart, the scalar field is represented by a function of the chart coordinates (an object of the type `ChartFunction` described above), which is accessible via the method `coord_function()`: # In[77]: f.coord_function(X_U) # In[78]: f.coord_function(X_U).display() # In[79]: f.coord_function(Y) # In[80]: f.coord_function(Y).display() # The "raw" symbolic expression (here a SymPy object) is returned by the method `expr()`: # In[81]: f.expr(X_U) # In[82]: f.expr(Y) # In[83]: f.expr(Y) is f.coord_function(Y).expr() #A scalar field can also be defined by some unspecified function of the coordinates:
# In[84]: h = U.scalar_field(function('H')(x, y, z), name='h') ; print(h) # In[85]: h.display() # In[86]: h.display(Y) # In[87]: h(p) # remember that p is the point of coordinates (1,2,-1) in the chart X_U #The parent of $f$ is the set $C^\infty(U)$ of all smooth scalar fields on $U$, which is a commutative algebra over $\mathbb{R}$:
# In[88]: CU = f.parent() ; CU # In[89]: print(CU) # In[90]: CU.category() #The base ring of the algebra is the field $\mathbb{R}$, which is represented here by SageMath's Symbolic Ring (SR):
#Arithmetic operations on scalar fields are defined through the algebra structure:
# In[92]: s = f + 2*h ; print(s) # In[93]: s.display() #The tangent vector space to the manifold at point $p$ is obtained as follows:
# In[94]: Tp = M.tangent_space(p) ; Tp # In[95]: print(Tp) #$T_p\, \mathcal{M}$ is a 2-dimensional vector space over $\mathbb{R}$ (represented here by SageMath's Symbolic Ring (SR)) :
# In[96]: print(Tp.category()) # In[97]: Tp.dim() #$T_p\, \mathcal{M}$ is automatically endowed with vector bases deduced from the vector frames defined around the point:
# In[98]: Tp.bases() #For the tangent space at the point $q$, on the contrary, there is only one pre-defined basis, since $q$ is not in the domain $U$ of the frame associated with coordinates $(r,\theta,\phi)$:
# In[99]: Tq = M.tangent_space(q) Tq.bases() #A random element:
# In[100]: v = Tp.an_element() ; print(v) # In[101]: v.display() # In[102]: u = Tq.an_element() ; print(u) # In[103]: u.display() #Note that, despite what the above simplified writing may suggest (the mention of the point $p$ or $q$ is omitted in the basis vectors), $u$ and $v$ are different vectors, for they belong to different vector spaces:
# In[104]: v.parent() # In[105]: u.parent() #In particular, it is not possible to add $u$ and $v$:
# In[106]: try: s = u + v except TypeError as exc: print("Error: " + str(exc)) #Each chart defines a vector frame on the chart domain: the so-called coordinate basis:
# In[107]: X.frame() # In[108]: X.frame().domain() # this frame is defined on the whole manifold # In[109]: Y.frame() # In[110]: Y.frame().domain() # this frame is defined only on U #The list of frames defined on a given open subset is returned by the method frames():
# In[111]: M.frames() # In[112]: U.frames() # In[113]: M.default_frame() #Unless otherwise specified (via the command set_default_frame()), the default frame is that associated with the default chart:
# In[114]: M.default_frame() is M.default_chart().frame() # In[115]: U.default_frame() is U.default_chart().frame() #Individual elements of a frame can be accessed by means of their indices:
# In[116]: e = U.default_frame() ; e2 = e[2] ; e2 # In[117]: print(e2) #We may define a new vector field as follows:
# In[118]: v = e[2] + 2*x*e[3] ; print(v) # In[119]: v.display() #A vector field can be defined by its components with respect to a given vector frame. When the latter is not specified, the open set's default frame is of course assumed:
# In[120]: v = U.vector_field(name='v') # vector field defined on the open set U v[1] = 1+y v[2] = -x v[3] = x*y*z v.display() #Vector fields on $U$ are Sage element objects, whose parent is the set $\mathfrak{X}(U)$ of vector fields defined on $U$:
# In[121]: v.parent() #The set $\mathfrak{X}(U)$ is a module over the commutative algebra $C^\infty(U)$ of scalar fields on $U$:
# In[122]: print(v.parent()) # In[123]: print(v.parent().category()) # In[124]: v.parent().base_ring() #A vector field acts on scalar fields:
# In[125]: f.display() # In[126]: s = v(f) ; print(s) # In[127]: s.display() # In[128]: e[3].display() # In[129]: e[3](f).display() #Unset components are assumed to be zero:
# In[130]: w = U.vector_field(name='w') w[2] = 3 w.display() #A vector field on $U$ can be expanded in the vector frame associated with the chart $(r,\theta,\phi)$:
# In[131]: v.display(Y.frame()) #By default, the components are expressed in terms of the default coordinates $(x,y,z)$. To express them in terms of the coordinates $(r,\theta,\phi)$, one should add the corresponding chart as the second argument of the method display():
# In[132]: v.display(Y.frame(), Y) # In[133]: for i in M.irange(): show(e[i].display(Y.frame(), Y)) #The components of a tensor field w.r.t. the default frame can also be obtained as a list, via the command [:]:
# In[134]: v[:] #An alternative is to use the method display_comp():
# In[135]: v.display_comp() #To obtain the components w.r.t. to another frame, one may go through the method comp() and specify the frame:
# In[136]: v.comp(Y.frame())[:] #However a shortcut is to provide the frame as the first argument of the square brackets:
# In[137]: v[Y.frame(), :] # In[138]: v.display_comp(Y.frame()) #Components are shown expressed in terms of the default's coordinates; to get them in terms of the coordinates $(r,\theta,\phi)$ instead, add the chart name as the last argument in the square brackets:
# In[139]: v[Y.frame(), :, Y] #or specify the chart in display_comp():
A vector field can be defined with components being unspecified functions of the coordinates:
# In[144]: u = U.vector_field(name='u') u[:] = [function('u_x')(x,y,z), function('u_y')(x,y,z), function('u_z')(x,y,z)] u.display() # In[145]: s = v + u ; s.set_name('s') ; s.display() #The value of a vector field at some point of the manifold is obtained via the method at():
# In[146]: vp = v.at(p) ; print(vp) # In[147]: vp.display() #Indeed, recall that, w.r.t. chart X_U=$(x,y,z)$, the coordinates of the point $p$ and the components of the vector field $v$ are
# In[148]: p.coord(X_U) # In[149]: v.display(X_U.frame(), X_U) #Note that to simplify the writing, the symbol used to denote the value of the vector field at point $p$ is the same as that of the vector field itself (namely $v$); this can be changed by the method set_name():
# In[150]: vp.set_name(latex_name='v|_p') vp.display() #Of course, $v|_p$ belongs to the tangent space at $p$:
# In[151]: vp.parent() # In[152]: vp in M.tangent_space(p) # In[153]: up = u.at(p) ; print(up) # In[154]: up.display() #A 1-form on $\mathcal{M}$ is a field of linear forms. For instance, it can be the differential of a scalar field:
# In[155]: df = f.differential() ; print(df) # In[156]: df.display() #In the above writing, the 1-form is expanded over the basis $(\mathrm{d}x, \mathrm{d}y, \mathrm{d}z)$ associated with the chart $(x,y,z)$. This basis can be accessed via the method coframe():
# In[157]: dX = X.coframe() ; dX #The list of all coframes defined on a given manifold open subset is returned by the method coframes():
# In[158]: M.coframes() #As for a vector field, the value of the differential form at some point on the manifold is obtained by the method at():
# In[159]: dfp = df.at(p) ; print(dfp) # In[160]: dfp.display() #Recall that
# In[161]: p.coord() #The linear form $\mathrm{d}f|_p$ belongs to the dual of the tangent vector space at $p$:
# In[162]: dfp.parent() # In[163]: dfp.parent() is M.tangent_space(p).dual() #As such, it is acting on vectors at $p$, yielding a real number:
# In[164]: print(vp) ; vp.display() # In[165]: dfp(vp) # In[166]: print(up) ; up.display() # In[167]: dfp(up) #The differential 1-form of the unspecified scalar field $h$:
# In[168]: h.display() # In[169]: h.display() ; dh = h.differential() ; dh.display() #A 1-form can also be defined from scratch:
# In[170]: om = U.one_form('omega', r'\omega') ; print(om) #It can be specified by providing its components in a given coframe:
# In[171]: om[:] = [x^2+y^2, z, x-z] # components in the default coframe (dx,dy,dz) om.display() #Of course, one may set the components in a frame different from the default one:
# In[172]: om[Y.frame(), :, Y] = [r*sin(th)*cos(ph), 0, r*sin(th)*sin(ph)] om.display(Y.frame(), Y) #The components in the coframe $(\mathrm{d}x,\mathrm{d}y,\mathrm{d}z)$ are updated automatically:
# In[173]: om.display() #Let us revert to the values set previously:
# In[174]: om[:] = [x^2+y^2, z, x-z] om.display() #This time, the components in the coframe $(\mathrm{d}r, \mathrm{d}\theta,\mathrm{d}\phi)$ are those that are updated:
# In[175]: om.display(Y.frame(), Y) #A 1-form acts on vector fields, resulting in a scalar field:
# In[176]: v.display() ; om.display() ; print(om(v)) ; om(v).display() # In[177]: df.display() ; print(df(v)) ; df(v).display() # In[178]: u.display() ; om(u).display() #In the case of a differential 1-form, the following identity holds:
# In[179]: df(v) == v(f) #1-forms are Sage element objects, whose parent is the $C^\infty(U)$-module $\Omega^{1}(U)$ of all 1-forms defined on $U$:
# In[180]: df.parent() # In[181]: print(df.parent()) # In[182]: print(om.parent()) #$\Omega^{1}(U)$ is actually the dual of the free module $\mathfrak{X}(U)$:
# In[183]: df.parent() is v.parent().dual() #The exterior product of two 1-forms is taken via the method wedge() and results in a 2-form:
# In[184]: a = om.wedge(df) ; print(a) ; a.display() #A matrix view of the components:
# In[185]: a[:] #Displaying only the non-vanishing components, skipping the redundant ones (i.e. those that can be deduced by antisymmetry):
# In[186]: a.display_comp(only_nonredundant=True) #The 2-form $\omega\wedge\mathrm{d}f$ can be expanded on the $(\mathrm{d}r,\mathrm{d}\theta,\mathrm{d}\phi)$ coframe:
# In[187]: a.display(Y.frame(), Y) #As a 2-form, $A:=\omega\wedge\mathrm{d}f$ can be applied to a pair of vectors and is antisymmetric:
# In[188]: a.set_name('A') print(a(u,v)) ; a(u,v).display() # In[189]: a(u,v) == - a(v,u) # In[190]: a.symmetries() #The exterior derivative of a differential form:
# In[191]: dom = om.exterior_derivative() ; print(dom) ; dom.display() #Instead of invoking the method exterior_derivative(), one can use the function xder, after having imported it from sage.manifolds.utilities:
# In[192]: from sage.manifolds.utilities import xder dom = xder(om) # In[193]: da = xder(a) ; print(da) ; da.display() #The exterior derivative is nilpotent:
# In[194]: ddf = xder(df) ; ddf.display() # In[195]: ddom = xder(dom) ; ddom.display() #The Lie derivative of any tensor field with respect to a vector field is computed by the method lie_derivative(), with the vector field as the argument:
# In[196]: lv_om = om.lie_derivative(v) ; print(lv_om) ; lv_om.display() # In[197]: # lu_dh = dh.lie_derivative(u) ; print(lu_dh) ; lu_dh.display() #Let us check Cartan identity on the 1-form $\omega$:
#$\mathcal{L}_v \omega = v\cdot \mathrm{d}\omega + \mathrm{d}\langle \omega, v\rangle$
#and on the 2-form $A$:
#$\mathcal{L}_v A = v\cdot \mathrm{d}A + \mathrm{d}(v\cdot A)$
# In[198]: om.lie_derivative(v) == v.contract(xder(om)) + xder(om(v)) # In[199]: a.lie_derivative(v) == v.contract(xder(a)) + xder(v.contract(a)) #The Lie derivative of a vector field along another one is the commutator of the two vectors fields:
# In[200]: # v.lie_derivative(u)(f) == u(v(f)) - v(u(f)) #Up to now, we have encountered tensor fields
#More generally, tensor fields of any type $(p,q)$ can be introduced in SageMath. For instance a tensor field of type (1,2) on the open subset $U$ is declared as follows:
# In[201]: t = U.tensor_field(1, 2, name='T') ; print(t) #As for vectors or 1-forms, the tensor's components with respect to the domain's default frame are set by means of square brackets:
# In[202]: t[1,2,1] = 1 + x^2 t[3,2,1] = x*y*z #Unset components are zero:
# In[203]: t.display() # In[204]: t[:] #Display of the nonzero components:
# In[205]: t.display_comp() #Double square brackets return the component (still w.r.t. the default frame) as a scalar field, while single square brackets return the expression of this scalar field in terms of the domain's default coordinates:
# In[206]: print(t[[1,2,1]]) ; t[[1,2,1]].display() # In[207]: print(t[1,2,1]) ; t[1,2,1] #A tensor field of type (1,2) maps a 3-tuple (1-form, vector field, vector field) to a scalar field:
# In[208]: print(t(om, u, v)) ; t(om, u, v).display() #As for vectors and differential forms, the tensor components can be taken in any frame defined on the manifold:
# In[209]: t[Y.frame(), 1,1,1, Y] #The tensor product $\otimes$ is denoted by `*`:
# In[210]: v.tensor_type() ; a.tensor_type() # In[211]: b = v*a ; print(b) ; b #The tensor product preserves the (anti)symmetries: since $A$ is a 2-form, it is antisymmetric with respect to its two arguments (positions 0 and 1); as a result, b is antisymmetric with respect to its last two arguments (positions 1 and 2):
# In[212]: a.symmetries() # In[213]: b.symmetries() #Standard tensor arithmetics is implemented:
# In[214]: s = - t + 2*f* b ; print(s) #Tensor contractions are dealt with by the methods trace() and contract(): for instance, let us contract the tensor $T$ w.r.t. its first two arguments (positions 0 and 1), i.e. let us form the tensor $c$ of components $c_i = T^k_{\ \, k i}$:
# In[215]: c = t.trace(0,1) print(c) #An alternative to the writing trace(0,1) is to use the index notation to denote the contraction: the indices are given in a string inside the [] operator, with '^' in front of the contravariant indices and '_' in front of the covariant ones:
# In[216]: c1 = t['^k_ki'] print(c1) c1 == c #The contraction is performed on the repeated index (here k); the letter denoting the remaining index (here i) is arbitrary:
# In[217]: t['^k_kj'] == c # In[218]: t['^b_ba'] == c #It can even be replaced by a dot:
# In[219]: t['^k_k.'] == c #LaTeX notations are allowed:
# In[220]: t['^{k}_{ki}'] == c #The contraction $T^i_{\ j k} v^k$ of the tensor fields $T$ and $v$ is taken as follows (2 refers to the last index position of $T$ and 0 to the only index position of v):
# In[221]: tv = t.contract(2, v, 0) print(tv) #Since 2 corresponds to the last index position of $T$ and 0 to the first index position of $v$, a shortcut for the above is
# In[222]: tv1 = t.contract(v) print(tv1) # In[223]: tv1 == tv #Instead of contract(), the index notation, combined with the * operator, can be used to denote the contraction:
# In[224]: t['^i_jk']*v['^k'] == tv #The non-repeated indices can be replaced by dots:
# In[225]: t['^._.k']*v['^k'] == tv #A Riemannian metric on the manifold $\mathcal{M}$ is declared as follows:
# In[226]: g = M.riemannian_metric('g') print(g) #It is a symmetric tensor field of type (0,2):
# In[227]: g.parent() # In[228]: print(g.parent()) # In[229]: g.symmetries() #The metric is initialized by its components with respect to some vector frame. For instance, using the default frame of $\mathcal{M}$:
# In[230]: g[1,1], g[2,2], g[3,3] = 1, 1, 1 g.display() #The components w.r.t. another vector frame are obtained as for any tensor field:
# In[231]: g.display(Y.frame(), Y) #Of course, the metric acts on vector pairs:
# In[232]: u.display() ; v.display(); print(g(u,v)) ; g(u,v).display() #The Levi-Civita connection associated to the metric $g$:
# In[233]: nabla = g.connection() print(nabla) ; nabla #The Christoffel symbols with respect to the manifold's default coordinates:
# In[234]: nabla.coef()[:] #The Christoffel symbols with respect to the coordinates $(r,\theta,\phi)$:
# In[235]: nabla.coef(Y.frame())[:, Y] # A nice view is obtained via the method `display()` (by default, only the nonzero connection coefficients are shown): # In[236]: nabla.display(frame=Y.frame(), chart=Y) #The connection acting as a covariant derivative:
# In[237]: nab_v = nabla(v) print(nab_v) # nab_v.display() # does not work with SymPy (trigsimp runs forever...) #Being a Levi-Civita connection, $\nabla_g$ is torsion.free:
# In[238]: print(nabla.torsion()) ; nabla.torsion().display() #In the present case, it is also flat:
# In[239]: print(nabla.riemann()) ; nabla.riemann().display() #Let us consider a non-flat metric, by changing $g_{rr}$ to $1/(1+r^2)$:
# In[240]: g[Y.frame(), 1,1, Y] = 1/(1+r^2) g.display(Y.frame(), Y) #For convenience, we change the default chart on the domain $U$ to Y=$(U,(r,\theta,\phi))$:
# In[241]: U.set_default_chart(Y) #In this way, we do not have to specify Y when asking for coordinate expressions in terms of $(r,\theta,\phi)$:
# In[242]: g.display(Y.frame()) #We recognize the metric of the hyperbolic space $\mathbb{H}^3$. Its expression in terms of the chart $(U,(x,y,z))$ is
# In[243]: g.display(X_U.frame(), X_U) #A matrix view of the components may be more appropriate:
# In[244]: g[X_U.frame(), :, X_U] #We extend these components, a priori defined only on $U$, to the whole manifold $\mathcal{M}$, by demanding the same coordinate expressions in the frame associated to the chart X=$(\mathcal{M},(x,y,z))$:
# In[245]: g.add_comp_by_continuation(X.frame(), U, X) g.display() #The Levi-Civita connection is automatically recomputed, after the change in $g$:
# In[246]: nabla = g.connection() #In particular, the Christoffel symbols are different:
# In[247]: nabla.display(only_nonredundant=True) # In[248]: nabla.display(frame=Y.frame(), chart=Y, only_nonredundant=True) #The Riemann tensor is now
# In[249]: Riem = nabla.riemann() print(Riem) ; Riem.display(Y.frame()) #Note that it can be accessed directely via the metric, without any explicit mention of the connection:
# In[250]: g.riemann() is nabla.riemann() #The Ricci tensor is
# In[251]: Ric = g.ricci() print(Ric) ; Ric.display(Y.frame()) #The Weyl tensor is:
# In[252]: C = g.weyl() print(C) ; C.display() #The Weyl tensor vanishes identically because the dimension of $\mathcal{M}$ is 3.
#Finally, the Ricci scalar is
# In[253]: R = g.ricci_scalar() print(R) ; R.display() #We recover the fact that $\mathbb{H}^3$ is a Riemannian manifold of constant negative curvature.
# #The most important tensor transformation induced by the metric $g$ is the so-called musical isomorphism, or index raising and index lowering:
# In[254]: print(t) # In[255]: t.display() # In[256]: t.display(X_U.frame(), X_U) #Raising the last index of $T$ with $g$:
# In[257]: s = t.up(g, 2) print(s) #Raising all the covariant indices of $T$ (i.e. those at the positions 1 and 2):
# In[258]: s = t.up(g) print(s) # In[259]: s = t.down(g) print(s) #The volume 3-form (Levi-Civita tensor) associated with the metric $g$ is
# In[260]: epsilon = g.volume_form() print(epsilon) ; epsilon.display() # In[261]: epsilon.display(Y.frame()) # In[262]: print(f) ; f.display() # In[263]: sf = f.hodge_dual(g) print(sf) ; sf.display() #We check the classical formula $\star f = f\, \epsilon_g$, or, more precisely, $\star f = f\, \epsilon_g|_U$ (for $f$ is defined on $U$ only):
# In[264]: sf == f * epsilon.restrict(U) #The Hodge dual of a 1-form is a 2-form:
# In[265]: print(om) ; om.display() # In[266]: som = om.hodge_dual(g) print(som) ; som.display() #The Hodge dual of a 2-form is a 1-form:
# In[267]: print(a) # In[268]: sa = a.hodge_dual(g) print(sa) ; sa.display() #Finally, the Hodge dual of a 3-form is a 0-form:
# In[269]: print(da) ; da.display() # In[270]: sda = da.hodge_dual(g) print(sda) ; sda.display() #In dimension 3 and for a Riemannian metric, the Hodge star is idempotent:
# In[271]: sf.hodge_dual(g) == f # In[272]: som.hodge_dual(g) == om # In[273]: sa.hodge_dual(g) == a # In[274]: sda.hodge_dual(g) == da #To get the list of functions (methods) that can be called on a object, type the name of the object, followed by a dot and the TAB key, e.g.
#sa.
To get information on an object or a method, use the question mark:
# In[275]: get_ipython().run_line_magic('pinfo', 'nabla') # In[276]: get_ipython().run_line_magic('pinfo', 'g.ricci_scalar') #Using a double question mark leads directly to the Python source code (SageMath is open source, isn't it?)
Have a look at the examples on SageManifolds page, especially the 2-dimensional sphere example for usage on a non-parallelizable manifold (each scalar field has to be defined in at least two coordinate charts, the module $\mathfrak{X}(\mathcal{M})$ is no longer free and each tensor field has to be defined in at least two vector frames).