GSoC Ideas : Tensor/Linear Algebra

110 views
Skip to first unread message

Matthew Rocklin

unread,
Feb 20, 2012, 4:19:54 PM2/20/12
to sy...@googlegroups.com
Hi Everyone, 

I would like to create a general tensor/linear algebra framework for SymPy. I'd like to hear ideas from the community about this. 

We already have a few linear algebraic projects within SymPy (i.e. Matrices, SparseMatrices, MatExprs, Indexed/IndexedBase code generation, Physics stuff, Geometric Algebra (sort of)) but they don't communicate well. It would be nice to create a general and abstract framework off of which these projects and others could hang and interact more naturally.

I'm writing the community about this for two reasons. 
Reason one: I'd like feedback as to whether or not this sort of undertaking is a good idea. If it is I'd welcome some thoughts on how it should be done and how it could be useful for future work. 

Reason two: I think I can separate this work into a few pieces, each of which would make for a good GSoC project for this year or next. Is this endeavor something into which the community would want to invest resources?

Here are some projects that interest me
  1. Framework design - we need a sufficiently general framework (this is hard and probably has to be half completed before GSoC time)
  2. Abstract Vector Spaces
  3. Tensor Math - Krastanov was talking about this and I think it's a great idea. There is a lot of good multilinear algebra out there that SymPy doesn't currently touch at all. 
  4. General storage - Efficient NDArray classes (dense, mutable, sparse, functional, numpy, external programs) - views of NDArrays (transpose, slices). 
  5. Theorem proving type system for tensors/matrices http://scicomp.stackexchange.com/questions/74/symbolic-software-packagbasic es-for-matrix-expressions
I've dumped some thoughts on the following wiki page

Comments or questions are welcome. 

-Matt

Aaron Meurer

unread,
Feb 20, 2012, 9:40:51 PM2/20/12
to sy...@googlegroups.com
On Mon, Feb 20, 2012 at 2:19 PM, Matthew Rocklin <mroc...@gmail.com> wrote:
> Hi Everyone,
>
> I would like to create a general tensor/linear algebra framework for SymPy.
> I'd like to hear ideas from the community about this.
>
> We already have a few linear algebraic projects within SymPy
> (i.e. Matrices, SparseMatrices, MatExprs, Indexed/IndexedBase code
> generation, Physics stuff, Geometric Algebra (sort of)) but they don't
> communicate well. It would be nice to create a general and abstract
> framework off of which these projects and others could hang and interact
> more naturally.
>
> I'm writing the community about this for two reasons.
> Reason one: I'd like feedback as to whether or not this sort of undertaking
> is a good idea. If it is I'd welcome some thoughts on how it should be done
> and how it could be useful for future work.

Definitely. One problem right now is that a lot of modules duplicate
work, because we don't really have a good centrailzed module for
things. For example, over GCI, a student merged together three
independent KroneckerDelta implementations (one in
sympy/physics/quantum, one in sympy/physics/secondquant, and one in
sympy/functions/special/tensor_functions.py). No doubt there are
other things still duplicated.

That's also why even if we are implementing some of these things to
help with physics, we should try to separate mathematical concepts
from physical concepts in the implementation.

>
> Reason two: I think I can separate this work into a few pieces, each of
> which would make for a good GSoC project for this year or next. Is this
> endeavor something into which the community would want to invest resources?

I think so. Some projects may depend on others (e.g., we're limited in
what we can do with slow matrices). But feel free to do this and add
the ideas to the GSoC ideas page. That page needs more ideas that
have more descriptions on them (like the ones at the bottom). Not
only will this help potential students, but it will help us a lot when
we apply.

Aaron Meurer

>
> Here are some projects that interest me
>
> Framework design - we need a sufficiently general framework (this is hard
> and probably has to be half completed before GSoC time)
> Abstract Vector Spaces
> Tensor Math - Krastanov was talking about this and I think it's a great
> idea. There is a lot of good multilinear algebra out there that SymPy
> doesn't currently touch at all.
> General storage - Efficient NDArray classes (dense, mutable, sparse,
> functional, numpy, external programs) - views of NDArrays (transpose,
> slices).
> Theorem proving type system for
> tensors/matrices http://scicomp.stackexchange.com/questions/74/symbolic-software-packagbasic es-for-matrix-expressions
>
> I've dumped some thoughts on the following wiki page
> https://github.com/sympy/sympy/wiki/Linear-Algebra-Vision
>
> Comments or questions are welcome.
>
> -Matt
>

> --
> You received this message because you are subscribed to the Google Groups
> "sympy" group.
> To post to this group, send email to sy...@googlegroups.com.
> To unsubscribe from this group, send email to
> sympy+un...@googlegroups.com.
> For more options, visit this group at
> http://groups.google.com/group/sympy?hl=en.

Aaron Meurer

unread,
Feb 20, 2012, 10:00:18 PM2/20/12
to sy...@googlegroups.com
I updated the section on that page about ground types to be a little
clearer about what we want there. Many things will actually involve
improvements to Poly() to get things to work (for example, the
addition of a Frac() class).

Aaron Meurer

Alan Bromborsky

unread,
Feb 20, 2012, 10:38:24 PM2/20/12
to sy...@googlegroups.com
I am rewriting the geometric algebra module using noncommuting symbols
to represent base vectors and base multivectors and let general
multivectors be linear combinations of the base multivectors and let
sympy do the heavy
lifting with regard to simplification and other operations. The current
GA module was written before I understood what sympy could do and is a
real mess from the point of view of being an extensible module.

Someone should employ all the inherent abilities of sympy to write a
module for Tensors using Penrose's abstract index
http://en.wikipedia.org/wiki/Abstract_index_notation notation especially
allowing for simplifications using tensor symmetries and also allow for
symbolic differentiation.


Matthew Rocklin

unread,
Feb 20, 2012, 11:20:29 PM2/20/12
to sy...@googlegroups.com
@Aaron Thanks for the feedback, the KroneckerDelta anecdote is hilarious. I suspect that a lot of current work could be refactored with a good linear algebra module. There is even a current GSoC project discussing units that could make use of it. Lots of things that SymPy finds difficult are easily representable as vector spaces. Thanks also for updating the page, ground types is not something I know much about. I'll add to the projects page when I get some time.

@Alan I completely agree about abstract index notation. I have a small section of the wiki page on tensor indexing. There I pose that we might be able to use indexing as a way to define most operations. The tricky thing with this project is how to describe all tensor-ish operations with a single syntax that can be applied to a wide variety of applications. Can we create a single system which can be cleanly extended to include both Riemannian geometry and dense matrix computations? My current thought is that indexed expressions are a powerful-yet-general description that these representations all share. Any work done on indexed expressions (such as your idea of finding symmetries) would apply for all application projects. 

I really appreciate the input so far. I think that if only one person thinks about this project then it will likely end up as just-another-linear-algebra-module. Extra perspectives greatly reduce this risk. I would be interested in what people could see coming out of a general linear algebra system. What projects would this facilitate?


For more options, visit this group at
http://groups.google.com/group/sympy?hl=en.
I am rewriting the geometric algebra module using noncommuting symbols to represent base vectors and base multivectors and let general multivectors be linear combinations of the base multivectors  and let sympy do the heavy
lifting with regard to simplification and other operations.  The current GA module was written before I understood what sympy could do and is a real mess from the point of view of being an extensible module.

Someone should employ all the inherent abilities of sympy to write a module for Tensors using Penrose's abstract index http://en.wikipedia.org/wiki/Abstract_index_notation notation especially allowing for simplifications using tensor symmetries and also allow for symbolic differentiation.
--
You received this message because you are subscribed to the Google Groups "sympy" group.
To post to this group, send email to sy...@googlegroups.com.
To unsubscribe from this group, send email to sympy+unsubscribe@googlegroups.com.

Alan Bromborsky

unread,
Feb 21, 2012, 8:27:06 AM2/21/12
to sy...@googlegroups.com

For more options, visit this group at
http://groups.google.com/group/sympy?hl=en.
I am rewriting the geometric algebra module using noncommuting symbols to represent base vectors and base multivectors and let general multivectors be linear combinations of the base multivectors  and let sympy do the heavy
lifting with regard to simplification and other operations.  The current GA module was written before I understood what sympy could do and is a real mess from the point of view of being an extensible module.

Someone should employ all the inherent abilities of sympy to write a module for Tensors using Penrose's abstract index http://en.wikipedia.org/wiki/Abstract_index_notation notation especially allowing for simplifications using tensor symmetries and also allow for symbolic differentiation.



--
You received this message because you are subscribed to the Google Groups "sympy" group.
To post to this group, send email to sy...@googlegroups.com.
To unsubscribe from this group, send email to sympy+un...@googlegroups.com.

For more options, visit this group at http://groups.google.com/group/sympy?hl=en.

--
You received this message because you are subscribed to the Google Groups "sympy" group.
To post to this group, send email to sy...@googlegroups.com.
To unsubscribe from this group, send email to sympy+un...@googlegroups.com.

For more options, visit this group at http://groups.google.com/group/sympy?hl=en.
I have one suggestion for concrete tensors as opposed to abstract tensors.  Represent concrete tensors with numpy arrays stuffed with sympy scalars and let numpy do the heavy lifting of concrete tensor operations.

krastano...@gmail.com

unread,
Feb 22, 2012, 7:51:47 PM2/22/12
to sy...@googlegroups.com
On 20 February 2012 23:20, Matthew Rocklin <mroc...@gmail.com> wrote:
> @Aaron Thanks for the feedback, the KroneckerDelta anecdote is hilarious. I
> suspect that a lot of current work could be refactored with a good linear
> algebra module. There is even a current GSoC project discussing units that
> could make use of it. Lots of things that SymPy finds difficult are
> easily representable as vector spaces. Thanks also for updating the page,
> ground types is not something I know much about. I'll add to the projects
> page when I get some time.
>
> @Alan I completely agree about abstract index notation. I have a small
> section of the wiki page on tensor indexing. There I pose that we might be
> able to use indexing as a way to define most operations. The tricky thing
> with this project is how to describe all tensor-ish operations with a single
> syntax that can be applied to a wide variety of applications. Can we create
> a single system which can be cleanly extended to include both Riemannian
> geometry and dense matrix computations? My current thought is that indexed
> expressions are a powerful-yet-general description that these
> representations all share.
I do not think that mixing tensors and matrices too much is a good
idea. The main problem with the proposals I have read on the wiki is
that matrices and linear operators are not distinguished sufficiently
well. Many of the proposals make sense for linear operator, but not
for matrices.

More precisely:

1. Matrices are only a representation of a linear operator.

2. A basis must be chosen before creating this representation. The
information about the basis is not contained in the matrix
representation.

3. The information about the vector space is lost when taking this
representation. For example a linear operator from R^n to R^m (vectors
formed of real numbers) and from R[n-1] to R[m-1] (polynomials with
real coefficients) will both be matrices in M(R, n, m) (n-by-m real
matrices)

4. MatrixSymbol is not enough. We need LinearOperator that contains
VectorSpace (one for the domain of definition and one for the Image).
Then we can generate matrices and matrix symbols from the
LinearOperator.

5. Tensor is an object closer to LinearOperator than to Matrix. It
contains the information about the TangentBundle. Only when we add a
basis and the indices stop being abstract we get actual
multidimensional arrays, but this is not the interesting part in their
usage.

6. Using the tensor as a base for matrix is a bad idea because:
6.1 Abstract tensors do not depend on the basis (they are akin to
LinearOperators, not matrices, see 5 and 2)
6.2 Even when the indices stop being abstract you will use tensors for
stuff that is much different than the main usage of matrices. I won't
use "sparse" tensors, but I may use "sparse" matrices. I don't care
about the ground types of the tensor, because I use it as a symbol,
not as a container. I care about the base type of matrices, because
they are used in the core of high performance algorithms just like the
poly module.
6.3 TangentBundle and VectorSpace are not the same at all in
implementation terms. (I can back down on this if I see a good
proposal).

So I imagine something like:
1. LinearOperator, VectorSpace, Vector (_not_ a container, just a
symbol, no shape, it may be a function or a polynomial if the vector
space is over them)

2. Tensor, TangentBundle, AbstractTensor (with abstract indices)
(Those are _not_ containers, but symbols)

3. Indexed (the class currently in the tensor module. It seems like
the way forward for interoperability with numpy) This is a container
with some other smart operation defined on it. It does not contain
information about the basis that is used. It has shape.

4. Matrices, MatrixSymbol, ImmutableMatrix. (2D container that rests
completely separate from Indexed for performance reasons. It has
shape. Does not contain information about the basis that is used).

5. When deciding on a basis, LinearOperator can generate a Matrix and
Tensor can generate an Idexed.

6. Maybe creation of IndexedSymbol akin to MatrixSymbol, and
ImmutableIndexed (I am not sure that Indexed is actually mutable)

7. Maybe the creation of common base for IndexedSymbol and
MatrixSymbol, but not for Indexed and Matrix

8. Ground types are import for Matrix, not for the others. Container
types (list, numpy arrays, etc.) are important for Matrix and Indexed,
not the others.

I am definitely not seeing the whole picture, but I consider the
distinction between LinearOperator over a VectorSpace and Matrix over
some basis a very important one.

The symmetries treating code should be implemented in LinearOperator
and Tensor and not in Matrix and Indexed in my opinion. There should
be some way to check the Matrix class for symmetries.

Stefan

>>>>> sympy+un...@googlegroups.com.


>>>>> For more options, visit this group at
>>>>> http://groups.google.com/group/sympy?hl=en.
>>
>> I am rewriting the geometric algebra module using noncommuting symbols to
>> represent base vectors and base multivectors and let general multivectors be
>> linear combinations of the base multivectors  and let sympy do the heavy
>> lifting with regard to simplification and other operations.  The current
>> GA module was written before I understood what sympy could do and is a real
>> mess from the point of view of being an extensible module.
>>
>> Someone should employ all the inherent abilities of sympy to write a
>> module for Tensors using Penrose's abstract index
>> http://en.wikipedia.org/wiki/Abstract_index_notation notation especially
>> allowing for simplifications using tensor symmetries and also allow for
>> symbolic differentiation.
>>
>>
>>
>> --
>> You received this message because you are subscribed to the Google Groups
>> "sympy" group.
>> To post to this group, send email to sy...@googlegroups.com.
>> To unsubscribe from this group, send email to

>> sympy+un...@googlegroups.com.


>> For more options, visit this group at
>> http://groups.google.com/group/sympy?hl=en.
>>
>

> --
> You received this message because you are subscribed to the Google Groups
> "sympy" group.
> To post to this group, send email to sy...@googlegroups.com.
> To unsubscribe from this group, send email to

> sympy+un...@googlegroups.com.

Matthew Rocklin

unread,
Feb 22, 2012, 8:30:39 PM2/22/12
to sy...@googlegroups.com
Hi Krastaonv

Thanks for the feedback. You bring up a number of important issues. I'd like to respond to the main one

I completely agree that the solution needs to clearly separate abstract tensors (multi-linear functions) from indexed collections of components (ndarrays). I think the difficult part is to find a formulation which arranges all of the classes/ideas you mention above (LinearOperator, Matrix, Tensor, Indexed, VectorSpace, etc...) in a simple and elegant manner. 

krastano...@gmail.com

unread,
Feb 22, 2012, 8:31:14 PM2/22/12
to sy...@googlegroups.com
In the last mail I completely disregarded the difference between
Tensor and TensorField (and Vector and VectorField). For example the
matrix describing a homogeneous VectorField in the Cartesian basis
will not depend on coordinates, but it will depend on coordinates in
spherical basis (because the basis of the VectorSpace in each point of
the VectorField is different).

On 22 February 2012 19:51, krastano...@gmail.com

krastano...@gmail.com

unread,
Feb 22, 2012, 8:50:51 PM2/22/12
to sy...@googlegroups.com
On 22 February 2012 20:30, Matthew Rocklin <mroc...@gmail.com> wrote:
> Hi Krastaonv
>
> Thanks for the feedback. You bring up a number of important issues. I'd like
> to respond to the main one
>
> I completely agree that the solution needs to clearly separate abstract
> tensors (multi-linear functions) from indexed collections of components
> (ndarrays). I think the difficult part is to find a formulation which
> arranges all of the classes/ideas you mention above (LinearOperator, Matrix,
> Tensor, Indexed, VectorSpace, etc...) in a simple and elegant manner.
I will try to write a proposal about the object tree that can be used
(I do not know when I will have the time). Something containing the
likes of

VectorSpace, Vector, Covector, AlgebraicDual, ContinuousDual,
LinearOperator, BilinearForm, VectorBundle, VectorField, TensorProduct
(of spaces and of tensors), TangentBundle, CotangentBundle,
TensorBundle, Tensor, TensorField

Then from those one should be able to create Indexed and Matrix
instances (also MatrixSymbol and IndexedSymbol probably).

Maybe all this can be started in a new module, something like abstract_vectors?

I suppose that before starting a tensor gsoc project, those should be in place.

A substantial part of the quantum mechanics module can be simplified
if those exist.

Stefan

Aaron Meurer

unread,
Feb 23, 2012, 4:44:18 AM2/23/12
to sy...@googlegroups.com
I agree with these points 100%. We don't want to limit the power of
one thing by relating it too much with another. The base classes
should be based on how things work instead of so much what exact
mathematical objects they are intended to represent. Of course, if
that concept coincides with a mathematical concept, it can use that
name.

This is a general thing that happens with implementing symbolic
mathematics In a computer. It forces you to think of things as actual
living objects, which have differences that are superficial, don't
make sense, or are just ignored in the classical mathematical way of
looking at things.

For example, mathematically the concept of immutablity is meaningless.
But it's very important in a computer implementation of a concept.

Another issue is that concepts have to be objectified to be used. To
take an example from this discussion, mathematically, a linear
operator on a finite dinensional vector space can always be considered
as equivalent to some matrix. But if you want to use that idea in
implementation, you have to actually generate the matrix, which is not
actually doable in full generality (or at lead not easily). A better
way is to implement things without using matrices at all. You don't
need a matrix representation to get things like L(a*x1 + b*x2) =
a*L(x1) + b*L(x2), even though technically one exists.

This also has a side effect of making things more general. In this
case, you could also work with infinite dimensional operators, even on
non-separable spaces, which matrices or even more advanced basis
representations would not allow.

Aaron Meurer

On Feb 22, 2012, at 5:52 PM, "krastano...@gmail.com"

Akshay Srinivasan

unread,
Feb 28, 2012, 10:55:14 AM2/28/12
to sympy


On Feb 23, 6:50 am, "krastanov.ste...@gmail.com"
<krastanov.ste...@gmail.com> wrote:
> On 22 February 2012 20:30, Matthew Rocklin <mrock...@gmail.com> wrote:> Hi Krastaonv
>
> > Thanks for the feedback. You bring up a number of important issues. I'd like
> > to respond to the main one
>
> > I completely agree that the solution needs to clearly separate abstract
> > tensors (multi-linear functions) from indexed collections of components
> > (ndarrays). I think the difficult part is to find a formulation which
> > arranges all of the classes/ideas you mention above (LinearOperator, Matrix,
> > Tensor, Indexed, VectorSpace, etc...) in a simple and elegant manner.
>
> I will try to write a proposal about the object tree that can be used
> (I do not know when I will have the time). Something containing the
> likes of
>
> VectorSpace, Vector, Covector, AlgebraicDual, ContinuousDual,
> LinearOperator, BilinearForm, VectorBundle, VectorField, TensorProduct
> (of spaces and of tensors), TangentBundle, CotangentBundle,
> TensorBundle, Tensor, TensorField
>
> Then from those one should be able to create Indexed and Matrix
> instances (also MatrixSymbol and IndexedSymbol probably).
>
> Maybe all this can be started in a new module, something like abstract_vectors?
>
> I suppose that before starting a tensor gsoc project, those should be in place.
>
> A substantial part of the quantum mechanics module can be simplified
> if those exist.
>
> Stefan
>
>
>
>
>
>
>
>
>
> > On Wed, Feb 22, 2012 at 6:51 PM, krastanov.ste...@gmail.com
> > <krastanov.ste...@gmail.com> wrote:
> >> > On Mon, Feb 20, 2012 at 9:38 PM, Alan Bromborsky <abro...@verizon.net>
> >> > wrote:
>
> >> >> On 02/20/2012 10:00 PM, Aaron Meurer wrote:
>
> >> >>> I updated the section on that page about ground types to be a little
> >> >>> clearer about what we want there.  Many things will actually involve
> >> >>> improvements to Poly() to get things to work (for example, the
> >> >>> addition of a Frac() class).
>
> >> >>> Aaron Meurer
>
> >> >>> On Mon, Feb 20, 2012 at 7:40 PM, Aaron Meurer<asmeu...@gmail.com>
> >> >>>  wrote:
>
> >> >>>> On Mon, Feb 20, 2012 at 2:19 PM, Matthew Rocklin<mrock...@gmail.com>
> >> >>>> only will this help potential...
>
> read more »

This'd be an interesting read given that you're interested with
objects on abstract manifolds.
http://groups.csail.mit.edu/mac/users/wisdom/AIM-2005-003.pdf
Reply all
Reply to author
Forward
0 new messages