Added extrow routine.

src/modules/psb_serial_mod.f90

@@ -288,6 +288,31 @@ module psb_serial_mod
      end subroutine psb_zspgtrow
   end interface
 
+  interface psb_sp_extrow
+     subroutine psb_dspextrow(irw,a,nz,ia,ja,val,info,iren,lrw)
+       use psb_spmat_type
+       type(psb_dspmat_type), intent(in) :: a
+       integer, intent(in)       :: irw
+       integer, intent(out)      :: nz
+       integer, intent(inout)    :: ia(:), ja(:)
+       real(kind(1.d0)),  intent(inout)    :: val(:)
+       integer, intent(in), target, optional :: iren(:)
+       integer, intent(in), optional :: lrw
+       integer, intent(out)  :: info
+     end subroutine psb_dspextrow
+     subroutine psb_zspextrow(irw,a,nz,ia,ja,val,info,iren,lrw)
+       use psb_spmat_type
+       type(psb_zspmat_type), intent(in) :: a
+       integer, intent(in)       :: irw
+       integer, intent(out)      :: nz
+       integer, intent(inout)    :: ia(:), ja(:)
+       complex(kind(1.d0)),  intent(inout)    :: val(:)
+       integer, intent(in), target, optional :: iren(:)
+       integer, intent(in), optional :: lrw
+       integer, intent(out)  :: info
+     end subroutine psb_zspextrow
+  end interface
+
   interface psb_neigh
      subroutine psb_dneigh(a,idx,neigh,n,info,lev)
        use psb_spmat_type

src/modules/psb_spmat_type.f90

@@ -125,6 +125,7 @@ contains
     type(psb_dspmat_type), intent(inout) :: mat
 
     nullify(mat%aspk,mat%ia1,mat%ia2,mat%pl,mat%pr)
+    mat%infoa(:) = 0
     mat%m=0
     mat%k=0
     mat%fida=''
@@ -614,6 +615,7 @@ contains
 
   end function psb_dspsizeof
 
+
   subroutine psb_dsp_free(a,info)
     implicit none
     !....Parameters...
@@ -649,6 +651,7 @@ contains
     type(psb_zspmat_type), intent(inout) :: mat
 
     nullify(mat%aspk,mat%ia1,mat%ia2,mat%pl,mat%pr)
+    mat%infoa(:) = 0
     mat%m=0
     mat%k=0
     mat%fida=''
@@ -1167,5 +1170,6 @@ contains
   End Subroutine psb_zsp_free
 
 
+
 end module psb_spmat_type
 

src/serial/Makefile

@@ -6,7 +6,7 @@ FOBJS = psb_cest.o psb_dcoins.o psb_dcsdp.o psb_dcsmm.o psb_dcsmv.o \
 	psb_dfixcoo.o psb_dipcoo2csr.o psb_dipcsr2coo.o psb_dneigh.o \
 	psb_dnumbmm.o psb_drwextd.o psb_dspgtdiag.o psb_dspgtrow.o \
 	psb_dspinfo.o psb_dspscal.o psb_dsymbmm.o psb_dtransp.o \
-	psb_dipcoo2csc.o lsame.o\
+	psb_dipcoo2csc.o psb_dspextrow.o lsame.o psb_zspextrow.o\
 	psb_zcsmm.o psb_zcsmv.o  psb_zspgtdiag.o psb_zspgtrow.o\
 	psb_zcsnmi.o psb_zcsrws.o psb_zcssm.o psb_zcssv.o psb_zcsdp.o\
 	psb_zfixcoo.o psb_zipcoo2csr.o psb_zipcsr2coo.o psb_zipcoo2csc.o \

src/serial/psb_dspextrow.f90

@@ -0,0 +1,162 @@
+!!$ 
+!!$              Parallel Sparse BLAS  v2.0
+!!$    (C) Copyright 2006 Salvatore Filippone    University of Rome Tor Vergata
+!!$                       Alfredo Buttari        University of Rome Tor Vergata
+!!$ 
+!!$  Redistribution and use in source and binary forms, with or without
+!!$  modification, are permitted provided that the following conditions
+!!$  are met:
+!!$    1. Redistributions of source code must retain the above copyright
+!!$       notice, this list of conditions and the following disclaimer.
+!!$    2. Redistributions in binary form must reproduce the above copyright
+!!$       notice, this list of conditions, and the following disclaimer in the
+!!$       documentation and/or other materials provided with the distribution.
+!!$    3. The name of the PSBLAS group or the names of its contributors may
+!!$       not be used to endorse or promote products derived from this
+!!$       software without specific written permission.
+!!$ 
+!!$  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
+!!$  ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+!!$  TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
+!!$  PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
+!!$  BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
+!!$  CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
+!!$  SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
+!!$  INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
+!!$  CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
+!!$  ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+!!$  POSSIBILITY OF SUCH DAMAGE.
+!!$ 
+!!$   
+! File:  psb_dspextrow.f90 
+! Subroutine: psb_dspextrow
+!    Gets one or more rows from a sparse matrix. 
+! Parameters:
+
+!*****************************************************************************
+!*                                                                           *
+!* Takes a specified row from matrix A and copies into NZ,IA,JA,VAL  in COO  *
+!* format.                                                                   *
+!*                                                                           *
+!*****************************************************************************
+subroutine psb_dspextrow(irw,a,nz,ia,ja,val,info,iren,lrw)
+  use psb_spmat_type
+  use psb_string_mod
+  type(psb_dspmat_type), intent(in) :: a
+  integer, intent(in)       :: irw
+  integer, intent(out)      :: nz
+  integer, intent(inout)    :: ia(:), ja(:)
+  real(kind(1.d0)),  intent(inout)    :: val(:)
+  integer, intent(in), target, optional :: iren(:)
+  integer, intent(in), optional :: lrw
+  integer, intent(out)  :: info
+  interface psb_spgtrow
+    subroutine psb_dspgtrow(irw,a,b,info,append,iren,lrw)
+      ! Output is always in  COO format  into B, irrespective of 
+      ! the input format 
+      use psb_spmat_type
+      use psb_const_mod
+      implicit none
+
+      type(psb_dspmat_type), intent(in)     :: a
+      integer, intent(in)                   :: irw
+      type(psb_dspmat_type), intent(inout)  :: b
+      integer,intent(out)                   :: info
+      logical, intent(in), optional         :: append
+      integer, intent(in), target, optional :: iren(:)
+      integer, intent(in), optional         :: lrw
+    end subroutine psb_dspgtrow
+  end interface
+
+  integer               :: lrw_, ierr(2), err_act
+  type(psb_dspmat_type) :: b
+  integer, pointer      :: iren_(:)
+  character(len=20)     :: name, ch_err
+
+
+  name='psb_sp_extrow'
+  info  = 0
+  call psb_erractionsave(err_act)
+  call psb_set_erraction(0)  
+
+  call psb_nullify_sp(b)
+
+  if (present(lrw)) then
+    lrw_ = lrw
+  else
+    lrw_ = irw
+  endif
+  if (lrw_ < irw) then
+    write(0,*) 'SPEXTROW input error: fixing lrw',irw,lrw_
+    lrw_ = irw
+  end if
+  call psb_sp_all(lrw_-irw+1,lrw_-irw+1,b,info)
+
+  if (present(iren)) then
+    call psb_spgtrow(irw,a,b,info,iren=iren,lrw=lrw_)
+  else 
+    call psb_spgtrow(irw,a,b,info,lrw=lrw_)
+  end if
+  if (info /= 0) then     
+    info=136
+    ch_err=a%fida(1:3)
+    call psb_errpush(info,name,a_err=ch_err)
+    goto 9999
+  end if
+
+  if (toupper(b%fida) /= 'COO') then 
+    info=4010
+    ch_err=a%fida(1:3)
+    call psb_errpush(info,name,a_err=ch_err)
+    goto 9999
+  endif
+
+  nz = b%infoa(psb_nnz_)
+
+  if (size(ia)>= nz) then 
+    ia(1:nz) = b%ia1(1:nz)
+  else
+    info    = 135
+    ierr(1) = 4
+    ierr(2) = size(ia)
+    call psb_errpush(info,name,i_err=ierr)
+    goto 9999
+  endif
+
+  if (size(ja)>= nz) then 
+    ja(1:nz) = b%ia2(1:nz)
+  else
+    info    = 135
+    ierr(1) = 5
+    ierr(2) = size(ja)
+    call psb_errpush(info,name,i_err=ierr)
+    goto 9999
+  endif
+
+  if (size(val)>= nz) then 
+    val(1:nz) = b%aspk(1:nz)
+  else
+    info    = 135
+    ierr(1) = 6
+    ierr(2) = size(val)
+    call psb_errpush(info,name,i_err=ierr)
+    goto 9999
+  endif
+
+
+  call psb_sp_free(b,info)
+
+  call psb_erractionrestore(err_act)
+  return
+
+9999 continue
+  call psb_erractionrestore(err_act)
+  if (err_act.eq.act_abort) then
+    call psb_error()
+    return
+  end if
+  return
+
+
+end subroutine psb_dspextrow
+

src/serial/psb_dspgtrow.f90

@@ -215,12 +215,11 @@ contains
         k=0
 
         do i=irw,lrw
-
           do j=a%ia2(i),a%ia2(i+1)-1
             k             = k + 1
-                b%aspk(nzb+k) = a%aspk(j)
             b%ia1(nzb+k)  = i
             b%ia2(nzb+k)  = a%ia1(j)
+            b%aspk(nzb+k) = a%aspk(j)
 !!$          write(0,*) 'csr_gtrow: in:',a%aspk(j),i,a%ia1(j)
           end do
         enddo

src/serial/psb_zspextrow.f90

@@ -0,0 +1,162 @@
+!!$ 
+!!$              Parallel Sparse BLAS  v2.0
+!!$    (C) Copyright 2006 Salvatore Filippone    University of Rome Tor Vergata
+!!$                       Alfredo Buttari        University of Rome Tor Vergata
+!!$ 
+!!$  Redistribution and use in source and binary forms, with or without
+!!$  modification, are permitted provided that the following conditions
+!!$  are met:
+!!$    1. Redistributions of source code must retain the above copyright
+!!$       notice, this list of conditions and the following disclaimer.
+!!$    2. Redistributions in binary form must reproduce the above copyright
+!!$       notice, this list of conditions, and the following disclaimer in the
+!!$       documentation and/or other materials provided with the distribution.
+!!$    3. The name of the PSBLAS group or the names of its contributors may
+!!$       not be used to endorse or promote products derived from this
+!!$       software without specific written permission.
+!!$ 
+!!$  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
+!!$  ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+!!$  TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
+!!$  PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
+!!$  BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
+!!$  CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
+!!$  SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
+!!$  INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
+!!$  CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
+!!$  ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+!!$  POSSIBILITY OF SUCH DAMAGE.
+!!$ 
+!!$   
+! File:  psb_zspextrow.f90 
+! Subroutine: psb_zspextrow
+!    Gets one or more rows from a sparse matrix. 
+! Parameters:
+
+!*****************************************************************************
+!*                                                                           *
+!* Takes a specified row from matrix A and copies into NZ,IA,JA,VAL  in COO  *
+!* format.                                                                   *
+!*                                                                           *
+!*****************************************************************************
+subroutine psb_zspextrow(irw,a,nz,ia,ja,val,info,iren,lrw)
+  use psb_spmat_type
+  use psb_string_mod
+  type(psb_zspmat_type), intent(in) :: a
+  integer, intent(in)       :: irw
+  integer, intent(out)      :: nz
+  integer, intent(inout)    :: ia(:), ja(:)
+  complex(kind(1.d0)),  intent(inout)    :: val(:)
+  integer, intent(in), target, optional :: iren(:)
+  integer, intent(in), optional :: lrw
+  integer, intent(out)  :: info
+  interface psb_spgtrow
+    subroutine psb_zspgtrow(irw,a,b,info,append,iren,lrw)
+      ! Output is always in  COO format  into B, irrespective of 
+      ! the input format 
+      use psb_spmat_type
+      use psb_const_mod
+      implicit none
+
+      type(psb_zspmat_type), intent(in)     :: a
+      integer, intent(in)                   :: irw
+      type(psb_zspmat_type), intent(inout)  :: b
+      integer,intent(out)                   :: info
+      logical, intent(in), optional         :: append
+      integer, intent(in), target, optional :: iren(:)
+      integer, intent(in), optional         :: lrw
+    end subroutine psb_zspgtrow
+  end interface
+
+  integer               :: lrw_, ierr(2), err_act
+  type(psb_zspmat_type) :: b
+  integer, pointer      :: iren_(:)
+  character(len=20)     :: name, ch_err
+
+
+  name='psb_sp_extrow'
+  info  = 0
+  call psb_erractionsave(err_act)
+  call psb_set_erraction(0)  
+
+  call psb_nullify_sp(b)
+
+  if (present(lrw)) then
+    lrw_ = lrw
+  else
+    lrw_ = irw
+  endif
+  if (lrw_ < irw) then
+    write(0,*) 'SPEXTROW input error: fixing lrw',irw,lrw_
+    lrw_ = irw
+  end if
+  call psb_sp_all(lrw_-irw+1,lrw_-irw+1,b,info)
+
+  if (present(iren)) then
+    call psb_spgtrow(irw,a,b,info,iren=iren,lrw=lrw_)
+  else 
+    call psb_spgtrow(irw,a,b,info,lrw=lrw_)
+  end if
+  if (info /= 0) then     
+    info=136
+    ch_err=a%fida(1:3)
+    call psb_errpush(info,name,a_err=ch_err)
+    goto 9999
+  end if
+
+  if (toupper(b%fida) /= 'COO') then 
+    info=4010
+    ch_err=a%fida(1:3)
+    call psb_errpush(info,name,a_err=ch_err)
+    goto 9999
+  endif
+
+  nz = b%infoa(psb_nnz_)
+
+  if (size(ia)>= nz) then 
+    ia(1:nz) = b%ia1(1:nz)
+  else
+    info    = 135
+    ierr(1) = 4
+    ierr(2) = size(ia)
+    call psb_errpush(info,name,i_err=ierr)
+    goto 9999
+  endif
+
+  if (size(ja)>= nz) then 
+    ja(1:nz) = b%ia2(1:nz)
+  else
+    info    = 135
+    ierr(1) = 5
+    ierr(2) = size(ja)
+    call psb_errpush(info,name,i_err=ierr)
+    goto 9999
+  endif
+
+  if (size(val)>= nz) then 
+    val(1:nz) = b%aspk(1:nz)
+  else
+    info    = 135
+    ierr(1) = 6
+    ierr(2) = size(val)
+    call psb_errpush(info,name,i_err=ierr)
+    goto 9999
+  endif
+
+
+  call psb_sp_free(b,info)
+
+  call psb_erractionrestore(err_act)
+  return
+
+9999 continue
+  call psb_erractionrestore(err_act)
+  if (err_act.eq.act_abort) then
+    call psb_error()
+    return
+  end if
+  return
+
+
+end subroutine psb_zspextrow
+

test/Fileread/mat_dist.f90

@@ -122,7 +122,7 @@ contains
     integer                     :: np, iam 
     integer                     :: ircode, length_row, i_count, j_count,&
          & k_count, blockdim, root, liwork, nrow, ncol, nnzero, nrhs,&
-         & i,j,k, ll, isize, iproc, nnr, err, err_act, int_err(5)
+         & i,j,k, ll, nz, isize, iproc, nnr, err, err_act, int_err(5)
     integer, pointer            :: iwork(:)
     character                   :: afmt*5, atyp*5
     integer, allocatable          :: irow(:),icol(:)
@@ -249,13 +249,17 @@ contains
         if (iam == root) then
 
           ll = 0
-          do j = i_count, j_count
+          do i= i_count, j_count-1
-            do k=a_glob%ia2(j),a_glob%ia2(j+1)-1
+            call psb_sp_extrow(i,a_glob,nz,&
-              ll       = ll+1
+                 & irow(ll+1:),icol(ll+1:),val(ll+1:), info)
-              irow(ll) = j
+            if (info /= 0) then            
-              icol(ll) = a_glob%ia1(k)
+              if (nz >min(size(irow(ll+1:)),size(icol(ll+1:)),size(val(ll+1:)))) then 
-              val(ll)  = a_glob%aspk(k)
+                write(0,*) 'Allocation failure? This should not happen!'
-            end do
+              end if
+              call psb_errpush(info,name,a_err=ch_err)
+              goto 9999
+            end if
+            ll = ll + nz
           end do
 
           if (iproc == iam) then
@@ -333,13 +337,17 @@ contains
           if (iam == root) then
 
             ll = 0
-            do j = i_count, i_count
+            do i= i_count, i_count
-              do k=a_glob%ia2(j),a_glob%ia2(j+1)-1
+              call psb_sp_extrow(i,a_glob,nz,&
-                ll       = ll+1
+                   & irow(ll+1:),icol(ll+1:),val(ll+1:), info)
-                irow(ll) = j
+              if (info /= 0) then            
-                icol(ll) = a_glob%ia1(k)
+                if (nz >min(size(irow(ll+1:)),size(icol(ll+1:)),size(val(ll+1:)))) then 
-                val(ll)  = a_glob%aspk(k)
+                  write(0,*) 'Allocation failure? This should not happen!'
-              end do
+                end if
+                call psb_errpush(info,name,a_err=ch_err)
+                goto 9999
+              end if
+              ll = ll + nz
             end do
 
             if (k_count == iam) then
@@ -550,7 +558,7 @@ contains
     integer                     :: np, iam
     integer                     :: ircode, length_row, i_count, j_count,&
          & k_count, blockdim, root, liwork, nrow, ncol, nnzero, nrhs,&
-         & i,j,k, ll, isize, iproc, nnr, err, err_act, int_err(5)
+         & i,j,k, ll, nz, isize, iproc, nnr, err, err_act, int_err(5)
     integer, pointer            :: iwork(:)
     character                   :: afmt*5, atyp*5
     integer, allocatable          :: irow(:),icol(:)
@@ -661,25 +669,19 @@ contains
       nnr = j_count - i_count
 
       if (iam == root) then
-        ll = a_glob%ia2(j_count)-a_glob%ia2(i_count)
+        
-        if (ll > size(val)) then 
+        ll = 0
-          deallocate(val,irow,icol)
+        do i= i_count, j_count-1
-          allocate(val(ll),irow(ll),icol(ll),stat=info)
+          call psb_sp_extrow(i,a_glob,nz,&
+               & irow(ll+1:),icol(ll+1:),val(ll+1:), info)
           if (info /= 0) then            
-            info=4010
+            if (nz >min(size(irow(ll+1:)),size(icol(ll+1:)),size(val(ll+1:)))) then 
-            ch_err='Allocate'
+              write(0,*) 'Allocation failure? This should not happen!'
+            end if
             call psb_errpush(info,name,a_err=ch_err)
             goto 9999
           end if
-
+          ll = ll + nz
-        endif
-        k = a_glob%ia2(i_count)
-        do i= i_count, j_count-1
-          do j = a_glob%ia2(i),a_glob%ia2(i+1)-1
-            irow(j-k+1)  = i
-            icol(j-k+1)  = a_glob%ia1(j) 
-            val(j-k+1) = a_glob%aspk(j)
-          end do
         end do
 
         if (iproc == iam) then
@@ -897,7 +899,7 @@ contains
     integer                     :: np, iam
     integer                     :: ircode, length_row, i_count, j_count,&
          & k_count, blockdim, root, liwork, nrow, ncol, nnzero, nrhs,&
-         & i,j,k, ll, isize, iproc, nnr, err, err_act, int_err(5)
+         & i,j,k, ll, nz, isize, iproc, nnr, err, err_act, int_err(5)
     integer, pointer            :: iwork(:)
     character                   :: afmt*5, atyp*5
     integer, allocatable          :: irow(:),icol(:)
@@ -1022,13 +1024,17 @@ contains
         if (iam == root) then
 
           ll = 0
-          do j = i_count, j_count
+          do i= i_count, j_count-1
-            do k=a_glob%ia2(j),a_glob%ia2(j+1)-1
+            call psb_sp_extrow(i,a_glob,nz,&
-              ll       = ll+1
+                 & irow(ll+1:),icol(ll+1:),val(ll+1:), info)
-              irow(ll) = j
+            if (info /= 0) then            
-              icol(ll) = a_glob%ia1(k)
+              if (nz >min(size(irow(ll+1:)),size(icol(ll+1:)),size(val(ll+1:)))) then 
-              val(ll)  = a_glob%aspk(k)
+                write(0,*) 'Allocation failure? This should not happen!'
-            end do
+              end if
+              call psb_errpush(info,name,a_err=ch_err)
+              goto 9999
+            end if
+            ll = ll + nz
           end do
 
           if (iproc == iam) then
@@ -1106,13 +1112,17 @@ contains
           if (iam == root) then
 
             ll = 0
-            do j = i_count, i_count
+            do i= i_count, i_count
-              do k=a_glob%ia2(j),a_glob%ia2(j+1)-1
+              call psb_sp_extrow(i,a_glob,nz,&
-                ll       = ll+1
+                   & irow(ll+1:),icol(ll+1:),val(ll+1:), info)
-                irow(ll) = j
+              if (info /= 0) then            
-                icol(ll) = a_glob%ia1(k)
+                if (nz >min(size(irow(ll+1:)),size(icol(ll+1:)),size(val(ll+1:)))) then 
-                val(ll)  = a_glob%aspk(k)
+                  write(0,*) 'Allocation failure? This should not happen!'
-              end do
+                end if
+                call psb_errpush(info,name,a_err=ch_err)
+                goto 9999
+              end if
+              ll = ll + nz
             end do
 
             if (k_count == iam) then
@@ -1322,7 +1332,7 @@ contains
     integer                     :: np, iam
     integer                     :: ircode, length_row, i_count, j_count,&
          & k_count, blockdim, root, liwork, nrow, ncol, nnzero, nrhs,&
-         & i,j,k, ll, isize, iproc, nnr, err, err_act, int_err(5)
+         & i,j,k, ll, nz, isize, iproc, nnr, err, err_act, int_err(5)
     integer, pointer            :: iwork(:)
     character                   :: afmt*5, atyp*5
     integer, allocatable          :: irow(:),icol(:)
@@ -1433,25 +1443,18 @@ contains
       nnr = j_count - i_count
 
       if (iam == root) then
-        ll = a_glob%ia2(j_count)-a_glob%ia2(i_count)
+        ll = 0
-        if (ll > size(val)) then 
+        do i= i_count, j_count-1
-          deallocate(val,irow,icol)
+          call psb_sp_extrow(i,a_glob,nz,&
-          allocate(val(ll),irow(ll),icol(ll),stat=info)
+               & irow(ll+1:),icol(ll+1:),val(ll+1:), info)
           if (info /= 0) then            
-            info=4010
+            if (nz >min(size(irow(ll+1:)),size(icol(ll+1:)),size(val(ll+1:)))) then 
-            ch_err='Allocate'
+              write(0,*) 'Allocation failure? This should not happen!'
+            end if
             call psb_errpush(info,name,a_err=ch_err)
             goto 9999
           end if
-
+          ll = ll + nz
-        endif
-        k = a_glob%ia2(i_count)
-        do i= i_count, j_count-1
-          do j = a_glob%ia2(i),a_glob%ia2(i+1)-1
-            irow(j-k+1)  = i
-            icol(j-k+1)  = a_glob%ia1(j) 
-            val(j-k+1) = a_glob%aspk(j)
-          end do
         end do
 
         if (iproc == iam) then