lapackpp-devel — Discussion of Lapack++ development

Re: [Lapackpp-devel] Error configuring with externally compiled Blas & Lapack

[Christian S. <sti...@tu...>]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

(copying to list as well)

Sonia Singhal schrieb:
> Ya I just saw that. HEre is what it reads.
> ----------------------------------------------------------------------
> configure:22164: checking for sgemm_ in
> /scratch2/sonia/packages/mkl/lib/32/libmkl_ia32.a
> configure:22197: gcc -o conftest -g -O2   conftest.c
> /scratch2/sonia/packages/mkl/lib/32/libmkl_ia32.a
> -L/afs/ece.cmu.edu/support/gcc/4.1.1/4.1.1-1/i386_suse93/image/usr/local/bin/../lib/gcc/i686-pc-linux-gnu/4.1.1
>
> -L/afs/ece.cmu.edu/support/gcc/4.1.1/4.1.1-1/i386_suse93/image/usr/local/bin/../lib/gcc
>
> -L/usr/local/lib/gcc/i686-pc-linux-gnu/4.1.1
> -L/afs/ece.cmu.edu/support/gcc/4.1.1/4.1.1-1/i386_suse93/image/usr/local/bin/../lib/gcc/i686-pc-linux-gnu/4.1.1/../../..
>
> -L/usr/local/lib/gcc/i686-pc-linux-gnu/4.1.1/../../.. -lgfortranbegin
> -lgfortran -lm -lgcc_s >&5
> /scratch2/sonia/packages/mkl/lib/32/libmkl_ia32.a(def_sgemm_omp.o): In
> function `mkl_blas_def_sgemm':
> tmp_sgemm_omp.c:(.text+0x1e): undefined reference to
> `__kmpc_global_thread_num'
> tmp_sgemm_omp.c:(.text+0xe0): undefined reference to `omp_in_parallel'
> ...

Then you're in (some) trouble. This means that if a program links to
MKL's libmkl_ia32.a, it *must* also link to various other libraries that
contain these functions. In the documentation of Intel's MKL library it
should say something along these lines:

"If a program wants to use the functions from MKL blas, it must link to
the libraris libfoo, libbla, and libpoo",

In which case you should write the full linker arguments into the
- --with-blas or --with-lapack parameters, like so:

     ./configure --with-blas="-L/my/mkl/dir/lib/32 -lmkl_ia32 -lfoo
- -lbla -lpoo"

(If you specify a library to the linker with the "-l" [ell] argument,
you say the library name without leading "lib" and without trailing
".so" or ".a".)

Keep adding more libraries until the full error message in config.log
doesn't complain anymore about "undefined reference" in any function of
the installed library.

Regards,

Christian
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.2.1 (MingW32)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iQCVAwUBRbjHEGXAi+BfhivFAQIrqwP/Tu+IxCKC9l63s8yRy1W+aCNwNDp8mblX
wxcYFiSFr5EjeyYNjsG3nAnDmxQ9soUynl81vBusbnP/FBlWfT1O/o0SfGIppmyn
B9whdDI5/Fhh1wFOpn4tB2TX004gGKP6I97UvjY1XrEvKWeI4guGzjnFcpoDK1FR
d6lUwhRNsLQ=
=jMjP
-----END PGP SIGNATURE-----

Re: [Lapackpp-devel] Error configuring with externally compiled Blas & Lapack

[Christian S. <sti...@tu...>]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

Sonia Singhal schrieb:
> I gave the full path, and the config.log looks like this
> 
> -----------------------------------------------------------------
> configure:22227: result: no
> configure:22237: checking for sgemm_
> /tmp/ccUC6VD8.o: In function `main':
> /ux0/soniasin/packages/lapackpp-2.5.0/conftest.c:81: undefined
> reference to `sgemm_'
> /tmp/ccUC6VD8.o:(.data+0x0): undefined reference to `sgemm_'
> collect2: ld returned 1 exit status
> configure:22308: $? = 1
> configure: failed program was:
> -------------------------------------------------------------------------------------

There are several checks for sgemm_. The first one doesn't specify any
blas library, so it usually fails on every system. The second check
("checking for sgemm_ in") then uses the blas library and (hopefully)
succeeds.

Christian

> Its not linking the blas lib that I have given in the config command line.
> ./configure --with-blas=/scratch2/sonia/packages/mkl/lib/32/libmkl_ia32.a
> --with-lapack=/scratch2/sonia/packages/mkl/lib/32/libmkl_lapack.a


-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.2.1 (MingW32)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iQCVAwUBRbjBbWXAi+BfhivFAQKy9QP8DSttKfpttsfFYYLFU9ZKcre0RAXtCvcJ
XVQ7b+pGUzpft8gk65ld8ePAegTakONO5PAGk3hh3PWvRc68K6MetqMfKGptJs5B
OPa0AA6gxOI5WhQIveRdt+CjG70IpLlNPh61iDZ81T4mZNsHrnXtJwQALqt4IwgG
F4JDWpAayt4=
=+TE8
-----END PGP SIGNATURE-----

Re: [Lapackpp-devel] Error configuring with externally compiled Blas & Lapack

[Sonia S. <son...@gm...>]

Hi,
I gave the full path, and the config.log looks like this

-----------------------------------------------------------------

| }
configure:22227: result: no
configure:22237: checking for sgemm_
configure:22302: gcc -o conftest -g -O2   conftest.c
-L/afs/ece.cmu.edu/support/gcc/4.1.1/4.1.1-1/i386_suse93/image/usr/local/bin/../lib/gcc/i686-pc-linux-gnu/4.1.1
-L/afs/ece.cmu.edu/support/gcc/4.1.1/4.1.1-1/i386_suse93/image/usr/local/bin/../lib/gcc
-L/usr/local/lib/gcc/i686-pc-linux-gnu/4.1.1
-L/afs/ece.cmu.edu/support/gcc/4.1.1/4.1.1-1/i386_suse93/image/usr/local/bin/../lib/gcc/i686-pc-linux-gnu/4.1.1/../../..
-L/usr/local/lib/gcc/i686-pc-linux-gnu/4.1.1/../../.. -lgfortranbegin
-lgfortran -lm -lgcc_s >&5
/tmp/ccUC6VD8.o: In function `main':
/ux0/soniasin/packages/lapackpp-2.5.0/conftest.c:81: undefined
reference to `sgemm_'
/tmp/ccUC6VD8.o:(.data+0x0): undefined reference to `sgemm_'
collect2: ld returned 1 exit status
configure:22308: $? = 1
configure: failed program was:
-------------------------------------------------------------------------------------

Its not linking the blas lib that I have given in the config command line.
./configure --with-blas=/scratch2/sonia/packages/mkl/lib/32/libmkl_ia32.a
--with-lapack=/scratch2/sonia/packages/mkl/lib/32/libmkl_lapack.a

~ Sonia

On 1/25/07, Christian Stimming <sti...@tu...> wrote:
> -----BEGIN PGP SIGNED MESSAGE-----
> Hash: SHA1
>
> Hi,
>
> Sonia Singhal schrieb:
> > I am trying to compile lapackpp-2.5.0 on Suse Linux.
> > I am using Intel's MKL library containing Blas and Lapack libraries.
> > When I run configure giving the path to these libraries, it gives the
> > error of not finding functions in the blas/lapack libraries.
> > ---------------------------------------------------------------------------------------------------------
> >> ./configure --with-blas=~/pacakages/mkl/lib/32/libmkl_ia32.a --with-lapack=~/pacakages/mkl/lib/32/libmkl_lapack.a
> > .
> > .
> > checking for dummy main to link with Fortran libraries... none
> > checking for Fortran name-mangling scheme... lower case, underscore,
> > no extra underscore
> > checking for sgemm_ in ~/pacakages/mkl/lib/32/libmkl_ia32.a... no
>
> This is the important error in the configure result -- according to the
> "nm" output below, this check should have succeeded instead:
>
> >> nm ~/pacakages/mkl/lib/32/libmkl_ia32.a | grep sgemm_
> > ...
> > _sgemm_fb.o:
> > 00000000 T sgemm_
>
> Can you open the file "config.log" that was created during configure,
> and look for the exact (compiler/linker) error message on the failed
> test? To find this place, search for the string "checking for sgemm_" in
> that file. Maybe the compiler doesn't like the tilde ("~"), in which
> case you should specify a complete absolute path name instead.
>
> Regards,
>
> Christian
> -----BEGIN PGP SIGNATURE-----
> Version: GnuPG v1.4.2.1 (MingW32)
> Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org
>
> iQCVAwUBRbi3mmXAi+BfhivFAQL5AwQAhD9+flxjEMyfZj3y42F248uDy94nVuD3
> r8y+sKzhYRWizuIGvHFeaibqHoZlE1gMTnT8K9V50nwYllc8NhU028Zm+K3z1yGM
> JTtNMt6IKCOBwQveOKQH1Rktxy6f1LVOCxyvGcQ/dvuDHqRc12h9McDKNNNZLUjh
> bPlHY0+KVlk=
> =6x2V
> -----END PGP SIGNATURE-----
>

Re: [Lapackpp-devel] Error configuring with externally compiled Blas & Lapack

[Christian S. <sti...@tu...>]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

Hi,

Sonia Singhal schrieb:
> I am trying to compile lapackpp-2.5.0 on Suse Linux.
> I am using Intel's MKL library containing Blas and Lapack libraries.
> When I run configure giving the path to these libraries, it gives the
> error of not finding functions in the blas/lapack libraries.
> ---------------------------------------------------------------------------------------------------------
>> ./configure --with-blas=~/pacakages/mkl/lib/32/libmkl_ia32.a --with-lapack=~/pacakages/mkl/lib/32/libmkl_lapack.a
> .
> .
> checking for dummy main to link with Fortran libraries... none
> checking for Fortran name-mangling scheme... lower case, underscore,
> no extra underscore
> checking for sgemm_ in ~/pacakages/mkl/lib/32/libmkl_ia32.a... no

This is the important error in the configure result -- according to the
"nm" output below, this check should have succeeded instead:

>> nm ~/pacakages/mkl/lib/32/libmkl_ia32.a | grep sgemm_
> ...
> _sgemm_fb.o:
> 00000000 T sgemm_

Can you open the file "config.log" that was created during configure,
and look for the exact (compiler/linker) error message on the failed
test? To find this place, search for the string "checking for sgemm_" in
that file. Maybe the compiler doesn't like the tilde ("~"), in which
case you should specify a complete absolute path name instead.

Regards,

Christian
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.2.1 (MingW32)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iQCVAwUBRbi3mmXAi+BfhivFAQL5AwQAhD9+flxjEMyfZj3y42F248uDy94nVuD3
r8y+sKzhYRWizuIGvHFeaibqHoZlE1gMTnT8K9V50nwYllc8NhU028Zm+K3z1yGM
JTtNMt6IKCOBwQveOKQH1Rktxy6f1LVOCxyvGcQ/dvuDHqRc12h9McDKNNNZLUjh
bPlHY0+KVlk=
=6x2V
-----END PGP SIGNATURE-----

[Lapackpp-devel] Error configuring with externally compiled Blas & Lapack

[Sonia S. <son...@gm...>]

Hi,

I am trying to compile lapackpp-2.5.0 on Suse Linux.
I am using Intel's MKL library containing Blas and Lapack libraries.
When I run configure giving the path to these libraries, it gives the
error of not finding functions in the blas/lapack libraries.
---------------------------------------------------------------------------------------------------------
> ./configure --with-blas=~/pacakages/mkl/lib/32/libmkl_ia32.a --with-lapack=~/pacakages/mkl/lib/32/libmkl_lapack.a
.
.
.
.
checking for dummy main to link with Fortran libraries... none
checking for Fortran name-mangling scheme... lower case, underscore,
no extra underscore
checking for sgemm_ in ~/pacakages/mkl/lib/32/libmkl_ia32.a... no
checking for sgemm_... no
checking if ATLAS should be used... yes
checking for ATLAS library... checking for f77blas library...
configure: WARNING: not found
checking for cblas library... configure: WARNING: not found
checking whether ATLAS is usable... configure: WARNING: no
checking for sgemm_ in -lblas... no
checking for sgemm_ in -lcxml... no
checking for sgemm_ in -ldxml... no
checking for sgemm_ in -lscs... no
checking for sgemm_ in -lcomplib.sgimath... no
checking for sgemm_ in -lblas... (cached) no
checking for sgemm_ in -lblas... (cached) no
checking for sgemm_ in -lblas32... no
checking for cheev_ in ~/pacakages/mkl/lib/32/libmkl_lapack.a... noblas
configure: error: Blas/Lapack was not found.
*** This means Lapack++ and matrix support cannot be compiled.
*** This makes this library unusable. Please get blas and lapack
*** installed. If you do have these installed, use the options
*** --with-blas=<libname> or --with-lapack=<libname> and/or set
*** the env variable LDFLAGS to include the appropriate linker
*** flags.
---------------------------------------------------------------------------------
But when I nm the libraries that I am linking to, I can find the
functions that configure was looking for.
--------------------------------------------------------------------------------
> nm ~/pacakages/mkl/lib/32/libmkl_ia32.a | grep sgemm_
def_sgemm_omp.o:
00000004 b ___kmpv_zeromkl_blas_def_sgemm_0
00000000 b ___kmpv_zeromkl_blas_def_sgemm_1
0000051a t L_mkl_blas_def_sgemm_276__par_loop0
0000081c t L_mkl_blas_def_sgemm_296__par_loop1
def_xsgemm_1.o:
00000000 T mkl_blas_def_xsgemm_1
         U mkl_blas_def_xsgemm_1
p3_sgemm_omp.o:
00000004 b ___kmpv_zeromkl_blas_p3_sgemm_0
00000000 b ___kmpv_zeromkl_blas_p3_sgemm_1
0000051a t L_mkl_blas_p3_sgemm_276__par_loop0
00000801 t L_mkl_blas_p3_sgemm_296__par_loop1
p4_sgemm_omp.o:
00000004 b ___kmpv_zeromkl_blas_p4_sgemm_0
00000000 b ___kmpv_zeromkl_blas_p4_sgemm_1
0000051a t L_mkl_blas_p4_sgemm_276__par_loop0
00000801 t L_mkl_blas_p4_sgemm_296__par_loop1
p4p_sgemm_omp.o:
00000004 b ___kmpv_zeromkl_blas_p4p_sgemm_0
00000000 b ___kmpv_zeromkl_blas_p4p_sgemm_1
0000051a t L_mkl_blas_p4p_sgemm_276__par_loop0
00000801 t L_mkl_blas_p4p_sgemm_296__par_loop1
p4m_sgemm_omp.o:
00000004 b ___kmpv_zeromkl_blas_p4m_sgemm_0
00000000 b ___kmpv_zeromkl_blas_p4m_sgemm_1
0000051a t L_mkl_blas_p4m_sgemm_276__par_loop0
00000801 t L_mkl_blas_p4m_sgemm_296__par_loop1
_sgemm_fb.o:
00000000 T sgemm_
-------------------------------------------------------
What could be the problem here ?

Thanks,
~ sonia

[Lapackpp-devel] Problem using lapackpp & an optimization package

[Sonia S. <son...@gm...>]

Hi,

I am doing development in Microsoft Visual Studio 2003. I have been using
lapackpp-2.5.0 for solving linear system of equations and doing SVD.
Recently I started using an open source optimization package IPOpt to solve
an optimization problem. This library uses certain routines in BLAS, LAPACK,
HSL.

The problem is that my code crashes in the lapackpp's call to LaSVD_IP(). It
crashes giving a "floating point underflow" error. This doesnt happen if I
remove ipopt's libs and code. I havnt really written any code for calling
this optimization library, just instantiated its class just to see link the
ipopt infrastructure into my existing code. I have used IPOpt earlier but
then wasnt using lapackpp and ipopt worked fine. I have posted this question
on the ipopt forum too. I am posting here too just in case anyone has any
idea of what the problem could be.

Thanks,
~ sonia

Re: [Lapackpp-devel] Problem finding Null space using LaSVD_IP()

[Christian S. <sti...@tu...>]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

Good to hear. Feel free to ask more questions if something appears to be
wrong. By the way you should be able to write std::cout<<U_; as well.

Regarding the mailing list address: This is a "mailing list", not a
"forum". You have been sending your emails to the address
"lapackpp-devel-owner (at) lists.sourceforge.net" which is wrong; the
correct email address to send your messages to is
lap...@li... (no -owner)! Thanks.

Christian

Sonia Singhal schrieb:
> hi,
> 
> I think I caught my mistake. I didnt make a copy of the matrix before I
> checked for the homogenous solution. :(
> ~ sonia
> 
> On 12/14/06, *Sonia Singhal * <son...@gm...
> <mailto:son...@gm...>> wrote:
> 
>     Hi,
> 
>     Thanks for the reply.
>     As for posting it on the forum. I had sent the mail to the forum, I
>     have no clue how you have received the mail. I initially emailed it
>     without joining the list and got an error message, maybe thats the
>     reason you got it.
> 
>     As for the problem, even I am not familiar with the algorithms used
>     for svd and am not sure if the null spaces would be the same, but I
>     assumed that since MATLAB also uses LAPACK the answers would be the
>     same.
>     Even if the answers are not the same and I assume that the NULL
>     spaces are returned, then multiplying the original matrix A with a
>     row ( or column ) of VT ( or V ) beyond the rank should be zero. I
>     did this and I got a very high value when I used Norm_Inf(). Whereas
>     in MATLAB I get a very small value ( 10^-14).
> 
>     Also right now I was making small examples ( 3x4, 3x5 ) and what I
>     found was for some examples I am getting right answers and for some
>     they are wrong. By right and wrong I mean comparing both with matlab
>     and checking if A*VT(r,:) is zero.
> 
>     For example :
>     This code / exmaple gives the correct answers
> 
>     //------------------------------------------------------------------------------------------------------
>             // Trying SVD
>             LaGenMatDouble ForSVD(3,5);
>             ForSVD(0,0) = 1;
>             ForSVD(0,1) = 2;
>             ForSVD(0,2) = 3;
>             ForSVD(0,3) = 4;
>             ForSVD(0,4) = 2;
>             ForSVD(1,0) = 2;
>             ForSVD(1,1) = 3;
>             ForSVD(1,2) = 9;
>             ForSVD(1,3) = 2;
>             ForSVD(1,4) = 0;
>             ForSVD(2,0) = 4;
>             ForSVD(2,1) = 5;
>             ForSVD(2,2) = 8;
>             ForSVD(2,3) = 7;
>             ForSVD(2,4) = 4;
>             LaVectorDouble S_(3);
>             LaGenMatDouble U_(3,3), VT_(5,5);
>             LaSVD_IP(ForSVD, S_, U_, VT_);
>             LaVectorDouble x_(5);
>            
>             for(i=0; i<5; ++i){
>                 x_(i) = VT_(4,i);
>             }
>             std::cout << "\tNorm_Inf(A*x)" << Norm_Inf(ForSVD*x_) <<
>     std::endl;
>             std::cout << "U:" << std::endl;
>             for(int i=0; i<3; ++i){
>                 for(int j=0; j<3; ++j){
>                     std::cout << U_(i,j) << " ";
>                 }
>                 std::cout << std::endl;
>             }       
>             std::cout << "V:" << std::endl;
>             for(int i=0; i<5; ++i){
>                 for(int j=0; j<5; ++j){
>                     std::cout << VT_(i,j) << " ";
>                 }
>                 std::cout << std::endl;
>             }
> 
>     //
>     -----------------------------------------------------------------------------------
>     and this doesnt.
>     //
>     -----------------------------------------------------------------------------------
> 
>             // Trying SVD
>             LaGenMatDouble ForSVD(3,5);
>             ForSVD(0,0) = 1;
>             ForSVD(0,1) = 2;
>             ForSVD(0,2) = 3;
>             ForSVD(0,3) = 4;
>             ForSVD(1,0) = 2;
>             ForSVD(1,1) = 3;
>             ForSVD(1,2) = 9;
>             ForSVD(1,3) = 2;
>             ForSVD(2,0) = 4;
>             ForSVD(2,1) = 5;
>             ForSVD(2,2) = 8;
>             ForSVD(2,3) = 7;
>             LaVectorDouble S_(3);
>             LaGenMatDouble U_(3,3), VT_(4,4);
>             LaSVD_IP(ForSVD, S_, U_, VT_);
>             LaVectorDouble x_(5);
>            
>             for(i=0; i<5; ++i){
>                 x_(i) = VT_(4,i);
>             }
>             std::cout << "\tNorm_Inf(A*x)" << Norm_Inf(ForSVD*x_) <<
>     std::endl;
>             std::cout << "U:" << std::endl;
>             for(int i=0; i<3; ++i){
>                 for(int j=0; j<3; ++j){
>                     std::cout << U_(i,j) << " ";
>                 }
>                 std::cout << std::endl;
>             }       
>             std::cout << "V:" << std::endl;
>             for(int i=0; i<4; ++i){
>                 for(int j=0; j<4; ++j){
>                     std::cout << VT_(i,j) << " ";
>                 }
>                 std::cout << std::endl;
>             }
>     //
>     ------------------------------------------------------------------------------------------------
>     Thanks for your time,
>     ~ Sonia

-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.2.1 (MingW32)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iQCVAwUBRYJoP2XAi+BfhivFAQKzCQP7BvwUeKqGlbLSYrh2rVbWoy7yAFjYD5AT
cBVEPecyYR1sX5ZN1Dj3CPRquGB1GcNexlnQfPJ9Tnawc4ZY4djjZa7BtOxXFzJq
ZEfXlXbRiR1emMi/ikW5DLvvWPrkUGz9uCeA1MTIppL6xtcyd4Lui9PygAL+rrhq
u19NPB9Le8k=
=aaop
-----END PGP SIGNATURE-----

Re: [Lapackpp-devel] =?GB2312?B?SGxlcKOh?=

[Christian S. <sti...@tu...>]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

(Please always write to the mailing list lapackpp-devel, never to
individual developers!)

settingyin schrieb:
> Stimming,
>=20
> Hello! I=A1=AFm using LAPACK++,but something comfused me.
>
> I want to get  A  .*  B=20
>=20
> A=3D 1  2      =20
> 3         4 =20
> B=3D3   2
> 4         4
> A  .*  B=3D 3    4
>          12  16
>=20
> But I can=A1=AFt get it direcly.

I don't fully understand your question, but if you mean you want to have
an "element-wise multiplication", then indeed this is *not* available in
lapackpp. The reason is that lapackpp is merely a wrapper to the
functionality that is offered by the numerical libraries LAPACK and
BLAS, and those don't have this operation (probably because in numerics
this operation is really never used). You will have to program it yoursel=
f.

Christian
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.2.1 (MingW32)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iQCVAwUBRYJnWmXAi+BfhivFAQLFyAP/Qb/m7jragdiss7pyLUyd/NzDjt1Ib2XZ
/MdjU0p1BQfML9ZLIN1GX2W0C7uL8PiDvoANeIUiNPp6zqaqfw3s/KJLjqJ/lDvc
tUSM2LtLgmbuH8X1yWSkaW/7SpSoonQaYdXV4cyK+7u/Akk1eOzF05nE00ukfPc7
vWCWmowRJ+E=3D
=3DAkCa
-----END PGP SIGNATURE-----

Re: [Lapackpp-devel] Problem finding Null space using LaSVD_IP()

[Christian S. <sti...@tu...>]

(Please send your message to "lapackpp-devel", not -owner. Thanks.)

Am Donnerstag, 14. Dezember 2006 21:18 schrieb Sonia Singhal:
> Hi,
> I am using the LAPACKPP on Windows with Visual Studio. I downloaded it a
> couple of days back and am using the 2.5.0 version.
>
> I am using LaSVD_IP() to find the null space of a matrix. My matrix is with
> m<n. I was getting wrong answers.
> I dumped the U, S, VT matrices generated by the above call and checked it
> with Matlab results.
> Apparently, the U and V( orVT) matrix columns are fine upto the rank of the
> matrix. Beyond that they are all messed up, hence my null space is
> incorrect.

I'm not so sure you are expecting the right thing here. As you said, the 
singular vectors for the non-null singular values are as expected. However, 
the rest of the U and VT matrices, as I understand SVD, are quite ambiguous, 
given that they only have to span the null space. So as long as these vectors 
span that space, they can be totally different. As far as I understand it, 
you would have to check whether the spanned space of all these vectors from 
lapackpp is identical to the space spanned by the matlab vectors (I don't 
have a good check for that at hand, but surely one can think of something 
appropriate). 

> Has this been observed by someone else ? Am I doing something wrong here ?
> the way I set up my matrices is
> U - m*m
> VT - n*n
> S - min(m,n) . I have tried n also as in my case m<n.
> A matrix is m*n

I think everything is working fine here - you would have gotten clear error 
messages if any of these matrices were of wrong dimension. 

Christian

Re: [Lapackpp-devel] error handling

[Christian S. <sti...@tu...>]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

Dominik Wagenfuehr schrieb:
>> As for the output - I think this is usually redirected to
> ~/.xsession-errors,
> 
> What's the Windows pendant?

I don't know. I don't know windows that much.

>> and at the next failed assert() you will get a thrown LaException
>> instead of a failed assert. Is this what you wanted?
> 
> I will take a look at it. Maybe I should describe how I do my error
> handling in libraries. :)
> 
> I have a special static class called MyMessage with different enums for
> each error that can occur. If in a function an error occures I do
> 
> MyMessage::Set(ERROR_DEFINEMENT, "Class::FunctionA");
> return(true);
> 
> The functionB that called this failed functionA sees in the return
> argument that something went wrong. It now can decide if to goes one or do
> 
> MyMessage::Append("Class::FunctionB");
> return(true);
> 
> In this way I will create an error stack in the static class. At top
> level (that means in the GUI in most cases) I finally now can get the
> error and the stack
> 
> MyMessage::GetMessage() and MyMessage::GetStack()
> 
> and give the user a correct output. This means in a console application
> a console output, in a GUI-app maybe a little window message.

I think you can do this very well by using the new modifications in
laexcp.h, and then you enclose each lapackpp call by

try{
 // ...
}
catch (LaException &e) {
  MyMessage::Set(BLA, e.what() );
}

and continue accordingly.

> Further sometimes the user can decide what he wants to do now. In most
> cases he should copy the message, save all his data (as far as this is
> possible) and send a bug report to me. But if the GUI  just crashes
> (maybe after an hour of computing something) without any output the user
> may become a little impatient. ;)
> 
> I know that this is complicated especially implementing the return
> values of all functions (it took me a dew days to change my code.) And
> now I could not return something else but do most things by call by
> reference. The advantage of it is that my library (or better the program
> of a user with my library) never crashes unless the programmer does not
> catch my error output and just go on as nothing happens.
> 
> You can read something about this indirect error control in
> Piegl/Tiller: "The NURBS Book", p. 607, chapter "B-spline Programming
> Concept / Error control". There my be better programming books for it
> but this one explained it clear enough. :)

This kind of error handling seems a bit twisted to me, and at least the
"one static class" concept isn't any different than the standard libc
approach with one static variable errno. Whatever. I guess it will work
anyway.

Have fun,

Christian
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.2.1 (MingW32)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iQCVAwUBRWcSK2XAi+BfhivFAQIUsAP+Phm/2i46JL2y66nkCvCs0U/pmEHHF7KA
2733vCd20WC8KReXdret1dbE8eWCO3tf4D3hpJ3hEJ9edFGPtZEgDrPuRynt8cXZ
oYSDBRok/fHCZFyje7TU5goWRywhOqWwdellG/zQW4N+UTv4NLBlxAgUh/xeMXbz
i10y18RRIuw=
=zNOb
-----END PGP SIGNATURE-----

Re: [Lapackpp-devel] error handling

[Dominik W. <dom...@ar...>]

> As for the output - I think this is usually redirected to
~/.xsession-errors,

What's the Windows pendant?

> and you could just as well redirect the stdout of your GUI-Tool to
somewhere else.

Unfortunately I'm not programming the GUI at all. I just programm
another library used by a GUI. Somewhat complicated here...

> and at the next failed assert() you will get a thrown LaException
instead of
> a failed assert. Is this what you wanted?

I will take a look at it. Maybe I should describe how I do my error
handling in libraries. :)

I have a special static class called MyMessage with different enums for
each error that can occur. If in a function an error occures I do

MyMessage::Set(ERROR_DEFINEMENT, "Class::FunctionA");
return(true);

The functionB that called this failed functionA sees in the return
argument that something went wrong. It now can decide if to goes one or do

MyMessage::Append("Class::FunctionB");
return(true);

In this way I will create an error stack in the static class. At top
level (that means in the GUI in most cases) I finally now can get the
error and the stack

MyMessage::GetMessage() and MyMessage::GetStack()

and give the user a correct output. This means in a console application
a console output, in a GUI-app maybe a little window message.

Further sometimes the user can decide what he wants to do now. In most
cases he should copy the message, save all his data (as far as this is
possible) and send a bug report to me. But if the GUI  just crashes
(maybe after an hour of computing something) without any output the user
may become a little impatient. ;)

I know that this is complicated especially implementing the return
values of all functions (it took me a dew days to change my code.) And
now I could not return something else but do most things by call by
reference. The advantage of it is that my library (or better the program
of a user with my library) never crashes unless the programmer does not
catch my error output and just go on as nothing happens.

You can read something about this indirect error control in
Piegl/Tiller: "The NURBS Book", p. 607, chapter "B-spline Programming
Concept / Error control". There my be better programming books for it
but this one explained it clear enough. :)

Greetings, Dominik

Re: [Lapackpp-devel] error handling

[Christian S. <sti...@tu...>]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

Dominik Wagenfuehr schrieb:
>> However, the assertions in lapackpp are intended to check for
> *programming* errors.
> 
>> a failed assertion should indeed print an error message including the
>> source code file name, line number, and the failed condition.
> 
> Without a console no output! :) As I said before if you use your in lib
> in a GUI-Tool the tool just quits without any hint why it does crash...
> Of course it's a programming error of the GUI-programmer because he has
> made a mistake in matrix allocation or something, but I do not know if a
> user cares about that. He only sees his GUI crash and loose all data. ;)

I still believe modifying the assertions is the wrong place to do
something about this. As for the output - I think this is usually
redirected to ~/.xsession-errors, and you could just as well redirect
the stdout of your GUI-Tool to somewhere else. Additionally, you can
install a different handler for abort()... but whatever. You've
challenged me to do something here.

By now I've modified the LaException class so that you can switch the
assert() behaviour from the normal abort() to throwing a LaException
exception. To do this, you will have to uncomment some lines in the
laexcp.h header file, which will effectively redefine the assert() macro
so that it throws an exception instead of aborting. So: Checkout the
latest CVS code; change the laexcp.h header; compiler and install; and
at the next failed assert() you will get a thrown LaException instead of
a failed assert. Is this what you wanted?

Christian
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.2.1 (MingW32)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iQCVAwUBRWcIe2XAi+BfhivFAQKpLQP/dtQ7x+58g+nqSPPSIR5azNkjOB0u2wQ4
3Iqr331itGsngCGs9z9ZFcBH8YM2R5xQ/AS+s5mkE4/PGFsOka/pwOVlDF2dh+gd
Eq1peR918lFLiiRnDm5fIsXgGdEocjJO+cwssZ1S6tzVIu5VYv8HWypW9YtlUkCg
rSwydZirhRw=
=/Bp9
-----END PGP SIGNATURE-----

Re: [Lapackpp-devel] Symmetric (packed) matrices

[Dominik W. <dom...@ar...>]

>> Do you know why LaSymmMatDouble needs two arguments in the constructor?
>> Unfortunately it's not even checked if they are the same, I mean "assert
>> ( m ==n );" ...

> Interesting :-) . No, I don't know. In the current form this is
> obviously wrong. The LaSymmMatDouble matrix, as well as the
> LowerTriangBlaBlaBla matrix, should have a one-argument constructor each.

I agree. At this time someone can allocate a "symmetric" matrix with
dimension 10x100 if he likes to do. I do not know what the result is if
you call a Lapack-function with it...

I don't use symmetric matrices that much so I do not know why it's
"derived" from a triangular matrix and why it uses two dimension. Maybe
I have time to review the code at my holidays at end of december.

> However, currently lapackpp doesn't use the "packed storage" matrix
> format at all. [...] You can now either add the packed storage with its
> related solver functions as an additional set of classes and functions
> (LaPSymmMatDouble, for example).

As I said above I do not use symmetric matrices and I do not know much
about the storage scheme for Lapack. So I'm afraid I won't be a help here...

> Or we can decide to switch the existing
> symmetric classes over to packed storage. Unfortunately I have no idea
> whether anyone relies on the existing symmetric classes.

Hm, good question. As in lapackpp-classes you have no choice to directly
access the "dead" part of the memory. So it may be no difference if we
change it to packed storage. But as you point out, not all
Lapack-functions are implemented for packed storange so a new class
would be much better.

Greetings, Dee

Re: [Lapackpp-devel] error handling

[Dominik W. <dom...@ar...>]

> However, the assertions in lapackpp are intended to check for
*programming* errors.

Okay. I have another aim here... In my libs I try to catch all errors
even programming errors. But it's okay if you do not want to add a real
error class.

> a failed assertion should indeed print an error message including the
> source code file name, line number, and the failed condition.

Without a console no output! :) As I said before if you use your in lib
in a GUI-Tool the tool just quits without any hint why it does crash...
Of course it's a programming error of the GUI-programmer because he has
made a mistake in matrix allocation or something, but I do not know if a
user cares about that. He only sees his GUI crash and loose all data. ;)

Greetings, Dominik

Re: [Lapackpp-devel] New resize-methods

[Christian S. <sti...@tu...>]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

Dominik Wagenfuehr schrieb:
> could you please add the resize-methods in the files attached (These are
> the changed CVS-files!). I need to resize a symmetric matrix and a
> vector and therefore had to implement it in these classes. I have
> changed nothing else, so the matrix classes are far from optimized.

Ok, applied. Thanks a lot. I've modified your proposed version in the
LaVectorDouble in order to have the same behaviour as the existing size
arguments in the constructor.

> Do you know why LaSymmMatDouble needs two arguments in the constructor?
> Unfortunately it's not even checked if they are the same, I mean "assert
> ( m ==n );" ...

Interesting :-) . No, I don't know. In the current form this is
obviously wrong. The LaSymmMatDouble matrix, as well as the
LowerTriangBlaBlaBla matrix, should have a one-argument constructor each.

Regards,

Christian
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.2.1 (MingW32)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iQCVAwUBRWQ882XAi+BfhivFAQK6HgQAqCJePmHhzxP19t4DGn3+vC/KuaJJ7FHM
yWvl2PLLOeyOG7L2/eOidRBrhz0JYer7qiF+YZYb8aYZPTvXNCtCThqgUyNTCSAG
f5Kv+rbOVM4D7c3tp5PbmQAJ/tMFHsTma8QQJexW6kgfnp6nXf1AofZgOK3tK0MB
y4Jo5IoYUrI=
=EDWd
-----END PGP SIGNATURE-----

[Lapackpp-devel] New resize-methods

[Dominik W. <dom...@ar...>]

Hello,

could you please add the resize-methods in the files attached (These are
the changed CVS-files!). I need to resize a symmetric matrix and a
vector and therefore had to implement it in these classes. I have
changed nothing else, so the matrix classes are far from optimized.

Do you know why LaSymmMatDouble needs two arguments in the constructor?
Unfortunately it's not even checked if they are the same, I mean "assert
( m ==n );" ...

Greetings, Dominik

Re: [Lapackpp-devel] error handling

[Christian S. <sti...@tu...>]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

Dominik Wagenfuehr schrieb:
> BTW: Have you thought about a better error handling? I think a library
> should never ever break a program. But your assertion does nothing else.
> The program quits without any message. I use Lapack++ in a GUI so I only
> see an output if I create console test case and try to reproduce the
> crash. Anyway: A backend should never crash a frontend.

Yes and no. For "normal" error events that can occur during runtime
("file not found" or similar) this is perfectly true. However, the
assertions in lapackpp are intended to check for *programming* errors.

The header file documentation clearly says (well, it *should* say) that
for certain operations, certain conditions have to be fulfilled. E.g.,
the matrix dimensions of Blas_Mat_Mat_Mult must fit to each other. This
is a condition which the application can easily check *before* it calls
these functions - and for this reason lapackpp can rely upon these
conditions being fulfilled. That's the meaning of such an assertion.
If those conditions are not met, then most probably the application
programmer just made a programming error. Having some exceptions instead
of assertions IMHO don't change this much - in any case a user cannot do
anything about this failed condition; the programmer must fix the
application anyway. Therefore I want to have these checks as easy as
possible, and simply writing assert(...) IMHO is the easiest way.

In principle each of these assertions should be mentioned in the
documentation, like "This function requires the following assertions to
be true: A.size(0) == bla bla bla". In practice, the assertion are a
very helpful tool to find out programming errors as early as possible.
As a rule of thumb I'd suggest in every (non-trivial) function, the
correct value range of the function arguments should be checked by
assertions.

Depending on your libc and the compile flags during compiling lapackpp,
a failed assertion should indeed print an error message including the
source code file name, line number, and the failed condition. If this
isn't the case, you should consider recompiling lapackpp with
- --enable-debug and/or CXXFLAGS="-g -O0" and so on.

Regards,

Christian
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.2.1 (MingW32)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iQCVAwUBRWMCvGXAi+BfhivFAQLgsAP9ElQxrNBXVeYjNMCgGEphCy0kyDzRua85
Vt8Wf6spWEqzN+dR5/5oXlLBl8wQIX80Ydga2N+xkSp40jODxz52yq5UZJvsP6wu
ihq6kOmSlGOcwUhUNyofvEbZfzrypL0dCIEjiyRzPdF4PZ2j0jeNVz1o50lioUhP
CZ2/rddqX3E=
=tYvX
-----END PGP SIGNATURE-----

Re: [Lapackpp-devel] committed files

[Christian S. <sti...@tu...>]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

Hi Dominik,

thanks a lot for your contributions. By now I've reviewed all of them
and commited them into CVS. In blas3pp.{h,cc} I've replaced the "char"
argument back to "bool" just as I proposed before.

Dominik Wagenfuehr schrieb:
> Attached you find a corrected version fo blas3pp.cc (There were some
> errors in it, sorry for that!) and a test case tsybmd.cc for
> LaSymmBandMatDouble and LaSymmBandFactDouble. The classes were correct
> only my example was wrong. ;)

In this mail from 11:06 Uhr you attached a file blas3pp.cc that was
identical to the CVS version of that file at that time, not one with
some errors fixed. I've committed your file of yesterday into CVS. If
there are still errors, feel free to send me an updated / corrected
version now.

> Another question: doxygen does not show the new comments in blas3pp.h
> because of the #ifdef before. What do we do here?

There doesn't seem to be an easy solution in terms of doxygen commands.
However, these whole functions don't need to have the full class
definition of LaSymmMatDouble available - it is sufficient for the
compiler to have a forward declaration of those classes only. Therefore
I've removed the ifdef that asks for the included header, and instead
I've added the appropriate forward declaration. As a nice side effect,
this solved the doxygen issues as well.

Christian

-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.2.1 (MingW32)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iQCVAwUBRWLu92XAi+BfhivFAQK5oQP/TkizqagPUKTooN1n44wgAREaMCK1s3i8
8omwkJVgujk54GiQsw+DTUun1YV1SYlfh9I+H+aolAE9bcWAQH7I7DBi8SoDO1h+
hmcakLyJb4hN+cV7yPrGzkMN/6nAUbwKisDlTkYnXL6KMDggX0FQJTEVi/rQOg3p
o86MA1+Utsc=
=1FGS
-----END PGP SIGNATURE-----

Re: [Lapackpp-devel] symmetry in packed storage, other questions

[Christian S. <sti...@tu...>]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

Dominik Wagenfuehr schrieb:
>> And I also cannot think of an easy solution here that would be all of
>> this: 1. easy
>> for non-mathematics, 2. give the correct results, and 3. is powerful
>> enough to actually achieve this memory optimization. Sigh.
> 
> Did you take a look at the PP-classes that means "real symmetric
> positive definite matrix stored in packed format", see
> http://docs.sun.com/source/806-7993/dpptrs.html. I haven't found a
> description what this matrices looks like but I would guess that only
> one half of the matrix is stored.

Thanks for pointing this out. The link to "dpptrs" and "dpptrf" can be
added to the documentation of LaSymmMatDouble, just for future
reference. The description of the format is in the description of the
input argument "A" and/or in the man-page of dpptrf, "man dpptrf". See
also http://www.netlib.org/lapack/lug/node121.html and its subsections
which is probably an even better explanation.

However, currently lapackpp doesn't use the "packed storage" matrix
format at all. That means none of the packed-storage LAPACK function are
declared in the C header files include/blas*.h. Obviously all the
"symmetric" classes and related LAPACK solver functions were only using
the normal storage scheme, where no storage savings were achived.

You can now either add the packed storage with its related solver
functions as an additional set of classes and functions
(LaPSymmMatDouble, for example). Or we can decide to switch the existing
symmetric classes over to packed storage. Unfortunately I have no idea
whether anyone relies on the existing symmetric classes. To be on the
safe side, I'd definitely prefer to have the packed storage as
additional classes...

Funny enough, LAPACK seems to have some more functionality for the
conventional storage than for the packed storage. For example, the
eigenvalue decompositions (dsyev,dsyevd, dsyevr, dsyevx) seems to be
unavailable for packed storage.

>> So if you're working on "dsyrk" of Blas_R1_Update, 
>> you can just as well add an alternative
>> function where the transpose setting is controlled by an additional
>> boolean argument. 
> 
> I will do so. But boolean is not enough because it could be N, T or C as
> argument.  Although I have never heard of Blas_R1_Update before I will
> use the same name and will set the last char parameter as default='N'.
> So we have backward compatibility and could use the new function. :)

What is the difference between 'T' or 'C'? I don't see any. Therefore
there are only two possibilities.

>> You're still invited to give me your sourceforge username and join the
>> lapackpp project as a developer :-)
> 
> Yes. The username is "deedw" without quotes of course. :) But I have
> never used CVS so it may be that I will send the new files to the
> mailing list.

You are added as a developer now. If the anonymous CVS access is too
slow, you can switch over to the developer access. Submitting changes by
email is still fine, though.

Christian
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.2.1 (MingW32)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iQCVAwUBRWLVHmXAi+BfhivFAQJiGgP8CxMVDuMpFsBFOYaE878YOvBYWGNb0DFQ
ATECzWqQOn/YwtS8OA4T9bRrraLKPlS0EhrNtyBGpS5MMaSuVOvmidoWedLyU6/y
yaSwRb3+5aplAorq7HxeUsGRymud4gj5Hiiw/HnGsZR9Co3IvFlDk8h6VIIaLgg/
ObWesLXvcq8=
=DADX
-----END PGP SIGNATURE-----

Re: [Lapackpp-devel] inline, symmetry, transposition, other questions

[Dominik W. <dom...@ar...>]

Attached you find a corrected version fo blas3pp.cc (There were some
errors in it, sorry for that!) and a test case tsybmd.cc for
LaSymmBandMatDouble and LaSymmBandFactDouble. The classes were correct
only my example was wrong. ;)

Another question: doxygen does not show the new comments in blas3pp.h
because of the #ifdef before. What do we do here?

Greetings, Dominik

Re: [Lapackpp-devel] inline, symmetry, transposition, other questions

[Dominik W. <dom...@ar...>]

Attached you find the new files. This time I used the cvs-ones. ;)

What I have changed:

blas3pp.h - add documentation and default value for trans
blas3pp.c - set default value and add assert

sybfd.h - delete old factorize class and use functions only
sybfd.cc - create new file with functions

sybmd.h - set correct class description and set size(0)=size(1)

> Look into matrix/test/tgc.cc which should give a good idea about what
to test
> for such a matrix class.

Unfortunately my test class does not what it should do.
LaSymmBandMatDouble is okay but LaSymmBandFactDouble returns a "wrong"
Cholesky factorization. I'll try to fix the error.

BTW: Have you thought about a better error handling? I think a library
should never ever break a program. But your assertion does nothing else.
The program quits without any message. I use Lapack++ in a GUI so I only
see an output if I create console test case and try to reproduce the
crash. Anyway: A backend should never crash a frontend.

Greetings, Dee

Re: [Lapackpp-devel] inline, symmetry, transposition, other questions

[Dominik W. <dom...@ar...>]

> Not having a class
> at all but instead offer several functions (which will then be wrapper
> functions for the LAPACK factorizations) is also one such possibility.

I think this makes more sense for sybfd.h than a a whole class you don't
use. So I will delete the class and will only use functions. The first
function will then be delete as well.

> And I also cannot think of an easy solution here that would be all of
this: 1. easy
> for non-mathematics, 2. give the correct results, and 3. is powerful
enough
> to actually achieve this memory optimization. Sigh.

Did you take a look at the PP-classes that means "real symmetric
positive definite matrix stored in packed format", see
http://docs.sun.com/source/806-7993/dpptrs.html. I haven't found a
description what this matrices looks like but I would guess that only
one half of the matrix is stored.

> So if you're
> working on "dsyrk" of Blas_R1_Update, you can just as well add an
alternative
> function where the transpose setting is controlled by an additional
boolean
> argument. And since you're adding a new function, you can even choose the
> name for yourself :-)

I will do so. But boolean is not enough because it could be N, T or C as
argument.  Although I have never heard of Blas_R1_Update before I will
use the same name and will set the last char parameter as default='N'.
So we have backward compatibility and could use the new function. :)

> You're still invited to give me your sourceforge username and join the
> lapackpp project as a developer :-)

Yes. The username is "deedw" without quotes of course. :) But I have
never used CVS so it may be that I will send the new files to the
mailing list.

Greetings, Dominik

Re: [Lapackpp-devel] inline, symmetry, transposition, other questions

[Christian S. <sti...@tu...>]

Am Sonntag, 19. November 2006 09:33 schrieb Dominik Wagenfuehr:
> > Please always use the code from CVS, now that you're actually working
> > with it. I already did some changes.
>
> You see, this is the first "project" I working with where several other
> people are workig with. Sorry again.

No problem. 

> > As for inline vs. non-inline: I want to reduce the "inline" functions
> > to a minimum, because they contradict the OO-paradigm of "hiding the
> > implementation" (and make debugging harder, among other things).
>
> Hm, maybe I'm more mathematican than programmer. I've never heard of
> these "inline" before and have left these entries as they were before. I
> have read what "inline" stands for but I cannot decide where it useful,
> because I would never use it on my own. But I try to guess the correct
> usage.

Ok, and it's of course a good approach to leave the stuff as it is as long as 
you  don't know the reasons for it. 

So for "inline" in general: Usually you don't need this and don't have to 
think about it - simply declare the functions in the .h file and define 
(implement) them in the .cc file. However, at runtime there is a certain 
fixed overhead for each time when a function is called. For some specific 
functions that are called a *lot* one wants to avoid that overhead, and the 
compiler can do that by directly inserting "the body" of the called function 
into the place where it has been called. This is called inlining. For this to 
work, the function's definition (not only the declaration) has to be in the 
header file, which contradicts the whole point of OOP "hiding the 
implementation" paradigms. Therefore it should be used only when one is 
really sure about the benefit. The obvious example is operator(): This is 
called a *lot* as soon as a client application uses a(i,j)=... inside of the 
O(n^2) loop. And its function definition contains very few statements - 
probably only one if-clause and one address calculation. In that case 
inlining that function may give a noticable performance improvement at 
runtime. In those classes that I've reviewed already, I made sure to keep 
only those methods inline. Everything else is just a normal function.

> > Also, please think about adding a test program (in
> > matrix/test/blablabla.cc) so that the correct functioning can be tested
> > automatically.
>
> I will take a look at the other files in there. I doesn't have a clue
> what I should test. ;)

Look into matrix/test/tgc.cc which should give a good idea about what to test 
for such a matrix class.

> > Actually I'm not so sure whether the function LaSymmBandMatFactorize()
> > is useful. If that functionality better lives inside some
> > LaSymmBandFactDouble class methods then I'd encourage you to remove the
> > function and move the code to suitable class methods.
>
> Me too. In my first changes I do not want to change too many things at
> once. But I haven't understood it either why not assigning / referencing
> a Matrix in the constructor and do the factorization with a class
> method. I will do this the next time.

I think there are several possibilities how such a factorization class can 
make sense. IMHO include/gfqrc.h is one such possibility. Not having a class 
at all but instead offer several functions (which will then be wrapper 
functions for the LAPACK factorizations) is also one such possibility. Your 
current sybfd.h is probably also fine. The current (very old) version of 
gfd.h, on the other hand, is no useful possibility.

> >> 2. What is the correct meaning of A.size(int d) ?
> >
> > No, I think size() should give the dimensions in terms of the
> > "mathematical" viewpoint, not in terms of the internal storage
>
> Okay, I will change this.

Fine.

> > Is sybmd.h really a *symmetric* banded matrix? Is this fact of symmetry
> > used anywhere inside this class?
> > Surely, in terms of optimized storage,
> > we would need to store only 1*p+1 instead of 2*p+1 elements...
>
> No and yes. ;) The memory space is not optimized for the same reason as
> LaSymmMatDouble uses the same memory as LaGenMatDouble. In each
> Lapack-function you can use the upper or lower part of a matrix to do
> the changes. Of course only one side is necessary in symmetric matrices
> but that is not interesting for lapack-functions. It must use the whole
> memory. Do you know how the memory is accessed by lapack or blas? Maybe
> we could change it.
>
> In other methods like copy or operator() the symmetry is used of course.

Ah! Oh! Ok! Suddenly I've understood the whole point of symmetric matrices in 
LAPACK! No, I think the rationale in LAPACK is that you will store two 
symmetric matrices U and L in the memory space of one non-symmetric matrix A 
- for U you only access the upper part of A, and for L you only access the 
lower part. That's why all these functions (e.g. dsysv) say that you can 
choose which part of A is being used as the symmetric matrix. 

However, that means that memory saving is only obtained when U references the 
same memory space as L. That can easily be done in Fortran or C, when the 
programmer allocates the matrix A by himself/herself anyway. However, in 
lapackpp when the LaSymmMatDouble class allocates the matrix A itself, it 
obviously cannot use the same allocated storage as another LaSymmMatDouble 
currently. So right now I'd say the classes are prepared to use the memory 
space of A for both one U and one L, but there are not class methods 
available to actually put both U and L into the same memory space A. And I 
also cannot think of an easy solution here that would be all of this: 1. easy 
for non-mathematics, 2. give the correct results, and 3. is powerful enough 
to actually achieve this memory optimization. Sigh.

> > Seems to me the matrix class itself rather implements a general band
> > matrix
>
> No. It is a banded matrix with same bandwith in subdiags and superdiags.
> :) So it's a type of special band matrix in terms of memory usage. But
> the class LaSymmBandMatDouble uses internally only half of the memory. I
> do not know much about Lapack storage but I think you could not and
> should not change it or Lapack-functions will fail. But I do not know it
> for sure.

Yes. Agreed.

> 1. In laslv.h you forget to
>
> #include LA_VECTOR_DOUBLE_H
>
> You will get an error including this file without including lavd.h
> yourself before.

Thanks. Applied to CVS.

> 2. Last time you suggested Blas_R1_Update for C := A' * A. Unfortunately
> this function is broken or does not work that way. The reason is that
> the lapack-function inside will alway perform
>
> C := alpha*A*A' + beta*C
>
> because "trans" is set to 'N'. You have no chance to do C := A' * A.
>
> The same with Blas_R2_Update. You only perform
>
> C := alpha*A*B' + alpha*B*A' + beta*C
>
> but alpha*A'*B + alpha*B'*A + beta*C is not acessable. You should add
> some new parameter for these functions or create different function as
> you have done with Blas_Mat_Mat_Mult and Blas_Mat_Trans_Mat_Mult.

Yes, I know the "transpose" argument isn't solved in a nice way in lapackpp. 
The current functions use the "Different function names for different 
arguments" approach. Alternatively, one could boolean argument(s) to control 
the transpositions. Actually I'd probably prefer the latter approach with 
explicit boolean arguments (with reasonable default values). So if you're 
working on "dsyrk" of Blas_R1_Update, you can just as well add an alternative 
function where the transpose setting is controlled by an additional boolean 
argument. And since you're adding a new function, you can even choose the 
name for yourself :-)  For the existing matrix multiplications I might add 
something like that in the future. (The existing functions will still be 
kept, of course.)

You're still invited to give me your sourceforge username and join the 
lapackpp project as a developer :-)

Regards,

Christian

Re: [Lapackpp-devel] Error in LaGenMatDouble::copy(...) and other questions

[Dominik W. <dom...@ar...>]

> You didn't send it to the mailing list, only to me directly. In the
> future please make sure you always send a CC to the mailing list

I'm sorry for that. Would be easier, if the mailing list would be
standard and your mail as CC. ;)

> Please always use the code from CVS, now that you're actually working
> with it. I already did some changes.

You see, this is the first "project" I working with where several other
people are workig with. Sorry again.

> As for inline vs. non-inline: I want to reduce the "inline" functions
> to a minimum, because they contradict the OO-paradigm of "hiding the
> implementation" (and make debugging harder, among other things).

Hm, maybe I'm more mathematican than programmer. I've never heard of
these "inline" before and have left these entries as they were before. I
have read what "inline" stands for but I cannot decide where it useful,
because I would never use it on my own. But I try to guess the correct
usage.

> If possible, please think about splitting the code to a .cc file so that
> only the useful functions stay inline and the others are defined in the
> .cc file.

I can do this the next week.

> Also, please think about adding a test program (in
> matrix/test/blablabla.cc) so that the correct functioning can be tested
> automatically.

I will take a look at the other files in there. I doesn't have a clue
what I should test. ;)

> Actually I'm not so sure whether the function LaSymmBandMatFactorize()
> is useful. If that functionality better lives inside some
> LaSymmBandFactDouble class methods then I'd encourage you to remove the
> function and move the code to suitable class methods.

Me too. In my first changes I do not want to change too many things at
once. But I haven't understood it either why not assigning / referencing
a Matrix in the constructor and do the factorization with a class
method. I will do this the next time.

>> 2. What is the correct meaning of A.size(int d) ?
> No, I think size() should give the dimensions in terms of the
> "mathematical" viewpoint, not in terms of the internal storage

Okay, I will change this.

I hope it's okay, if I paste the other mail below for answering...

> Are you subscribed to the lapackpp-cvs list? You should.

I wilI! :)

> Is sybmd.h really a *symmetric* banded matrix? Is this fact of symmetry
> used anywhere inside this class?
> Surely, in terms of optimized storage,
> we would need to store only 1*p+1 instead of 2*p+1 elements...

No and yes. ;) The memory space is not optimized for the same reason as
LaSymmMatDouble uses the same memory as LaGenMatDouble. In each
Lapack-function you can use the upper or lower part of a matrix to do
the changes. Of course only one side is necessary in symmetric matrices
but that is not interesting for lapack-functions. It must use the whole
memory. Do you know how the memory is accessed by lapack or blas? Maybe
we could change it.

In other methods like copy or operator() the symmetry is used of course.

> Seems to me the matrix class itself rather implements a general band
matrix

No. It is a banded matrix with same bandwith in subdiags and superdiags.
:) So it's a type of special band matrix in terms of memory usage. But
the class LaSymmBandMatDouble uses internally only half of the memory. I
do not know much about Lapack storage but I think you could not and
should not change it or Lapack-functions will fail. But I do not know it
for sure.

Again two other hints:

1. In laslv.h you forget to

#include LA_VECTOR_DOUBLE_H

You will get an error including this file without including lavd.h
yourself before.

2. Last time you suggested Blas_R1_Update for C := A' * A. Unfortunately
this function is broken or does not work that way. The reason is that
the lapack-function inside will alway perform

C := alpha*A*A' + beta*C

because "trans" is set to 'N'. You have no chance to do C := A' * A.

The same with Blas_R2_Update. You only perform

C := alpha*A*B' + alpha*B*A' + beta*C

but alpha*A'*B + alpha*B'*A + beta*C is not acessable. You should add
some new parameter for these functions or create different function as
you have done with Blas_Mat_Mat_Mult and Blas_Mat_Trans_Mat_Mult.

Greetings, Dominik

Re: [Lapackpp-devel] Error in LaGenMatDouble::copy(...) and other questions

[Christian S. <sti...@tu...>]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

Dominik Wagenfuehr schrieb:
> Please take a close look at the documentation. I'm not that much
> experienced with doxygen syntax... Most parts I have copied from gmd.h.
> I hope this is okay. :)

Review finished. Feel free to work with the latest CVS code now. Are you
subscribed to the lapackpp-cvs list? You should.

But after I read the wikipedia page
http://en.wikipedia.org/wiki/Band_matrix I actually asked myself: Is
sybmd.h really a *symmetric* banded matrix? Is this fact of symmetry
used anywhere inside this class? Surely, in terms of optimized storage,
we would need to store only 1*p+1 instead of 2*p+1 elements... Seems to
me the matrix class itself rather implements a general band matrix, and
only the particular LAPACK functions will use this only as symmetric
matrix. Maybe in that case we might consider renaming this class to
LaBandMatDouble and have a different class make up a
LaSymmBandMatDouble? Whatever. I don't really care.

Christian


-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.2.1 (MingW32)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iQCVAwUBRV2Yv2XAi+BfhivFAQKItAP8DOd0a01RlOPPDpaaqQtkMrjBOc00kR7c
6BApgPQtuNeXn4vKSbxNh2ZTZE0iw9wnq9xIln9x7+gSvgwqlz6MmiTLxBctADno
DS/3roxSQOe63F8srQtq21oxrMgikXeAtRG3qRLzb0/R9+zE2mGwvNGCWblG2FTH
9iNIyr4WJKI=
=9FlK
-----END PGP SIGNATURE-----