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

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

1. I havent' found the exact error yet, but if I debug my program I will
get a segmentation fault in

LaGenMatDouble B(A);

or

LaGenMatDouble B;
B.copy(A);

A is a (4 x 2)-matrix with some values I have calculated before.
Unfortunately I cannot debug the lapack-libary so I don't know where the
error could be. I try to extract the relevant code parts from my program
to reconstruct it.

2. Is it possible that you forgot to create the copy-function for
LaSymmBandMatDouble? I only get

/usr/local/include/lapackpp/sybmd.h:51: warning: inline function
'LaSymmBandMatDouble& LaSymmBandMatDouble::copy(const
LaSymmBandMatDouble&)' used but never defined

I have created my own copy-function yet but as I see you have used some
templates for this. I don't understand much about templates but it would
be great if you can fix it.

3. Is there such a simple thing like A := A + alpha*B ? I have only
found Blas_Mat_Mat_Mult and actually I perform

Blas_Mat_Mat_Mult(B, I, A, alpha);

where I as the Identity Matrix. Especially for LaSymmBandMatDouble you
could save a lot lot lot of time if there would be such a simple
operation for different types of matrices. I don't know if there exists
a BLAS-function for it (I haven't found it yet), but I have named my
function

void Blas_Mat_Mat_Add(LaSymmBandMatDouble& A, const LaSymmBandMatDouble&
B, const double& alpha = 1.0)
{
  assert(A.size(0) == B.size(0));
  assert(A.size(1) == B.size(1));

  int n=A.size(0), p=A.subdiags(), n2;
  int i,j;
 
  for (i=0; i<n; i++)
  {
    n2 = (i+p+1 < n) ? i+p+1 : n;
    for (j=0, j<n2; j++)
    {
      A(i,j)+=alpha*B(i,j);
    }
  }
}

4. Same with the operation A = B^T * B. There seems to be no
BLAS-function for it even if you know that A is symmetric. You could
save half of the time. And further if you know that B has a special
structure and A will be banded symmetric...

Greetings, Dominik

PS: I have tested the functions in sybmd.h and sybfd.h and they seem to
work fine. (Ignoring the little copy-incident above!) I will try to
comment the functions...