diff -N -r -c HLib-1.3/Library/h2virtual.c HLib-1.3new/Library/h2virtual.c
*** HLib-1.3/Library/h2virtual.c	2004-11-23 18:16:20.000000000 +0100
--- HLib-1.3new/Library/h2virtual.c	2005-06-30 01:00:50.809099600 +0200
***************
*** 2567,2573 ****
    double *X;
    int ldX;
    double *Q;
-   double *Q2;
    double *lambda, *tau, *work;
    double h, fmax, sqeps;
    int k, k2, smax, kt, kmax;
--- 2567,2572 ----
***************
*** 2583,2589 ****
  
    sqeps = h2r->sqeps;
    kmax = h2r->kmax;
-   Q = h2r->Q;
  
    /* Get rank of original cluster basis */
  
--- 2582,2587 ----
***************
*** 2735,2740 ****
--- 2733,2741 ----
  
    /* Expand row or column matrix */
  
+   provide_array(m * k, &h2r->Q, &h2r->Qsize);
+   Q = h2r->Q;
+ 
    dgemm_("Transposed", "Not Transposed",
  	 &m, &k, &k,
  	 deins_,
***************
*** 2745,2803 ****
  
    /* Compute SVD */
  
!   if(m <= k) {
!     ldwork = 8 * k;
!     rank = m;
!     provide_array(rank + ldwork, &h2r->tmp, &h2r->tmpsize);
!     lambda = h2r->tmp;
!     work = lambda + rank;
      
!     dgesvd_("Overwrite Q with left SVs", "No right SVs",
! 	    &m, &k,
! 	    Q, &m,
! 	    lambda,
! 	    NULL, eins_,
! 	    NULL, eins_,
! 	    work, &ldwork,
! 	    &info);
!     assert(info == 0);
!   }
!   else {
!     ldwork = 8 * m;
!     rank = k;
!     provide_array(m * k + 2 * rank + ldwork, &h2r->tmp, &h2r->tmpsize);
!     Q2 = h2r->tmp;
!     lambda = Q2 + m * k;
!     tau = lambda + rank;
!     work = tau + rank;
! 
!     dgeqrf_(&m, &k, Q, &m, tau, work, &ldwork, &info);
!     assert(info == 0);
! 
!     copyblock_lapack(Q, m, Q2, m, m, k);
! 
!     for(j=0; j<k; j++)
!       for(i=j+1; i<m; i++)
! 	Q[i+j*m] = 0.0;
! 
!     dgesvd_("Overwrite Q with left SVs", "No right SVs",
! 	    &k, &k,
! 	    Q, &m,
! 	    lambda,
! 	    NULL, eins_,
! 	    NULL, eins_,
! 	    work, &ldwork,
! 	    &info);
!     assert(info == 0);
! 
!     dormqr_("Left side", "Transposed",
! 	    &m, &k, &k,
! 	    Q2, &m, tau,
! 	    Q, &m,
! 	    work, &ldwork,
! 	    &info);
!     assert(info == 0);
!   }
  
    /* Determine rank of new cluster basis */
  
--- 2746,2766 ----
  
    /* Compute SVD */
  
!   ldwork = 8 * k;
!   rank = imin(k, m);
!   provide_array(rank + ldwork, &h2r->tmp, &h2r->tmpsize);
!   lambda = h2r->tmp;
!   work = lambda + rank;
      
!   dgesvd_("Overwrite Q with left SVs", "No right SVs",
! 	  &m, &k,
! 	  Q, &m,
! 	  lambda,
! 	  NULL, eins_,
! 	  NULL, eins_,
! 	  work, &ldwork,
! 	  &info);
!   assert(info == 0);
  
    /* Determine rank of new cluster basis */