Feeds:
Posts

SAGE tip: GRmodule. Day 07.

This code is supposed to be (if some one does the work in the future) located in sage.tensor.differential_form_element.

The code presented below is a slight modification of Joris code for differential forms manipulation on SAGE.

Needed modules

from sage.symbolic.ring import SymbolicRing, SR
from sage.rings.ring_element import RingElement
from sage.algebras.algebra_element import AlgebraElement
from sage.rings.integer import Integer
from sage.combinat.permutation import Permutation

The advantage of using this is that, tensors defined here are an Algebra element, not just a python object as in the previous code.

The sage.combinat.permutation won’t be used (yet), but could be useful if tensor symmetries are defined.

TensorFormatter

class TensorFormatter:
r"""
This class contains all the functionality to print a tensor in a
graphically pleasing way.  This class is called by the _latex_ and
_repr_ methods of the Tensor class.
"""
def __init__(self, space):
r"""
Construct a tensor formatter.  See
TensorFormatter for more information.
"""
self._space = space

def repr(self, comp, fun):
r"""
String representation of a primitive tensor, i.e. a function
times a tensor product of d's of the coordinate functions.

INPUT:

- comp -- a subscript of a differential form.

- fun -- the component function of this form.

EXAMPLES::

sage: from sage.tensor.tensor_element import TensorFormatter
sage: x, y, z = var('x, y, z')
sage: U = CoordinatePatch((x, y, z))
sage: D = TensorFormatter(U)
sage: D.repr((0, 1), z^3)
'z^3*dx@dy'

"""

str = "@".join( \
[('d%s' % self._space.coordinate(c).__repr__()) for c in comp])

if fun == 1 and len(comp) > 0:
# We have a non-trivial form whose component function is 1,
# so we just return the formatted form part and ignore the 1.
return str
else:
funstr = fun._repr_()

if not self._is_atomic(funstr):
funstr = '(' + funstr + ')'

if len(str) > 0:
return funstr + "*" + str
else:
return funstr

def latex(self, comp, fun):
r"""
Latex representation of a primitive differential form, i.e. a function
times a tensor product of d's of the coordinate functions.

INPUT:

- comp -- a subscript of a differential form.

- fun -- the component function of this form.

EXAMPLES::

sage: from sage.tensor.tensor_element import TensorFormatter
sage: x, y, z = var('x, y, z')
sage: U = CoordinatePatch((x, y, z))
sage: D = TensorFormatter(U)
sage: D.latex((0, 1), z^3)
'z^{3} d x \otimes d y'

"""

from sage.misc.latex import latex

str = " \otimes ".join( \
[('d %s' % latex(self._space.coordinate(c))) for c in comp])

if fun == 1 and len(comp) > 0:
return str
else:
funstr = latex(fun)

if not self._is_atomic(funstr):
funstr = '(' + funstr + ')'

return funstr + " " + str

def _is_atomic(self, str):
r"""
Helper function to check whether a given string expression
is atomic.

EXAMPLES::

sage: x, y, z = var('x, y, z')
sage: U = CoordinatePatch((x, y, z))
sage: from sage.tensor.tensor_element import TensorFormatter
sage: D = TensorFormatter(U)
sage: D._is_atomic('a + b')
False
sage: D._is_atomic('(a + b)')
True
"""
level = 0
for n, c in enumerate(str):
if c == '(':
level += 1
elif c == ')':
level -= 1

if c == '+' or c == '-':
if level == 0 and n > 0:
return False
return True

The only I’ve changed here is “DifferentialForm” by “Tensor” and “\wedge” by “\otimes”

The above code allows to write the tensor product in a basis, $dx^1\otimes\cdots\otimes dx^n$. The chosen symbol for denoting the tensor product was @.

Tensor Class

This code is incomplete due to:

• I’ve not defined the TensorsAlgebra, which should be done in parallel.
• There are a lot of attributes not presented in this class.
• class Tensor(AlgebraElement):
r"""
Tensor class.
"""

def __init__(self, parent, degree, fun = None):
r"""
Construct a tensor.

INPUT:

- parent -- Parent algebra of tensors.

- degree -- Degree of the tensor.

- fun (default: None) -- Initialize this differential form with the given function.  If the degree is not zero, this argument is silently ignored.

EXAMPLES::

sage: x, y, z = var('x, y, z')
sage: F = Tensors(); F
Algebra of tensors in the variables x, y, z
sage: f = Tensor(F, 0, sin(z)); f
sin(z)

"""

from sage.tensor.tensorss import Tensors
if not isinstance(parent, Tensors):
raise TypeError, "Parent not an algebra of tensors."

RingElement.__init__(self, parent)

self._degree = degree
self._components = {}

if degree == 0 and fun is not None:
self.__setitem__([], fun)

def __getitem__(self, subscript):
r"""
Return a given component of the tensor.

INPUT:

- subscript: subscript of the component.  Must be an integer
or a list of integers.

EXAMPLES::

sage: x, y, z = var('x, y, z')
sage: F = Tensors(); F
Algebra of tensors in the variables x, y, z
sage: f = Tensor(F, 0, sin(x*y)); f
sin(x*y)
sage: f[()]
sin(x*y)
"""

if isinstance(subscript, (Integer, int)):
subscript = (subscript, )
else:
subscript = tuple(subscript)

dim = self.parent().base_space().dim()
if any([s >= dim for s in subscript]):
raise ValueError, "Index out of bounds."

if len(subscript) != self._degree:
raise TypeError, "%s is not a subscript of degree %s" %\
(subscript, self._degree)

"""sign, subscript = sort_subscript(subscript)"""

if subscript in self._components:
return sign*self._components[subscript]
else:
return 0

def __setitem__(self, subscript, fun):
r"""
Modify a given component of the tensor.

INPUT:

- subscript: subscript of the component.  Must be an integer or a list of integers.

EXAMPLES::

sage: F = Tensors(); F
Algebra of tensors in the variables x, y, z
sage: f = Tensor(F, 2)
sage: f[1, 2] = x; f
x*dy@dz
"""

if isinstance(subscript, (Integer, int)):
subscript = (subscript, )
else:
subscript = tuple(subscript)

dim = self.parent().base_space().dim()
if any([s >= dim for s in subscript]):
raise ValueError, "Index out of bounds."

if len(subscript) != self._degree:
raise TypeError, "%s is not a subscript of degree %s" %\
(subscript, self._degree)

"""sign, subscript = sort_subscript(subscript)"""
self._components[subscript] = SR(fun)

Ok, so again I’ve changed “DifferentialForm(s)” by “Tensor(s)”, drop the permutation of indices (’cause tensors do not need to be neither symmetric nor anti-symmetric.

I’ll keep working with this code… It’s all by now. Oh! I’ll post next week the rest of the code based in Sergey’s GRPy.

Enjoy.

Dox

Read Full Post »

SAGE tip: GRmodule. Day 06.

I’ve just updated the SAGE worksheet which uses the definitions described in the previous posts.

• There are some explanations in text format
• The code has been hidden… because is long.
• Moreover… I’ve discover something really amazing! Joris Vankerschaver‘s code of the differential form package in SAGE. Thus, I could use some ideas from Joris’ code to improve GRmodule. Nice, Isn’t it?

Let’s hope I could so something nice this weekend!

Don’t forget check the worksheet, and post some comment for feedback! 😉

Enjoy!

Dox

Read Full Post »

SAGE tip: GRmodule. Day 05

Hi everyone!

This time the Christoffel connection will be defined.

The code

As usual, here is the code:

class Christoffel(Tensor):
'''The class to represent Christoffel Symbols of the second kind. Please
note that while it inherits from Tensor, Christoffel symbols are
NOT tensors'''

def __init__(self,metr,symbol='C',rank=(1,2),sh=(1,-1,-1)):

# The metric
self.g_down = metr

# Since we have a metric we do indeed have a coordinate system
self.rep  = self.g_down.rep

self.g_up = metr.inverse

# Please note that this call will trigger a call to allocate in
# the Tensor class, but the allocate will actually be the allocate
# defined below
super(Christoffel,self).__init__(symbol,rank,sh,coords=metr.coord)

def allocate(self,rank):
Tensor.allocate(self,rank)
# Now that we have allocated things, time to actually compute things
for i in xrange(self.dim):
for k in xrange(self.dim):
for l in xrange(self.dim):
sum = 0
for m in xrange(self.dim):
term1 = diff(self.g_down[m,k],self.g_down.coords[l])
term2 = diff(self.g_down[m,l],self.g_down.coords[k])
term3 = diff(self.g_down[k,l],self.g_down.coords[m])

tr = self.g_up[i,m] * (term1+term2-term3)

sum += tr
res = sum/2
self.components[i,k,l] = res
self.getNonZero()

This code is almost a copy of Sergey’s one, except for the use of xrange instead of np.arange, and the fact that I’ve dropped the minus signs denoting the shape of the tensors.

Sage implementation

This time I won’t present a Python file, but a SAGE file, GRmodule.

Enjoy!

Dox

Read Full Post »

Hi everyone, today I start to write a new class for dealing with non-Abelian forms, i.e., forms with values in a certain Lie-algebra. This could be really useful when computing Yang-Mills theories in physics.

So, I started by defining a new object which have two entries, a differential form and a matrix, and call it nAform. The code I wrote was,

class nAform(object):
def __init__(self, a, b):
self._form = a
self._matrix = b
if isinstance(other, nAform):
if (self._matrix == other._matrix):
return nAform(self._form + other._form, self._matrix)
else:
return NotImplemented
return NotImplemented
def __mul__(self, other):
if isinstance(other, nAform):
return nAform(self._form.wedge(other._form), self._matrix.commutator(other._matrix))
return NotImplemented
def diff(self):
return nAform(self._form.diff(), self._matrix)
def __repr__(self):
return str((self._form, self._matrix))
def __str__(self):
return self.__repr__()

Explanation

One should enter a couple of arguments when defining the nAform object. The __init__ attribute recognize them.

Then an addition attribute is defined, this is incomplete!!!

Another attribute is the multiplication, which take the wedge product of the forms and the commutator of the matrices.

I also implement the exterior derivative on nAform objects.

Finally the __repr__ and __str__ are attributes for returning the data.

TO-DO

• I couldn’t implement the addition of nAform’s if the matrices are different. As
Nicolas M. Thiery note, this objects should define a Monoid (or something quit close to it). But my programming skills are not so developed.
• It would be great if one could define the multiplication by a constant or function.
• The show attribute is not implemented, but one can show either forms of matrices by themselves… I think that is someone knows how they work, would be easy to do that! 😉
• One could try to implement the simplification attributes on differential forms.
• I don’t remember exactly why I was looking for it, but I think could be useful to define and attribute on forms which show the generator, like .gens(), but for a given differential form, say
sage: x, y, z = var('x, y, z')
sage: U = CoordinatePatch((x, y, z)); U
Open subset of R^3 with coordinates x, y, z
sage: F = DifferentialForms(U); F
Algebra of differential forms in the variables x, y, z
sage: F.gens()
(dx, dy, dz)
sage: F.ngens()
3

I’d like one which do like this,

sage: form1 = DifferentialForm(F, 1); form1.ngens()
3
sage: form1 = DifferentialForm(F, 1); form1.gens()
(dx, dy, dz)
form2 = DifferentialForm(F, 2); form2.ngens()
3
sage: form1 = DifferentialForm(F, 1); form1.gens()
(dx/\dy, dy/\dz, dz/\dx)
• Simple SAGE sheet

Worksheet

Hope you can help me with these plenty tasks.

Enjoy!

Dox

Read Full Post »

SAGE tip: GRmodule. Day 04.

Hello everyone. I’ll do it quick, ’cause I’m really tired, and by the way, today… Feb. 8th, I’m turning 30 years old 🙂

First of all, I create a page in my google site for the GRmodule files. Note that the name is GRmodule and not GR-module, because the dash is not admissible in a module name (beginner mistake!)

In the last post…

The Metric class was defined, and the first steps of the implementation were walked.

Inverting the metric

Once the components of the metric are given, one can call the invert method for assigning the inverse of the metric. Don’t know why, but I’d like to keep track of what’s going on! so, I change a little bit the code from Sergei, and my invert method returns the inverse metric tensor,

  def invert(self):
'''Find the inverse of the metric and store the result in a
Metric object self.inverse'''

'''Create a unit matrix of dimension dim'''
temp = sympy.eye(self.dim)

'''Assign the values of the metric to temp'''
for key in self.components.keys():
id = tuple(np.abs(key))
temp[id] = self.components[key]

'''invert the matrix with inv() from sympy'''
inv = temp.inv()
'''convert the matrix in a dictionary'''
inverse = self._dictkeycopy(self.components)
for i in range(self.dim):
for j in range(self.dim):
inverse[i,j] = inv[i,j]
self.inverse = Metric(self.coord,rank=(2,0),sh=(1,1),symbol='g_inv')
self.inverse.components = inverse
return self.inverse

First a temporal matrix is created, and just for assuring it’s invertible, one creates a unit matrix of dimension dim. That’s what the sympy.eye does!

Secondly, the values of the metric are assigned to temp.

The temp matrix is inverted. The sympy command to inver a matrix is inv.

Next, the inverse matrix is converted in a dictionary (just like in previous cases).

Finally, the characteristic of tensor is given to the inverse matrix object… and it’s returned!, i.e., it can be assigned.

In the implementation file the components of the metric are assigned, the metric is invert… and the result of the inversion is printed.

How to run it?!

• Download the latest files and store them in the same folder, e.g., GR.
• Open a terminal and move to the GR folder  cd path/to/GR
• In the terminal type  $python Proof-GR-module.py > Result.txt • Open the Result.txt file, i.e., $ emacs Result.txt
• Enjoy life!

DOX

Happy B-day to me!!! 🙂

Read Full Post »

SAGE tip: GR-module. Day 03.

In the last two post I introduce the classes formalTensor (which has been not that useful) and Tensor, the latter include a series of attributes.

Now is time to get started with GR 😛

Metric Class

I left this class just as Sergei defined, it looked weird to me at the beginning, but again I decipher it 😉 with the help of Reinteract, b.t.w., those who want to play around with Reinteract, it’s in Ubuntu repositories, so just type in your console,

\$ sudo apt-get install reinteract

Ok. Continuing… here is the code,

class Metric(Tensor):
'''Represents a metric. Note that coordinates now MUST be provided'''

def __init__(self,coord,rank=(0,2),sh=(-1,-1),symbol='g'):
self.coord=coord
super(Metric,self).__init__(symbol,rank,sh,coords=coord)

def invert(self):
'''Find the inverse of the metric and store the result in a
Metric object self.inverse'''

'''Create a unit matrix of dimension dim'''
temp = sympy.eye(self.dim)

'''Assign the values of the metric to temp'''
for key in self.components.keys():
id = tuple(np.abs(key))
temp[id] = self.components[key]

'''invert the matrix with inv() from sympy'''
inv = temp.inv()
'''convert the matrix in a dictionary'''
inverse = self._dictkeycopy(self.components)
for i in range(self.dim):
for j in range(self.dim):
inverse[i,j] = inv[i,j]
self.inverse = Metric(self.coord,rank=(2,0),sh=(1,1),symbol='g_inv')
self.inverse.components = inverse

Inherit…

The first I notice, as a non-expert programmer (or non-programmer at all 😛 ) was the word Tensor inside the brackets of the Metric class. This means that Metric is a Tensor… just as a Tensor is an Object. Thus, all methods defined in the Tensor class are applicable to the Metric class.

Coordinate system needed

The __init__ method define the rank, shape and symbol of the metric, which are always the same. However, the set of coordinates must be given.

Example of implementation

Coordinated should be entered like a tuple, and the name of coordinates MUST be “declared” as sympy.Symbol

 t = Symbol('t')
r = Symbol('r')
th = Symbol('theta')
ph = Symbol('phi')
g = Metric((t,r,th,ph))

Given the this data, the Metric class calls the Tensor class and create a covariant rank 2 tensor. And we should give the components. For example, the easiest Schwarzschild metric (the metric is easier than the name!!!),

g[0,0] = -(1 - 2*M/r)
g[1,1] = 1/(1 - 2*M/r)
g[2,2] = r**2
g[3,3] = r**2*(sin(th))**2

WAIT!!!! Don’t forget to declare the mass parameter, so, again…

M = Symbol('M')
g[0,0] = -(1 - 2*M/r)
g[1,1] = 1/(1 - 2*M/r)
g[2,2] = r**2
g[3,3] = r**2*(sin(th))**2

… Sorry guys, it’s again too late! and is already Monday!

Ok, tomorrow we’ll continue with the inverse metric 😉

Get GRmodule.py
and the proof-file.py.

Cya tomorrow!

Dox

Read Full Post »

SAGE tip: GR-module. Day 02.

Ok guys! time to continue 🙂

In the last post I stop just before explaining how to allocate a tensor.

Allocate method

In the allocation method there is a line where the itertools.product command is used. It was the first time I saw such a command, so I opened reinteract (it’s like ipython) to test it.

import itertools
rank = (0,2)
n = rank[0] + rank[1]
indc = list(itertools.product(range(4),repeat=n))
indc

gives the output,
[(0, 0), (0, 1), (0, 2), (0, 3), (1, 0), (1, 1), (1, 2), (1, 3), (2, 0), (2, 1), (2, 2), (2, 3), (3, 0), (3, 1), (3, 2), (3, 3)]

Now I know!!! It gives you a list of the Cartesian product, i.e., a list of all possible values of the indices of a tensor of rank $n$ 🙂 Very clever Sergei!

Next there is a loop, which ends with an assignment. In the loop, every possible combination of indices is stored in a list called mastr and finally a dictionary is created by zipping this mastr list with a zeros list.

This loop assignment complication looks like a work around for differentiate upper and lower components. However, for my mental health… I’ll avoid (at least by now) the shape of the tensor, and the whole loop-assignment can be changed by,

self.components = dict(zip(indc,[0 for i in xrange(len(indc))]))

right after the Cartesian product definition.

_dictkeycopy Attribute

This attribute takes a dictionary (code name hay), extract the keys of the dictionary, and create a new dictionary with those keys but assigning zero values to every key.

:-S Got confused?

it’s like having a copy of the possible indices values of the tensor, zipped with zeros.

getNonZero attribute

This is another important attribute!!! We all know that not all components of tensor in GR are different than zero, so this attribute “print” only those that are non-vanishing 😉

I think it’s more or less clear what it does. It sweeps the values on the dictionary, and if the value is non-zero, it stores both key and value in a new dictionary 🙂

__str__ attribute

It prints the result of getNonZero. 😛

Now what?!

On the next post I’ll start by defining a metric tensor. However, it could be useful to define instead a set of vielbein and a signature… Perhaps later!

Get GR-module.py (version 2)

Got to work right now!

Enjoy life!

DOX

Read Full Post »

Older Posts »