Re: norms and inner products for vectors of matrices

classic Classic list List threaded Threaded
10 messages Options
Reply | Threaded
Open this post in threaded view
|

Re: norms and inner products for vectors of matrices

Karl Meerbergen
Peter Gottschling wrote:

>Concerning the norms, I would like to limit the result type to real
>values as I have seen it in all mathematical definitions.  Does
>somebody sees a reason to vectors or matrices as results?
>
>The definitions of vector<vector> norms are straight forward:
>- norm_1 := sum (norm_1(x[i])
>- norm_2 := sqrt(sum(norm_2(x[i])^2))            // there might be more
>efficient ways to compute this
>- norm_inf := max(norm_inf(x[i])
>
>Best,
>Peter
>
>
>  
>

I am afraid that we cannot define norm_1, etc in this way. (I too made
this mistake in a previous mail) But we can define functions
generalized_norm_1 := sum (generalized_norm_1(x[i])
generalized_norm_2 :=sqrt(  sum (generalized_norm_2(x[i])^2 )
generalized_norm_inf := max (generalized_norm_inf(x[i])

norm_1, norm_2, norm_inf are well defined in linear algebra and these
definitions should be respected:
norm_1 := sum( abs(x[i]) )
etc.
These should only be used for vector and matrix objects with 'scalar'
value_type's.


Karl

_______________________________________________
glas mailing list
[hidden email]
http://lists.boost.org/mailman/listinfo.cgi/glas
Reply | Threaded
Open this post in threaded view
|

Re: norms and inner products for vectors of matrices

Peter Gottschling

On 13.01.2006, at 05:02, Karl Meerbergen wrote:

> Peter Gottschling wrote:
>
>> Concerning the norms, I would like to limit the result type to real
>> values as I have seen it in all mathematical definitions.  Does
>> somebody sees a reason to vectors or matrices as results?
>>
>> The definitions of vector<vector> norms are straight forward:
>> - norm_1 := sum (norm_1(x[i])
>> - norm_2 := sqrt(sum(norm_2(x[i])^2))            // there might be
>> more
>> efficient ways to compute this
>> - norm_inf := max(norm_inf(x[i])
>>
>> Best,
>> Peter
>>
>>
>>
>>
>
> I am afraid that we cannot define norm_1, etc in this way. (I too made
> this mistake in a previous mail) But we can define functions
> generalized_norm_1 := sum (generalized_norm_1(x[i])
> generalized_norm_2 :=sqrt(  sum (generalized_norm_2(x[i])^2 )
> generalized_norm_inf := max (generalized_norm_inf(x[i])
>
> norm_1, norm_2, norm_inf are well defined in linear algebra and these
> definitions should be respected:
> norm_1 := sum( abs(x[i]) )
> etc.
> These should only be used for vector and matrix objects with 'scalar'
> value_type's.
>
>
> Karl
>
Karl, I agree that it is not a good idea to call always the recursive
definition.

(Even if I think it is possible by defining norm_1 of scalars as abs
etc.  Anyway, it's ugly.)

The way I would like to define it is to have one function and dispatch
first between scalar value types and nested value types. For instance

template <typename Vector>
double norm_1(const Vector& v, type::scalar)
{
     return sum(abs(x[i]));       // well this is only pseudo-code
}

template <typename Vector>
double norm_1(const Vector& v, type::nested)
{
     return sum(norm_1(x[i]));
}

template <typename Vector>
double norm_1(const Vector& v)
{
     nesting_traits<typename Vector::value_type>::type nesting;
     return norm_1(v, nesting);
}

With this code we have the classical defintion if the value_type is
scalar and the recursive definition is the value_type itself is a
vector.  Note that this also handles higher levels of nesting than 2.

What remains is the question how we define norm_1(vector<matrix<T> >) ?
  I'm not convinced that my definition is perfect for this.  OTOH, I
wonder if there is a reasonable norm definition for vector<matrix> ????

May be we should discuss first what we consider the norm of a matrix.  
Matrices are typically used as linear operators and operator norms can
be defined in terms of vector norms norm_xx(op)=sup( norm_xx( op(v) ) /
norm_xx (v) ).  I vote for defining matrix norms as operator norms.  
OTOH, matrices model the concept of VectorSpace and it is
mathematically perfectly correct to apply the vector norm directly.  To
do this I suggest to use views for it, e.g.

vector_view<matrix<...> > vv(m);
alpha= norm_xx(vv);

One complication I see is that nested matrices has to be considered as
nested vectors, in other words the view must be applied recursively.

Coming back to the question what is norm_xx( vector<matrix> ) it would
be a mixture of vector and operator norms, which doesn't look
reasonable to me.  Probably we see it clearer if we know what usage
vector<matrix> have (if there is any).

Best,
Peter
------------
Peter Gottschling
Research Associate
Open Systems Laboratory
Indiana University
301i Lindley Hall
Bloomington, IN 47405
Tel.: +1 812 855-8898   Fax: +1 812 856 0853
http://www.osl.iu.edu/~pgottsch

_______________________________________________
glas mailing list
[hidden email]
http://lists.boost.org/mailman/listinfo.cgi/glas
Reply | Threaded
Open this post in threaded view
|

Re: norms and inner products for vectors of matrices

Karl Meerbergen
Peter Gottschling wrote:

>On 13.01.2006, at 05:02, Karl Meerbergen wrote:
>
>  
>
>>Peter Gottschling wrote:
>>
>>    
>>
>>>Concerning the norms, I would like to limit the result type to real
>>>values as I have seen it in all mathematical definitions.  Does
>>>somebody sees a reason to vectors or matrices as results?
>>>
>>>The definitions of vector<vector> norms are straight forward:
>>>- norm_1 := sum (norm_1(x[i])
>>>- norm_2 := sqrt(sum(norm_2(x[i])^2))            // there might be
>>>more
>>>efficient ways to compute this
>>>- norm_inf := max(norm_inf(x[i])
>>>
>>>Best,
>>>Peter
>>>
>>>
>>>
>>>
>>>      
>>>
>>I am afraid that we cannot define norm_1, etc in this way. (I too made
>>this mistake in a previous mail) But we can define functions
>>generalized_norm_1 := sum (generalized_norm_1(x[i])
>>generalized_norm_2 :=sqrt(  sum (generalized_norm_2(x[i])^2 )
>>generalized_norm_inf := max (generalized_norm_inf(x[i])
>>
>>norm_1, norm_2, norm_inf are well defined in linear algebra and these
>>definitions should be respected:
>>norm_1 := sum( abs(x[i]) )
>>etc.
>>These should only be used for vector and matrix objects with 'scalar'
>>value_type's.
>>
>>
>>Karl
>>
>>    
>>
>Karl, I agree that it is not a good idea to call always the recursive
>definition.
>
>(Even if I think it is possible by defining norm_1 of scalars as abs
>etc.  Anyway, it's ugly.)
>
>The way I would like to define it is to have one function and dispatch
>first between scalar value types and nested value types. For instance
>
>template <typename Vector>
>double norm_1(const Vector& v, type::scalar)
>{
>     return sum(abs(x[i]));       // well this is only pseudo-code
>}
>
>template <typename Vector>
>double norm_1(const Vector& v, type::nested)
>{
>     return sum(norm_1(x[i]));
>}
>
>template <typename Vector>
>double norm_1(const Vector& v)
>{
>     nesting_traits<typename Vector::value_type>::type nesting;
>     return norm_1(v, nesting);
>}
>

Agreed. Actually, the current implementation defines norm_1=abs for a
scalar. The nesting_traits is not required.

>
>With this code we have the classical defintion if the value_type is
>scalar and the recursive definition is the value_type itself is a
>vector.  Note that this also handles higher levels of nesting than 2.
>
Agreed. The current implementation uses a similar mechanism

>
>What remains is the question how we define norm_1(vector<matrix<T> >) ?
>  I'm not convinced that my definition is perfect for this.  OTOH, I
>wonder if there is a reasonable norm definition for vector<matrix> ????
>

sum(norm_1(v[i]) is a norm, there is no doubt about that. But there is
no correspondance with the algebraic matrix 1-norm.

>
>May be we should discuss first what we consider the norm of a matrix.  
>Matrices are typically used as linear operators and operator norms can
>be defined in terms of vector norms norm_xx(op)=sup( norm_xx( op(v) ) /
>norm_xx (v) ).  I vote for defining matrix norms as operator norms.  
>

Agreed. That is the way matrix norms are defined. We should stay as
close as possible to the operator norms.

>OTOH, matrices model the concept of VectorSpace and it is
>mathematically perfectly correct to apply the vector norm directly.  To
>do this I suggest to use views for it, e.g.
>
>vector_view<matrix<...> > vv(m);
>alpha= norm_xx(vv);
>
>One complication I see is that nested matrices has to be considered as
>nested vectors, in other words the view must be applied recursively.
>

I do not see a problem with that. The same holds for all functions, e.g.
trans, herm, conj, etc

>
>Coming back to the question what is norm_xx( vector<matrix> ) it would
>be a mixture of vector and operator norms, which doesn't look
>reasonable to me.  Probably we see it clearer if we know what usage
>vector<matrix> have (if there is any).
>

It depends what the user wants to do with the nested matrix. If there is
an operator behind, the user can create a norm for that. Unfortunately,
for such interpretations, GLAS cannot provides the functions, but GLAS
can support the creation of such functions in an easy way.

Conclusion;
* we keep the current implementation norm_1=sum(norm_1(v[i]) ?
* If the user wants a more appropriate norm for his application, he will
have to write it. As Toon mentioned, we cannot support all possible
interpretations of nested collections. (But we can provide tools to
develop appropriate functions). Recall that the collection attributes
are very powerful for that.


Karl


>
>Best,
>Peter
>------------
>Peter Gottschling
>Research Associate
>Open Systems Laboratory
>Indiana University
>301i Lindley Hall
>Bloomington, IN 47405
>Tel.: +1 812 855-8898   Fax: +1 812 856 0853
>http://www.osl.iu.edu/~pgottsch
>
>_______________________________________________
>glas mailing list
>[hidden email]
>http://lists.boost.org/mailman/listinfo.cgi/glas
>
>  
>


--
==============================================
Look at our unique training program and
Register on-line at http://www.fft.be/?id=35
----------------------------------------------
Karl Meerbergen
Free Field Technologies
16 place de l'Université
B-1348 Louvain-la-Neuve - BELGIUM
Company Phone:  +32 10 45 12 26
Company Fax:    +32 10 45 46 26
Mobile Phone:   +32 474 26 66 59
Home Phone:     +32 2 306 38 10
mailto:[hidden email]
http://www.fft.be
============================================
NEW ADDRESS FROM FEBRUARY 1st ONWARD:

Axis Park Louvain-la-Neuve
rue Emile Francqui, 1
B-1435 Mont-Saint Guibert - BELGIUM
Same telephone, same fax, same e-mail
============================================



_______________________________________________
glas mailing list
[hidden email]
http://lists.boost.org/mailman/listinfo.cgi/glas
Reply | Threaded
Open this post in threaded view
|

Re: norms and inner products for vectors of matrices

Karl Meerbergen
In reply to this post by Peter Gottschling
Hi,

Here is a question about nested vectors.

Suppose I have
vector< vector< double > > v ;
vector< double > w ;
double d ;

d*v : scalar times vector;
trans(w)*w: dot product of vector of vector with vector of vector: sum_i
trans(w[i])*w[i]
trans(v)*v: dot product of vector with vector: sum_i v[i]*v[i]

How should trans(v)*w be interpreted?
* result is a vector<double> which is the sum: sum_i v[i]*w[i] ?
* result is a vector<double> with element i being: trans(v)*w[i] ?

Best,

Karl

_______________________________________________
glas mailing list
[hidden email]
http://lists.boost.org/mailman/listinfo.cgi/glas
Reply | Threaded
Open this post in threaded view
|

Re: norms and inner products for vectors of matrices

Andrew Lumsdaine
I think the answer is: "it depends."

In a generic library, the operator * is defined between two concepts  
for which multiplication makes sense, i.e., a scalar times a scalar,  
a scalar times a member of a vector space, a linear operator times a  
member of a vector space.

So, if you have some particular types, the meaning of the operator *  
between those types depends on the manner in which those types model  
the concepts for which * is defined.

For instance, a vector<vector<double>> could be a row matrix, a  
column matrix, a diagonal matrix, a banded matrix, etc.  Each of  
these different interpretations of that type (different ways the type  
can model the concept) will have different interpretations of the  
multiply operator.

I think it is very important to keep the mathematical concepts and  
the concrete implementations distinct from each other.

On Jan 17, 2006, at 4:15 AM, Karl Meerbergen wrote:

> Hi,
>
> Here is a question about nested vectors.
>
> Suppose I have
> vector< vector< double > > v ;
> vector< double > w ;
> double d ;
>
> d*v : scalar times vector;
> trans(w)*w: dot product of vector of vector with vector of vector:  
> sum_i
> trans(w[i])*w[i]
> trans(v)*v: dot product of vector with vector: sum_i v[i]*v[i]
>
> How should trans(v)*w be interpreted?
> * result is a vector<double> which is the sum: sum_i v[i]*w[i] ?
> * result is a vector<double> with element i being: trans(v)*w[i] ?
>
> Best,
>
> Karl
>
> _______________________________________________
> glas mailing list
> [hidden email]
> http://lists.boost.org/mailman/listinfo.cgi/glas

_______________________________________________
glas mailing list
[hidden email]
http://lists.boost.org/mailman/listinfo.cgi/glas
Reply | Threaded
Open this post in threaded view
|

Re: norms and inner products for vectors of matrices

Matthias Troyer
I fully agree with Andrew and would like to propose that an operator*  
for vector< vector< double > > should not be defined by default. The  
users can always define the appropriate operator* for these types if  
they need them. The library should define an operator* only where it  
makes clear mathematical sense.

Matthias



On Jan 17, 2006, at 2:51 PM, Andrew Lumsdaine wrote:

> I think the answer is: "it depends."
>
> In a generic library, the operator * is defined between two concepts
> for which multiplication makes sense, i.e., a scalar times a scalar,
> a scalar times a member of a vector space, a linear operator times a
> member of a vector space.
>
> So, if you have some particular types, the meaning of the operator *
> between those types depends on the manner in which those types model
> the concepts for which * is defined.
>
> For instance, a vector<vector<double>> could be a row matrix, a
> column matrix, a diagonal matrix, a banded matrix, etc.  Each of
> these different interpretations of that type (different ways the type
> can model the concept) will have different interpretations of the
> multiply operator.
>
> I think it is very important to keep the mathematical concepts and
> the concrete implementations distinct from each other.
>
> On Jan 17, 2006, at 4:15 AM, Karl Meerbergen wrote:
>
>> Hi,
>>
>> Here is a question about nested vectors.
>>
>> Suppose I have
>> vector< vector< double > > v ;
>> vector< double > w ;
>> double d ;
>>
>> d*v : scalar times vector;
>> trans(w)*w: dot product of vector of vector with vector of vector:
>> sum_i
>> trans(w[i])*w[i]
>> trans(v)*v: dot product of vector with vector: sum_i v[i]*v[i]
>>
>> How should trans(v)*w be interpreted?
>> * result is a vector<double> which is the sum: sum_i v[i]*w[i] ?
>> * result is a vector<double> with element i being: trans(v)*w[i] ?
>>
>> Best,
>>
>> Karl
>>
>> _______________________________________________
>> glas mailing list
>> [hidden email]
>> http://lists.boost.org/mailman/listinfo.cgi/glas
>
> _______________________________________________
> glas mailing list
> [hidden email]
> http://lists.boost.org/mailman/listinfo.cgi/glas

_______________________________________________
glas mailing list
[hidden email]
http://lists.boost.org/mailman/listinfo.cgi/glas
Reply | Threaded
Open this post in threaded view
|

Re: norms and inner products for vectors of matrices

Russell Smiley
In reply to this post by Karl Meerbergen
While I can see that it's true that vector<vector<double> > is ambiguous
in terms of how operators might be applied to it (ie. Is this a
representation of a matrix, etc) , I don't quite understand why we have
to maintain that ambiguity in a numerical library.

What is wrong with having a concrete type 'matrix' (along with others)
such that you can still use vector<vector<double> > if necessary, but
resolve the ambiguity when needed?

For example:

vector<vector<double> > v;
// do some stuff with v...
Y = v*v;
// perhaps v*v is scalar multiple of elements of v;

// now assign v to the matrix so we can do some operations on it,
// making it clear that this is a diagonal matrix, somehow.
matrix<double> m(v, matrix::diagonal);
// etc, or

matrix<double, matrix::diagonal> m(v);

Z = m*m;
// now m*m (implicitly v*v) is matrix multiply of some form.

I assume there is some penalty for the assignment (perhaps not a
significant one if it's a reference copy) but at least the ambiguity is
resolved and syntactically it's fairly transparent/simple.

Russell.

> -----Original Message-----
> From: [hidden email]
> [mailto:[hidden email]] On Behalf Of Andrew Lumsdaine
> Sent: Tuesday, January 17, 2006 8:51 AM
> To: Generic Linear Algorithm Software
> Subject: Re: [glas] norms and inner products for vectors of matrices
>
>
> I think the answer is: "it depends."
>
> In a generic library, the operator * is defined between two concepts  
> for which multiplication makes sense, i.e., a scalar times a scalar,  
> a scalar times a member of a vector space, a linear operator times a  
> member of a vector space.
>
> So, if you have some particular types, the meaning of the operator *  
> between those types depends on the manner in which those types model  
> the concepts for which * is defined.
>
> For instance, a vector<vector<double>> could be a row matrix, a  
> column matrix, a diagonal matrix, a banded matrix, etc.  Each of  
> these different interpretations of that type (different ways
> the type  
> can model the concept) will have different interpretations of the  
> multiply operator.
>
> I think it is very important to keep the mathematical concepts and  
> the concrete implementations distinct from each other.
>

_______________________________________________
glas mailing list
[hidden email]
http://lists.boost.org/mailman/listinfo.cgi/glas
Reply | Threaded
Open this post in threaded view
|

Re: norms and inner products for vectors of matrices

Wolfgang Bangerth-3
In reply to this post by Karl Meerbergen

> How should trans(v)*w be interpreted?
> * result is a vector<double> which is the sum: sum_i v[i]*w[i] ?
> * result is a vector<double> with element i being: trans(v)*w[i] ?

I don't think there is a universal interpretation of what that expression
should mean, so simply don't implement it. If someone wants to attach a
particular meaning to it, they can always overload operators in their
application programs.

W.

-------------------------------------------------------------------------
Wolfgang Bangerth                email:            [hidden email]
                                 www: http://www.math.tamu.edu/~bangerth/

_______________________________________________
glas mailing list
[hidden email]
http://lists.boost.org/mailman/listinfo.cgi/glas
Reply | Threaded
Open this post in threaded view
|

Re: norms and inner products for vectors of matrices

Karl Meerbergen
Can we conclude that

1) 'Nested' collections are supported in containers and views (matrix
row, column, etc) only, but not in computational algorithms?

2) Those who want to use vector<vector<T>> etc in computations have to
implement their own.

3) GLAS should be relatively easy to extend, i.e. it should be possible
to write
glas::norm_1(v), where norm_1 is dispatched to a user defined function.


Karl



Wolfgang Bangerth wrote:

>>How should trans(v)*w be interpreted?
>>* result is a vector<double> which is the sum: sum_i v[i]*w[i] ?
>>* result is a vector<double> with element i being: trans(v)*w[i] ?
>>    
>>
>
>I don't think there is a universal interpretation of what that expression
>should mean, so simply don't implement it. If someone wants to attach a
>particular meaning to it, they can always overload operators in their
>application programs.
>
>W.
>
>-------------------------------------------------------------------------
>Wolfgang Bangerth                email:            [hidden email]
>                                 www: http://www.math.tamu.edu/~bangerth/
>
>_______________________________________________
>glas mailing list
>[hidden email]
>http://lists.boost.org/mailman/listinfo.cgi/glas
>
>  
>



_______________________________________________
glas mailing list
[hidden email]
http://lists.boost.org/mailman/listinfo.cgi/glas
Reply | Threaded
Open this post in threaded view
|

Re: norms and inner products for vectors of matrices

Matthias Troyer

On Jan 18, 2006, at 1:57 PM, Karl Meerbergen wrote:

> Can we conclude that
>
> 1) 'Nested' collections are supported in containers and views (matrix
> row, column, etc) only, but not in computational algorithms?
>
> 2) Those who want to use vector<vector<T>> etc in computations have to
> implement their own.

It should not be necessary to write your own algorithm. It should be  
enough to make sure that your data type models the required concepts,  
i.e. that it satisfies all the requirements. This might mean  
providing an operator* or a norm_1 function.


>
> 3) GLAS should be relatively easy to extend, i.e. it should be  
> possible
> to write
> glas::norm_1(v), where norm_1 is dispatched to a user defined  
> function.

Of course. Isn't that the whole idea behind the project?

Matthias

_______________________________________________
glas mailing list
[hidden email]
http://lists.boost.org/mailman/listinfo.cgi/glas