ufl team mailing list archive

["Martin Sandve Alnæs" <martinal@xxxxxxxxx>]

2008/4/7, Kent-Andre Mardal <kent-and@xxxxxxxxx>:
>
>  Difficult subject, in particular when it comes to matrices.

I'd say it's most difficult when it comes to _vectors_, not matrices.
Actually, we're not dealing with vectors and matrices, we're dealing
with tensors. A rank 1 tensor is not the same as a vector in linear
algebra, there is no concept of rows or columns.

>  I think we should let * and dot be the same, I mean dot is just a word
>  for *, right ?

Not really, that's the whole point. Did you look at the examples? "*"
means multiplication, but what does multiplication mean for tensors of
rank > 0? It can be interpreted differently, and that's why I'm asking
here.

>  The meaning of *should be the same as in Matlab, i.e.
>  A*B is not the contraction, but the matrix-matrix product.

That's one of the trivial examples in my previous email. A*B and A*v
work out fine, but v*A and u*v would be invalid products in
Matlab/linear algebra notation. If you write (with tensors) "u v", it
(often, usually, always?) means outer(u, v). Translating from linear
algebra notation to UFL tensor notation, "u v" is invalid, "u^T v" is
dot(u,v) and "u v^T" is outer(u, v).

"Paper notation" contains a lot of special cases that doesn't
generalize, and would lead to a tangled mess of code. That's why we
must make our own precise definitions, that are always well defined
and translates well to code. And since we're dealing with tensor
algebra, we must use tensor algebra as our basis and not linear
algebra where they differ.

-- 
Martin


>  Kent
>
>
>  man, 07.04.2008 kl. 11.51 +0200, skrev Martin Sandve Alnæs:
>
> > I'd like to hear opinions of how * should behave in different contexts.
>  >
>  > To begin with, we will of course agree on the meaning of * between
>  > scalar expressions, which may or may not have free indices. Repeated
>  > free indices will be interpreted as implicit summation. Single free
>  > indices are combined to define a larger set of free indices. This is
>  > already implemented. (Some of the implicit summation is not, but I
>  > have it on paper).
>  >
>  > We may perhaps also agree that multiplying a matrix (rank 2 tensor)
>  > with a vector (rank 1 tensor) means contraction over the last index of
>  > the matrix and the index of the vector, which is equivalent with the
>  > dot product. This will allow f.ex. Poissons equation with a conduction
>  > tensor (M Du : Dv) to be written like inner(M*grad(u), grad(v)). Ok?
>  >
>  > But if so, is this to be a special case, or do we simply define * to be dot?
>  > This will eventually lead to (let u,v be vectors, A,B,C matrices, a scalar):
>  >
>  > Very natural notation:
>  > u = A*v = dot(A, v) = "A v" in linear algebra notation
>  > C = A*B = dot(A, B) = "A B" in linear algebra notation
>  >
>  > Feels natural to me:
>  > u = v*A = dot(v, A) = "v^T A" in linear algebra notation
>  > a = u*A*v = (u*A)*v = u*(A*v)
>  >   = dot(u, A*v) = dot(u*A, v)
>  >   = dot(dot(u, A), v) = dot(u, dot(A, v)) = "u^T A v" in linear algebra notation
>  >
>  > However, keeping * and dot separate may make things simpler in other
>  > ways in the implementation. And easy code understanding isn't always
>  > the same as shorter expressions. So what do you think?
>  >
>  > Note that dot(f, g) is defined differently in UFL than in FFC:
>  >     dot(f, g) == sum(f[...,k] * g[k, ...] for k in range(dim))
>  > where ... means the obvious for a tensor of rank > 1.
>  > Using ... to "skip indices" is also implemented.
>  > The contraction A:B is instead implemented as inner(A, B).
>  >
>  > Another example, Navier-Stokes with explicit convection term (freely
>  > from memory):
>  > a = u - w + dt * ( w * grad(u) * v - nu*inner(grad(u), grad(v)) +
>  > p*div(v) - f*v )
>  > vs.
>  > a = u - w + dt * ( dot(dot(w, grad(u)), v) - nu*inner(grad(u),
>  > grad(v)) + p*div(v) - dot(f, v) )
>  >
>  > --
>  > Martin
>
> > _______________________________________________
>  > UFL-dev mailing list
>  > UFL-dev@xxxxxxxxxx
>  > http://fenics.org/mailman/listinfo/ufl-dev
>
>