Merge branch 'development'

Changelog

@@ -1,5 +1,10 @@
 Changelog. A lot less detailed than usual, at least for past
 	history.
+2018/10/10: new ICTXT argument to prec%init method.
+2018/09/04: Modified vector class get_vect method
+2018/08/10: Optional arguments in GETROW method. 	
+2018/07/30: Improved TRIL/TRIU implementations. 	
+2018/06/14: New FCG code. 	
 2018/04/23: Change default for CDALL with VL. New GLOBAL argument for
 	reductions. 	
 2018/04/15: Fixed pargen benchmark programs. Made MOLD mandatory.

base/modules/psb_const_mod.F90

@@ -91,10 +91,10 @@ module psb_const_mod
   ! 
   ! Version
   !
-  character(len=*), parameter    :: psb_version_string_ = "3.5.0"
+  character(len=*), parameter    :: psb_version_string_ = "3.6.0"
   integer(psb_ipk_), parameter   :: psb_version_major_  = 3
-  integer(psb_ipk_), parameter   :: psb_version_minor_  = 5
+  integer(psb_ipk_), parameter   :: psb_version_minor_  = 6
-  integer(psb_ipk_), parameter   :: psb_patchlevel_     = 1
+  integer(psb_ipk_), parameter   :: psb_patchlevel_     = 0
 
   !
   !     Handy & miscellaneous constants

base/modules/serial/psb_c_base_mat_mod.f90

@@ -316,7 +316,7 @@ module psb_c_base_mat_mod
   !
   interface 
     subroutine psb_c_base_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_c_base_sparse_mat, psb_spk_
       class(psb_c_base_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -327,7 +327,7 @@ module psb_c_base_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale, chksz
     end subroutine psb_c_base_csgetrow
   end interface
   
@@ -355,7 +355,7 @@ module psb_c_base_mat_mod
   !
   interface 
     subroutine psb_c_base_csgetblk(imin,imax,a,b,info,&
-         & jmin,jmax,iren,append,rscale,cscale)
+         & jmin,jmax,iren,append,rscale,cscale,chksz)
       import :: psb_ipk_, psb_c_base_sparse_mat, psb_c_coo_sparse_mat, psb_spk_
       class(psb_c_base_sparse_mat), intent(in) :: a
       class(psb_c_coo_sparse_mat), intent(inout) :: b
@@ -364,7 +364,7 @@ module psb_c_base_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_c_base_csgetblk
   end interface
   
@@ -1550,7 +1550,7 @@ module psb_c_base_mat_mod
   !! \see psb_c_base_mat_mod::psb_c_base_csgetrow
   interface 
     subroutine psb_c_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_c_coo_sparse_mat, psb_spk_
       class(psb_c_coo_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -1561,7 +1561,7 @@ module psb_c_base_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_c_coo_csgetrow
   end interface
   

base/modules/serial/psb_c_base_vect_mod.f90

@@ -688,19 +688,24 @@ contains
   !! \brief  Extract a copy of the contents
   !!
   !    
-  function  c_base_get_vect(x) result(res)
+  function  c_base_get_vect(x,n) result(res)
     class(psb_c_base_vect_type), intent(inout) :: x
     complex(psb_spk_), allocatable                 :: res(:)
     integer(psb_ipk_) :: info
+    integer(psb_ipk_), optional :: n
+    ! Local variables
+    integer(psb_ipk_) :: isz
     
     if (.not.allocated(x%v)) return 
     if (.not.x%is_host()) call x%sync()
-    allocate(res(x%get_nrows()),stat=info) 
+    isz = x%get_nrows()
+    if (present(n)) isz = max(0,min(isz,n))
+    allocate(res(isz),stat=info) 
     if (info /= 0) then 
       call psb_errpush(psb_err_alloc_dealloc_,'base_get_vect')
       return
     end if
-    res(:) = x%v(:)
+    res(1:isz) = x%v(1:isz)
   end function c_base_get_vect
     
   !
@@ -744,9 +749,10 @@ contains
         
     integer(psb_ipk_) :: info, first_, last_, nr
 
-    first_=1
+
-    last_=min(psb_size(x%v),size(val))
+    first_                     = 1
     if (present(first)) first_ = max(1,first)
+    last_                      = min(psb_size(x%v),first_+size(val)-1)
     if (present(last))  last_  = min(last,last_)
 
     if (allocated(x%v)) then 

base/modules/serial/psb_c_csc_mat_mod.f90

@@ -313,7 +313,7 @@ module psb_c_csc_mat_mod
   !! \see psb_c_base_mat_mod::psb_c_base_csgetrow
   interface 
     subroutine psb_c_csc_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_c_csc_sparse_mat, psb_spk_
       class(psb_c_csc_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -324,27 +324,27 @@ module psb_c_csc_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_c_csc_csgetrow
   end interface
 
-  !> \memberof psb_c_csc_sparse_mat
+!!$  !> \memberof psb_c_csc_sparse_mat
-  !! \see psb_c_base_mat_mod::psb_c_base_csgetblk
+!!$  !! \see psb_c_base_mat_mod::psb_c_base_csgetblk
-  interface 
+!!$  interface 
-    subroutine psb_c_csc_csgetblk(imin,imax,a,b,info,&
+!!$    subroutine psb_c_csc_csgetblk(imin,imax,a,b,info,&
-       & jmin,jmax,iren,append,rscale,cscale)
+!!$       & jmin,jmax,iren,append,rscale,cscale,chksz)
-      import :: psb_ipk_, psb_c_csc_sparse_mat, psb_spk_, psb_c_coo_sparse_mat
+!!$      import :: psb_ipk_, psb_c_csc_sparse_mat, psb_spk_, psb_c_coo_sparse_mat
-      class(psb_c_csc_sparse_mat), intent(in) :: a
+!!$      class(psb_c_csc_sparse_mat), intent(in) :: a
-      class(psb_c_coo_sparse_mat), intent(inout) :: b
+!!$      class(psb_c_coo_sparse_mat), intent(inout) :: b
-      integer(psb_ipk_), intent(in)                  :: imin,imax
+!!$      integer(psb_ipk_), intent(in)                  :: imin,imax
-      integer(psb_ipk_),intent(out)                  :: info
+!!$      integer(psb_ipk_),intent(out)                  :: info
-      logical, intent(in), optional        :: append
+!!$      logical, intent(in), optional        :: append
-      integer(psb_ipk_), intent(in), optional        :: iren(:)
+!!$      integer(psb_ipk_), intent(in), optional        :: iren(:)
-      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+!!$      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-      logical, intent(in), optional        :: rscale,cscale
+!!$      logical, intent(in), optional        :: rscale,cscale,chksz
-    end subroutine psb_c_csc_csgetblk
+!!$    end subroutine psb_c_csc_csgetblk
-  end interface
+!!$  end interface
-    
+!!$    
   !> \memberof psb_c_csc_sparse_mat
   !! \see psb_c_base_mat_mod::psb_c_base_cssv
   interface 

base/modules/serial/psb_c_csr_mat_mod.f90

@@ -80,6 +80,8 @@ module psb_c_csr_mat_mod
     procedure, pass(a) :: aclsum      => psb_c_csr_aclsum
     procedure, pass(a) :: reallocate_nz => psb_c_csr_reallocate_nz
     procedure, pass(a) :: allocate_mnnz => psb_c_csr_allocate_mnnz
+    procedure, pass(a) :: tril          => psb_c_csr_tril
+    procedure, pass(a) :: triu          => psb_c_csr_triu
     procedure, pass(a) :: cp_to_coo   => psb_c_cp_csr_to_coo
     procedure, pass(a) :: cp_from_coo => psb_c_cp_csr_from_coo
     procedure, pass(a) :: cp_to_fmt   => psb_c_cp_csr_to_fmt
@@ -170,6 +172,93 @@ module psb_c_csr_mat_mod
       integer(psb_ipk_), intent(in), optional     :: ivr(:), ivc(:)
     end subroutine psb_c_csr_print
   end interface
+  !
+  !> Function  tril:
+  !! \memberof  psb_c_base_sparse_mat
+  !! \brief  Copy the lower triangle, i.e. all entries
+  !!         A(I,J) such that J-I <= DIAG
+  !!         default value is DIAG=0, i.e. lower triangle up to
+  !!         the main diagonal.
+  !!         DIAG=-1 means copy the strictly lower triangle
+  !!         DIAG= 1 means copy the lower triangle plus the first diagonal
+  !!                 of the upper triangle.
+  !!         Moreover, apply a clipping by copying entries A(I,J) only if
+  !!         IMIN<=I<=IMAX
+  !!         JMIN<=J<=JMAX
+  !!         
+  !!  \param l     the output (sub)matrix
+  !!  \param info  return code
+  !!  \param diag [0] the last diagonal (J-I) to be considered.
+  !!  \param imin [1] the minimum row index we are interested in 
+  !!  \param imax [a\%get_nrows()] the minimum row index we are interested in 
+  !!  \param jmin [1] minimum col index 
+  !!  \param jmax [a\%get_ncols()] maximum col index 
+  !!  \param iren(:) [none] an array to return renumbered indices (iren(ia(:)),iren(ja(:))
+  !!  \param rscale [false] map [min(ia(:)):max(ia(:))] onto [1:max(ia(:))-min(ia(:))+1]
+  !!  \param cscale [false] map [min(ja(:)):max(ja(:))] onto [1:max(ja(:))-min(ja(:))+1]
+  !!          ( iren cannot be specified with rscale/cscale)
+  !!  \param append [false] append to ia,ja 
+  !!  \param nzin [none]  if append, then first new entry should go in entry nzin+1
+  !!  \param u  [none]  copy of the complementary triangle
+  !!           
+  !
+  interface 
+    subroutine psb_c_csr_tril(a,l,info,diag,imin,imax,&
+         & jmin,jmax,rscale,cscale,u)
+      import :: psb_ipk_, psb_c_csr_sparse_mat, psb_c_coo_sparse_mat, psb_spk_
+      class(psb_c_csr_sparse_mat), intent(in) :: a
+      class(psb_c_coo_sparse_mat), intent(out) :: l
+      integer(psb_ipk_),intent(out)              :: info
+      integer(psb_ipk_), intent(in), optional    :: diag,imin,imax,jmin,jmax
+      logical, intent(in), optional              :: rscale,cscale
+      class(psb_c_coo_sparse_mat), optional, intent(out) :: u
+    end subroutine psb_c_csr_tril
+  end interface
+  
+  !
+  !> Function  triu:
+  !! \memberof  psb_c_csr_sparse_mat
+  !! \brief  Copy the upper triangle, i.e. all entries
+  !!         A(I,J) such that DIAG <= J-I
+  !!         default value is DIAG=0, i.e. upper triangle from 
+  !!         the main diagonal up.
+  !!         DIAG= 1 means copy the strictly upper triangle
+  !!         DIAG=-1 means copy the upper triangle plus the first diagonal
+  !!                 of the lower triangle.
+  !!         Moreover, apply a clipping by copying entries A(I,J) only if
+  !!         IMIN<=I<=IMAX
+  !!         JMIN<=J<=JMAX
+  !!         Optionally copies the lower triangle at the same time
+  !!         
+  !!  \param u     the output (sub)matrix
+  !!  \param info  return code
+  !!  \param diag [0] the last diagonal (J-I) to be considered.
+  !!  \param imin [1] the minimum row index we are interested in 
+  !!  \param imax [a\%get_nrows()] the minimum row index we are interested in 
+  !!  \param jmin [1] minimum col index 
+  !!  \param jmax [a\%get_ncols()] maximum col index 
+  !!  \param iren(:) [none] an array to return renumbered indices (iren(ia(:)),iren(ja(:))
+  !!  \param rscale [false] map [min(ia(:)):max(ia(:))] onto [1:max(ia(:))-min(ia(:))+1]
+  !!  \param cscale [false] map [min(ja(:)):max(ja(:))] onto [1:max(ja(:))-min(ja(:))+1]
+  !!          ( iren cannot be specified with rscale/cscale)
+  !!  \param append [false] append to ia,ja 
+  !!  \param nzin [none]  if append, then first new entry should go in entry nzin+1
+  !!  \param l  [none]  copy of the complementary triangle
+  !!           
+  !
+  interface 
+    subroutine psb_c_csr_triu(a,u,info,diag,imin,imax,&
+         & jmin,jmax,rscale,cscale,l)
+      import :: psb_ipk_, psb_c_csr_sparse_mat, psb_c_coo_sparse_mat, psb_spk_
+      class(psb_c_csr_sparse_mat), intent(in) :: a
+      class(psb_c_coo_sparse_mat), intent(out) :: u
+      integer(psb_ipk_),intent(out)              :: info
+      integer(psb_ipk_), intent(in), optional    :: diag,imin,imax,jmin,jmax
+      logical, intent(in), optional              :: rscale,cscale
+      class(psb_c_coo_sparse_mat), optional, intent(out) :: l
+    end subroutine psb_c_csr_triu
+  end interface
+  
   
   !> \memberof psb_c_csr_sparse_mat
   !! \see psb_c_base_mat_mod::psb_c_base_cp_to_coo
@@ -316,7 +405,7 @@ module psb_c_csr_mat_mod
   !! \see psb_c_base_mat_mod::psb_c_base_csgetrow
   interface 
     subroutine psb_c_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_c_csr_sparse_mat, psb_spk_
       class(psb_c_csr_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -327,26 +416,26 @@ module psb_c_csr_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_c_csr_csgetrow
   end interface
-
+!!$
-  !> \memberof psb_c_csr_sparse_mat
+!!$  !> \memberof psb_c_csr_sparse_mat
-  !! \see psb_c_base_mat_mod::psb_c_base_csgetblk
+!!$  !! \see psb_c_base_mat_mod::psb_c_base_csgetblk
-  interface 
+!!$  interface 
-    subroutine psb_c_csr_csgetblk(imin,imax,a,b,info,&
+!!$    subroutine psb_c_csr_csgetblk(imin,imax,a,b,info,&
-       & jmin,jmax,iren,append,rscale,cscale)
+!!$       & jmin,jmax,iren,append,rscale,cscale)
-      import :: psb_ipk_, psb_c_csr_sparse_mat, psb_spk_, psb_c_coo_sparse_mat
+!!$      import :: psb_ipk_, psb_c_csr_sparse_mat, psb_spk_, psb_c_coo_sparse_mat
-      class(psb_c_csr_sparse_mat), intent(in) :: a
+!!$      class(psb_c_csr_sparse_mat), intent(in) :: a
-      class(psb_c_coo_sparse_mat), intent(inout) :: b
+!!$      class(psb_c_coo_sparse_mat), intent(inout) :: b
-      integer(psb_ipk_), intent(in)                  :: imin,imax
+!!$      integer(psb_ipk_), intent(in)                  :: imin,imax
-      integer(psb_ipk_),intent(out)                  :: info
+!!$      integer(psb_ipk_),intent(out)                  :: info
-      logical, intent(in), optional        :: append
+!!$      logical, intent(in), optional        :: append
-      integer(psb_ipk_), intent(in), optional        :: iren(:)
+!!$      integer(psb_ipk_), intent(in), optional        :: iren(:)
-      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+!!$      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-      logical, intent(in), optional        :: rscale,cscale
+!!$      logical, intent(in), optional        :: rscale,cscale
-    end subroutine psb_c_csr_csgetblk
+!!$    end subroutine psb_c_csr_csgetblk
-  end interface
+!!$  end interface
     
   !> \memberof psb_c_csr_sparse_mat
   !! \see psb_c_base_mat_mod::psb_c_base_cssv

base/modules/serial/psb_c_mat_mod.f90

@@ -450,7 +450,7 @@ module psb_c_mat_mod
   
   interface 
     subroutine psb_c_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_cspmat_type, psb_spk_
       class(psb_cspmat_type), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -461,7 +461,7 @@ module psb_c_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_c_csgetrow
   end interface
   

base/modules/serial/psb_c_vect_mod.F90

@@ -227,13 +227,14 @@ contains
 
   end subroutine c_vect_bld_n
 
-  function  c_vect_get_vect(x) result(res)
+  function  c_vect_get_vect(x,n) result(res)
     class(psb_c_vect_type), intent(inout)  :: x
     complex(psb_spk_), allocatable                 :: res(:)
     integer(psb_ipk_) :: info
+    integer(psb_ipk_), optional :: n
 
     if (allocated(x%v)) then
-      res = x%v%get_vect()
+      res = x%v%get_vect(n)
     end if
   end function c_vect_get_vect
 

base/modules/serial/psb_d_base_mat_mod.f90

@@ -316,7 +316,7 @@ module psb_d_base_mat_mod
   !
   interface 
     subroutine psb_d_base_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_d_base_sparse_mat, psb_dpk_
       class(psb_d_base_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -327,7 +327,7 @@ module psb_d_base_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale, chksz
     end subroutine psb_d_base_csgetrow
   end interface
   
@@ -355,7 +355,7 @@ module psb_d_base_mat_mod
   !
   interface 
     subroutine psb_d_base_csgetblk(imin,imax,a,b,info,&
-         & jmin,jmax,iren,append,rscale,cscale)
+         & jmin,jmax,iren,append,rscale,cscale,chksz)
       import :: psb_ipk_, psb_d_base_sparse_mat, psb_d_coo_sparse_mat, psb_dpk_
       class(psb_d_base_sparse_mat), intent(in) :: a
       class(psb_d_coo_sparse_mat), intent(inout) :: b
@@ -364,7 +364,7 @@ module psb_d_base_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_d_base_csgetblk
   end interface
   
@@ -1550,7 +1550,7 @@ module psb_d_base_mat_mod
   !! \see psb_d_base_mat_mod::psb_d_base_csgetrow
   interface 
     subroutine psb_d_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_d_coo_sparse_mat, psb_dpk_
       class(psb_d_coo_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -1561,7 +1561,7 @@ module psb_d_base_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_d_coo_csgetrow
   end interface
   

base/modules/serial/psb_d_base_vect_mod.f90

@@ -688,19 +688,24 @@ contains
   !! \brief  Extract a copy of the contents
   !!
   !    
-  function  d_base_get_vect(x) result(res)
+  function  d_base_get_vect(x,n) result(res)
     class(psb_d_base_vect_type), intent(inout) :: x
     real(psb_dpk_), allocatable                 :: res(:)
     integer(psb_ipk_) :: info
+    integer(psb_ipk_), optional :: n
+    ! Local variables
+    integer(psb_ipk_) :: isz
     
     if (.not.allocated(x%v)) return 
     if (.not.x%is_host()) call x%sync()
-    allocate(res(x%get_nrows()),stat=info) 
+    isz = x%get_nrows()
+    if (present(n)) isz = max(0,min(isz,n))
+    allocate(res(isz),stat=info) 
     if (info /= 0) then 
       call psb_errpush(psb_err_alloc_dealloc_,'base_get_vect')
       return
     end if
-    res(:) = x%v(:)
+    res(1:isz) = x%v(1:isz)
   end function d_base_get_vect
     
   !
@@ -744,9 +749,10 @@ contains
         
     integer(psb_ipk_) :: info, first_, last_, nr
 
-    first_=1
+
-    last_=min(psb_size(x%v),size(val))
+    first_                     = 1
     if (present(first)) first_ = max(1,first)
+    last_                      = min(psb_size(x%v),first_+size(val)-1)
     if (present(last))  last_  = min(last,last_)
 
     if (allocated(x%v)) then 

base/modules/serial/psb_d_csc_mat_mod.f90

@@ -313,7 +313,7 @@ module psb_d_csc_mat_mod
   !! \see psb_d_base_mat_mod::psb_d_base_csgetrow
   interface 
     subroutine psb_d_csc_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_d_csc_sparse_mat, psb_dpk_
       class(psb_d_csc_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -324,27 +324,27 @@ module psb_d_csc_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_d_csc_csgetrow
   end interface
 
-  !> \memberof psb_d_csc_sparse_mat
+!!$  !> \memberof psb_d_csc_sparse_mat
-  !! \see psb_d_base_mat_mod::psb_d_base_csgetblk
+!!$  !! \see psb_d_base_mat_mod::psb_d_base_csgetblk
-  interface 
+!!$  interface 
-    subroutine psb_d_csc_csgetblk(imin,imax,a,b,info,&
+!!$    subroutine psb_d_csc_csgetblk(imin,imax,a,b,info,&
-       & jmin,jmax,iren,append,rscale,cscale)
+!!$       & jmin,jmax,iren,append,rscale,cscale,chksz)
-      import :: psb_ipk_, psb_d_csc_sparse_mat, psb_dpk_, psb_d_coo_sparse_mat
+!!$      import :: psb_ipk_, psb_d_csc_sparse_mat, psb_dpk_, psb_d_coo_sparse_mat
-      class(psb_d_csc_sparse_mat), intent(in) :: a
+!!$      class(psb_d_csc_sparse_mat), intent(in) :: a
-      class(psb_d_coo_sparse_mat), intent(inout) :: b
+!!$      class(psb_d_coo_sparse_mat), intent(inout) :: b
-      integer(psb_ipk_), intent(in)                  :: imin,imax
+!!$      integer(psb_ipk_), intent(in)                  :: imin,imax
-      integer(psb_ipk_),intent(out)                  :: info
+!!$      integer(psb_ipk_),intent(out)                  :: info
-      logical, intent(in), optional        :: append
+!!$      logical, intent(in), optional        :: append
-      integer(psb_ipk_), intent(in), optional        :: iren(:)
+!!$      integer(psb_ipk_), intent(in), optional        :: iren(:)
-      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+!!$      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-      logical, intent(in), optional        :: rscale,cscale
+!!$      logical, intent(in), optional        :: rscale,cscale,chksz
-    end subroutine psb_d_csc_csgetblk
+!!$    end subroutine psb_d_csc_csgetblk
-  end interface
+!!$  end interface
-    
+!!$    
   !> \memberof psb_d_csc_sparse_mat
   !! \see psb_d_base_mat_mod::psb_d_base_cssv
   interface 

base/modules/serial/psb_d_csr_mat_mod.f90

@@ -80,6 +80,8 @@ module psb_d_csr_mat_mod
     procedure, pass(a) :: aclsum      => psb_d_csr_aclsum
     procedure, pass(a) :: reallocate_nz => psb_d_csr_reallocate_nz
     procedure, pass(a) :: allocate_mnnz => psb_d_csr_allocate_mnnz
+    procedure, pass(a) :: tril          => psb_d_csr_tril
+    procedure, pass(a) :: triu          => psb_d_csr_triu
     procedure, pass(a) :: cp_to_coo   => psb_d_cp_csr_to_coo
     procedure, pass(a) :: cp_from_coo => psb_d_cp_csr_from_coo
     procedure, pass(a) :: cp_to_fmt   => psb_d_cp_csr_to_fmt
@@ -170,6 +172,93 @@ module psb_d_csr_mat_mod
       integer(psb_ipk_), intent(in), optional     :: ivr(:), ivc(:)
     end subroutine psb_d_csr_print
   end interface
+  !
+  !> Function  tril:
+  !! \memberof  psb_d_base_sparse_mat
+  !! \brief  Copy the lower triangle, i.e. all entries
+  !!         A(I,J) such that J-I <= DIAG
+  !!         default value is DIAG=0, i.e. lower triangle up to
+  !!         the main diagonal.
+  !!         DIAG=-1 means copy the strictly lower triangle
+  !!         DIAG= 1 means copy the lower triangle plus the first diagonal
+  !!                 of the upper triangle.
+  !!         Moreover, apply a clipping by copying entries A(I,J) only if
+  !!         IMIN<=I<=IMAX
+  !!         JMIN<=J<=JMAX
+  !!         
+  !!  \param l     the output (sub)matrix
+  !!  \param info  return code
+  !!  \param diag [0] the last diagonal (J-I) to be considered.
+  !!  \param imin [1] the minimum row index we are interested in 
+  !!  \param imax [a\%get_nrows()] the minimum row index we are interested in 
+  !!  \param jmin [1] minimum col index 
+  !!  \param jmax [a\%get_ncols()] maximum col index 
+  !!  \param iren(:) [none] an array to return renumbered indices (iren(ia(:)),iren(ja(:))
+  !!  \param rscale [false] map [min(ia(:)):max(ia(:))] onto [1:max(ia(:))-min(ia(:))+1]
+  !!  \param cscale [false] map [min(ja(:)):max(ja(:))] onto [1:max(ja(:))-min(ja(:))+1]
+  !!          ( iren cannot be specified with rscale/cscale)
+  !!  \param append [false] append to ia,ja 
+  !!  \param nzin [none]  if append, then first new entry should go in entry nzin+1
+  !!  \param u  [none]  copy of the complementary triangle
+  !!           
+  !
+  interface 
+    subroutine psb_d_csr_tril(a,l,info,diag,imin,imax,&
+         & jmin,jmax,rscale,cscale,u)
+      import :: psb_ipk_, psb_d_csr_sparse_mat, psb_d_coo_sparse_mat, psb_dpk_
+      class(psb_d_csr_sparse_mat), intent(in) :: a
+      class(psb_d_coo_sparse_mat), intent(out) :: l
+      integer(psb_ipk_),intent(out)              :: info
+      integer(psb_ipk_), intent(in), optional    :: diag,imin,imax,jmin,jmax
+      logical, intent(in), optional              :: rscale,cscale
+      class(psb_d_coo_sparse_mat), optional, intent(out) :: u
+    end subroutine psb_d_csr_tril
+  end interface
+  
+  !
+  !> Function  triu:
+  !! \memberof  psb_d_csr_sparse_mat
+  !! \brief  Copy the upper triangle, i.e. all entries
+  !!         A(I,J) such that DIAG <= J-I
+  !!         default value is DIAG=0, i.e. upper triangle from 
+  !!         the main diagonal up.
+  !!         DIAG= 1 means copy the strictly upper triangle
+  !!         DIAG=-1 means copy the upper triangle plus the first diagonal
+  !!                 of the lower triangle.
+  !!         Moreover, apply a clipping by copying entries A(I,J) only if
+  !!         IMIN<=I<=IMAX
+  !!         JMIN<=J<=JMAX
+  !!         Optionally copies the lower triangle at the same time
+  !!         
+  !!  \param u     the output (sub)matrix
+  !!  \param info  return code
+  !!  \param diag [0] the last diagonal (J-I) to be considered.
+  !!  \param imin [1] the minimum row index we are interested in 
+  !!  \param imax [a\%get_nrows()] the minimum row index we are interested in 
+  !!  \param jmin [1] minimum col index 
+  !!  \param jmax [a\%get_ncols()] maximum col index 
+  !!  \param iren(:) [none] an array to return renumbered indices (iren(ia(:)),iren(ja(:))
+  !!  \param rscale [false] map [min(ia(:)):max(ia(:))] onto [1:max(ia(:))-min(ia(:))+1]
+  !!  \param cscale [false] map [min(ja(:)):max(ja(:))] onto [1:max(ja(:))-min(ja(:))+1]
+  !!          ( iren cannot be specified with rscale/cscale)
+  !!  \param append [false] append to ia,ja 
+  !!  \param nzin [none]  if append, then first new entry should go in entry nzin+1
+  !!  \param l  [none]  copy of the complementary triangle
+  !!           
+  !
+  interface 
+    subroutine psb_d_csr_triu(a,u,info,diag,imin,imax,&
+         & jmin,jmax,rscale,cscale,l)
+      import :: psb_ipk_, psb_d_csr_sparse_mat, psb_d_coo_sparse_mat, psb_dpk_
+      class(psb_d_csr_sparse_mat), intent(in) :: a
+      class(psb_d_coo_sparse_mat), intent(out) :: u
+      integer(psb_ipk_),intent(out)              :: info
+      integer(psb_ipk_), intent(in), optional    :: diag,imin,imax,jmin,jmax
+      logical, intent(in), optional              :: rscale,cscale
+      class(psb_d_coo_sparse_mat), optional, intent(out) :: l
+    end subroutine psb_d_csr_triu
+  end interface
+  
   
   !> \memberof psb_d_csr_sparse_mat
   !! \see psb_d_base_mat_mod::psb_d_base_cp_to_coo
@@ -316,7 +405,7 @@ module psb_d_csr_mat_mod
   !! \see psb_d_base_mat_mod::psb_d_base_csgetrow
   interface 
     subroutine psb_d_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_d_csr_sparse_mat, psb_dpk_
       class(psb_d_csr_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -327,26 +416,26 @@ module psb_d_csr_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_d_csr_csgetrow
   end interface
-
+!!$
-  !> \memberof psb_d_csr_sparse_mat
+!!$  !> \memberof psb_d_csr_sparse_mat
-  !! \see psb_d_base_mat_mod::psb_d_base_csgetblk
+!!$  !! \see psb_d_base_mat_mod::psb_d_base_csgetblk
-  interface 
+!!$  interface 
-    subroutine psb_d_csr_csgetblk(imin,imax,a,b,info,&
+!!$    subroutine psb_d_csr_csgetblk(imin,imax,a,b,info,&
-       & jmin,jmax,iren,append,rscale,cscale)
+!!$       & jmin,jmax,iren,append,rscale,cscale)
-      import :: psb_ipk_, psb_d_csr_sparse_mat, psb_dpk_, psb_d_coo_sparse_mat
+!!$      import :: psb_ipk_, psb_d_csr_sparse_mat, psb_dpk_, psb_d_coo_sparse_mat
-      class(psb_d_csr_sparse_mat), intent(in) :: a
+!!$      class(psb_d_csr_sparse_mat), intent(in) :: a
-      class(psb_d_coo_sparse_mat), intent(inout) :: b
+!!$      class(psb_d_coo_sparse_mat), intent(inout) :: b
-      integer(psb_ipk_), intent(in)                  :: imin,imax
+!!$      integer(psb_ipk_), intent(in)                  :: imin,imax
-      integer(psb_ipk_),intent(out)                  :: info
+!!$      integer(psb_ipk_),intent(out)                  :: info
-      logical, intent(in), optional        :: append
+!!$      logical, intent(in), optional        :: append
-      integer(psb_ipk_), intent(in), optional        :: iren(:)
+!!$      integer(psb_ipk_), intent(in), optional        :: iren(:)
-      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+!!$      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-      logical, intent(in), optional        :: rscale,cscale
+!!$      logical, intent(in), optional        :: rscale,cscale
-    end subroutine psb_d_csr_csgetblk
+!!$    end subroutine psb_d_csr_csgetblk
-  end interface
+!!$  end interface
     
   !> \memberof psb_d_csr_sparse_mat
   !! \see psb_d_base_mat_mod::psb_d_base_cssv

base/modules/serial/psb_d_mat_mod.f90

@@ -450,7 +450,7 @@ module psb_d_mat_mod
   
   interface 
     subroutine psb_d_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_dspmat_type, psb_dpk_
       class(psb_dspmat_type), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -461,7 +461,7 @@ module psb_d_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_d_csgetrow
   end interface
   

base/modules/serial/psb_d_vect_mod.F90

@@ -227,13 +227,14 @@ contains
 
   end subroutine d_vect_bld_n
 
-  function  d_vect_get_vect(x) result(res)
+  function  d_vect_get_vect(x,n) result(res)
     class(psb_d_vect_type), intent(inout)  :: x
     real(psb_dpk_), allocatable                 :: res(:)
     integer(psb_ipk_) :: info
+    integer(psb_ipk_), optional :: n
 
     if (allocated(x%v)) then
-      res = x%v%get_vect()
+      res = x%v%get_vect(n)
     end if
   end function d_vect_get_vect
 

base/modules/serial/psb_i_base_vect_mod.f90

@@ -657,19 +657,24 @@ contains
   !! \brief  Extract a copy of the contents
   !!
   !    
-  function  i_base_get_vect(x) result(res)
+  function  i_base_get_vect(x,n) result(res)
     class(psb_i_base_vect_type), intent(inout) :: x
     integer(psb_ipk_), allocatable                 :: res(:)
     integer(psb_ipk_) :: info
+    integer(psb_ipk_), optional :: n
+    ! Local variables
+    integer(psb_ipk_) :: isz
     
     if (.not.allocated(x%v)) return 
     if (.not.x%is_host()) call x%sync()
-    allocate(res(x%get_nrows()),stat=info) 
+    isz = x%get_nrows()
+    if (present(n)) isz = max(0,min(isz,n))
+    allocate(res(isz),stat=info) 
     if (info /= 0) then 
       call psb_errpush(psb_err_alloc_dealloc_,'base_get_vect')
       return
     end if
-    res(:) = x%v(:)
+    res(1:isz) = x%v(1:isz)
   end function i_base_get_vect
     
   !
@@ -713,9 +718,10 @@ contains
         
     integer(psb_ipk_) :: info, first_, last_, nr
 
-    first_=1
+
-    last_=min(psb_size(x%v),size(val))
+    first_                     = 1
     if (present(first)) first_ = max(1,first)
+    last_                      = min(psb_size(x%v),first_+size(val)-1)
     if (present(last))  last_  = min(last,last_)
 
     if (allocated(x%v)) then 

base/modules/serial/psb_i_vect_mod.F90

@@ -200,13 +200,14 @@ contains
 
   end subroutine i_vect_bld_n
 
-  function  i_vect_get_vect(x) result(res)
+  function  i_vect_get_vect(x,n) result(res)
     class(psb_i_vect_type), intent(inout)  :: x
     integer(psb_ipk_), allocatable                 :: res(:)
     integer(psb_ipk_) :: info
+    integer(psb_ipk_), optional :: n
 
     if (allocated(x%v)) then
-      res = x%v%get_vect()
+      res = x%v%get_vect(n)
     end if
   end function i_vect_get_vect
 

base/modules/serial/psb_s_base_mat_mod.f90

@@ -316,7 +316,7 @@ module psb_s_base_mat_mod
   !
   interface 
     subroutine psb_s_base_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_s_base_sparse_mat, psb_spk_
       class(psb_s_base_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -327,7 +327,7 @@ module psb_s_base_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale, chksz
     end subroutine psb_s_base_csgetrow
   end interface
   
@@ -355,7 +355,7 @@ module psb_s_base_mat_mod
   !
   interface 
     subroutine psb_s_base_csgetblk(imin,imax,a,b,info,&
-         & jmin,jmax,iren,append,rscale,cscale)
+         & jmin,jmax,iren,append,rscale,cscale,chksz)
       import :: psb_ipk_, psb_s_base_sparse_mat, psb_s_coo_sparse_mat, psb_spk_
       class(psb_s_base_sparse_mat), intent(in) :: a
       class(psb_s_coo_sparse_mat), intent(inout) :: b
@@ -364,7 +364,7 @@ module psb_s_base_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_s_base_csgetblk
   end interface
   
@@ -1550,7 +1550,7 @@ module psb_s_base_mat_mod
   !! \see psb_s_base_mat_mod::psb_s_base_csgetrow
   interface 
     subroutine psb_s_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_s_coo_sparse_mat, psb_spk_
       class(psb_s_coo_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -1561,7 +1561,7 @@ module psb_s_base_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_s_coo_csgetrow
   end interface
   

base/modules/serial/psb_s_base_vect_mod.f90

@@ -688,19 +688,24 @@ contains
   !! \brief  Extract a copy of the contents
   !!
   !    
-  function  s_base_get_vect(x) result(res)
+  function  s_base_get_vect(x,n) result(res)
     class(psb_s_base_vect_type), intent(inout) :: x
     real(psb_spk_), allocatable                 :: res(:)
     integer(psb_ipk_) :: info
+    integer(psb_ipk_), optional :: n
+    ! Local variables
+    integer(psb_ipk_) :: isz
     
     if (.not.allocated(x%v)) return 
     if (.not.x%is_host()) call x%sync()
-    allocate(res(x%get_nrows()),stat=info) 
+    isz = x%get_nrows()
+    if (present(n)) isz = max(0,min(isz,n))
+    allocate(res(isz),stat=info) 
     if (info /= 0) then 
       call psb_errpush(psb_err_alloc_dealloc_,'base_get_vect')
       return
     end if
-    res(:) = x%v(:)
+    res(1:isz) = x%v(1:isz)
   end function s_base_get_vect
     
   !
@@ -744,9 +749,10 @@ contains
         
     integer(psb_ipk_) :: info, first_, last_, nr
 
-    first_=1
+
-    last_=min(psb_size(x%v),size(val))
+    first_                     = 1
     if (present(first)) first_ = max(1,first)
+    last_                      = min(psb_size(x%v),first_+size(val)-1)
     if (present(last))  last_  = min(last,last_)
 
     if (allocated(x%v)) then 

base/modules/serial/psb_s_csc_mat_mod.f90

@@ -313,7 +313,7 @@ module psb_s_csc_mat_mod
   !! \see psb_s_base_mat_mod::psb_s_base_csgetrow
   interface 
     subroutine psb_s_csc_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_s_csc_sparse_mat, psb_spk_
       class(psb_s_csc_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -324,27 +324,27 @@ module psb_s_csc_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_s_csc_csgetrow
   end interface
 
-  !> \memberof psb_s_csc_sparse_mat
+!!$  !> \memberof psb_s_csc_sparse_mat
-  !! \see psb_s_base_mat_mod::psb_s_base_csgetblk
+!!$  !! \see psb_s_base_mat_mod::psb_s_base_csgetblk
-  interface 
+!!$  interface 
-    subroutine psb_s_csc_csgetblk(imin,imax,a,b,info,&
+!!$    subroutine psb_s_csc_csgetblk(imin,imax,a,b,info,&
-       & jmin,jmax,iren,append,rscale,cscale)
+!!$       & jmin,jmax,iren,append,rscale,cscale,chksz)
-      import :: psb_ipk_, psb_s_csc_sparse_mat, psb_spk_, psb_s_coo_sparse_mat
+!!$      import :: psb_ipk_, psb_s_csc_sparse_mat, psb_spk_, psb_s_coo_sparse_mat
-      class(psb_s_csc_sparse_mat), intent(in) :: a
+!!$      class(psb_s_csc_sparse_mat), intent(in) :: a
-      class(psb_s_coo_sparse_mat), intent(inout) :: b
+!!$      class(psb_s_coo_sparse_mat), intent(inout) :: b
-      integer(psb_ipk_), intent(in)                  :: imin,imax
+!!$      integer(psb_ipk_), intent(in)                  :: imin,imax
-      integer(psb_ipk_),intent(out)                  :: info
+!!$      integer(psb_ipk_),intent(out)                  :: info
-      logical, intent(in), optional        :: append
+!!$      logical, intent(in), optional        :: append
-      integer(psb_ipk_), intent(in), optional        :: iren(:)
+!!$      integer(psb_ipk_), intent(in), optional        :: iren(:)
-      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+!!$      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-      logical, intent(in), optional        :: rscale,cscale
+!!$      logical, intent(in), optional        :: rscale,cscale,chksz
-    end subroutine psb_s_csc_csgetblk
+!!$    end subroutine psb_s_csc_csgetblk
-  end interface
+!!$  end interface
-    
+!!$    
   !> \memberof psb_s_csc_sparse_mat
   !! \see psb_s_base_mat_mod::psb_s_base_cssv
   interface 

base/modules/serial/psb_s_csr_mat_mod.f90

@@ -80,6 +80,8 @@ module psb_s_csr_mat_mod
     procedure, pass(a) :: aclsum      => psb_s_csr_aclsum
     procedure, pass(a) :: reallocate_nz => psb_s_csr_reallocate_nz
     procedure, pass(a) :: allocate_mnnz => psb_s_csr_allocate_mnnz
+    procedure, pass(a) :: tril          => psb_s_csr_tril
+    procedure, pass(a) :: triu          => psb_s_csr_triu
     procedure, pass(a) :: cp_to_coo   => psb_s_cp_csr_to_coo
     procedure, pass(a) :: cp_from_coo => psb_s_cp_csr_from_coo
     procedure, pass(a) :: cp_to_fmt   => psb_s_cp_csr_to_fmt
@@ -170,6 +172,93 @@ module psb_s_csr_mat_mod
       integer(psb_ipk_), intent(in), optional     :: ivr(:), ivc(:)
     end subroutine psb_s_csr_print
   end interface
+  !
+  !> Function  tril:
+  !! \memberof  psb_s_base_sparse_mat
+  !! \brief  Copy the lower triangle, i.e. all entries
+  !!         A(I,J) such that J-I <= DIAG
+  !!         default value is DIAG=0, i.e. lower triangle up to
+  !!         the main diagonal.
+  !!         DIAG=-1 means copy the strictly lower triangle
+  !!         DIAG= 1 means copy the lower triangle plus the first diagonal
+  !!                 of the upper triangle.
+  !!         Moreover, apply a clipping by copying entries A(I,J) only if
+  !!         IMIN<=I<=IMAX
+  !!         JMIN<=J<=JMAX
+  !!         
+  !!  \param l     the output (sub)matrix
+  !!  \param info  return code
+  !!  \param diag [0] the last diagonal (J-I) to be considered.
+  !!  \param imin [1] the minimum row index we are interested in 
+  !!  \param imax [a\%get_nrows()] the minimum row index we are interested in 
+  !!  \param jmin [1] minimum col index 
+  !!  \param jmax [a\%get_ncols()] maximum col index 
+  !!  \param iren(:) [none] an array to return renumbered indices (iren(ia(:)),iren(ja(:))
+  !!  \param rscale [false] map [min(ia(:)):max(ia(:))] onto [1:max(ia(:))-min(ia(:))+1]
+  !!  \param cscale [false] map [min(ja(:)):max(ja(:))] onto [1:max(ja(:))-min(ja(:))+1]
+  !!          ( iren cannot be specified with rscale/cscale)
+  !!  \param append [false] append to ia,ja 
+  !!  \param nzin [none]  if append, then first new entry should go in entry nzin+1
+  !!  \param u  [none]  copy of the complementary triangle
+  !!           
+  !
+  interface 
+    subroutine psb_s_csr_tril(a,l,info,diag,imin,imax,&
+         & jmin,jmax,rscale,cscale,u)
+      import :: psb_ipk_, psb_s_csr_sparse_mat, psb_s_coo_sparse_mat, psb_spk_
+      class(psb_s_csr_sparse_mat), intent(in) :: a
+      class(psb_s_coo_sparse_mat), intent(out) :: l
+      integer(psb_ipk_),intent(out)              :: info
+      integer(psb_ipk_), intent(in), optional    :: diag,imin,imax,jmin,jmax
+      logical, intent(in), optional              :: rscale,cscale
+      class(psb_s_coo_sparse_mat), optional, intent(out) :: u
+    end subroutine psb_s_csr_tril
+  end interface
+  
+  !
+  !> Function  triu:
+  !! \memberof  psb_s_csr_sparse_mat
+  !! \brief  Copy the upper triangle, i.e. all entries
+  !!         A(I,J) such that DIAG <= J-I
+  !!         default value is DIAG=0, i.e. upper triangle from 
+  !!         the main diagonal up.
+  !!         DIAG= 1 means copy the strictly upper triangle
+  !!         DIAG=-1 means copy the upper triangle plus the first diagonal
+  !!                 of the lower triangle.
+  !!         Moreover, apply a clipping by copying entries A(I,J) only if
+  !!         IMIN<=I<=IMAX
+  !!         JMIN<=J<=JMAX
+  !!         Optionally copies the lower triangle at the same time
+  !!         
+  !!  \param u     the output (sub)matrix
+  !!  \param info  return code
+  !!  \param diag [0] the last diagonal (J-I) to be considered.
+  !!  \param imin [1] the minimum row index we are interested in 
+  !!  \param imax [a\%get_nrows()] the minimum row index we are interested in 
+  !!  \param jmin [1] minimum col index 
+  !!  \param jmax [a\%get_ncols()] maximum col index 
+  !!  \param iren(:) [none] an array to return renumbered indices (iren(ia(:)),iren(ja(:))
+  !!  \param rscale [false] map [min(ia(:)):max(ia(:))] onto [1:max(ia(:))-min(ia(:))+1]
+  !!  \param cscale [false] map [min(ja(:)):max(ja(:))] onto [1:max(ja(:))-min(ja(:))+1]
+  !!          ( iren cannot be specified with rscale/cscale)
+  !!  \param append [false] append to ia,ja 
+  !!  \param nzin [none]  if append, then first new entry should go in entry nzin+1
+  !!  \param l  [none]  copy of the complementary triangle
+  !!           
+  !
+  interface 
+    subroutine psb_s_csr_triu(a,u,info,diag,imin,imax,&
+         & jmin,jmax,rscale,cscale,l)
+      import :: psb_ipk_, psb_s_csr_sparse_mat, psb_s_coo_sparse_mat, psb_spk_
+      class(psb_s_csr_sparse_mat), intent(in) :: a
+      class(psb_s_coo_sparse_mat), intent(out) :: u
+      integer(psb_ipk_),intent(out)              :: info
+      integer(psb_ipk_), intent(in), optional    :: diag,imin,imax,jmin,jmax
+      logical, intent(in), optional              :: rscale,cscale
+      class(psb_s_coo_sparse_mat), optional, intent(out) :: l
+    end subroutine psb_s_csr_triu
+  end interface
+  
   
   !> \memberof psb_s_csr_sparse_mat
   !! \see psb_s_base_mat_mod::psb_s_base_cp_to_coo
@@ -316,7 +405,7 @@ module psb_s_csr_mat_mod
   !! \see psb_s_base_mat_mod::psb_s_base_csgetrow
   interface 
     subroutine psb_s_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_s_csr_sparse_mat, psb_spk_
       class(psb_s_csr_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -327,26 +416,26 @@ module psb_s_csr_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_s_csr_csgetrow
   end interface
-
+!!$
-  !> \memberof psb_s_csr_sparse_mat
+!!$  !> \memberof psb_s_csr_sparse_mat
-  !! \see psb_s_base_mat_mod::psb_s_base_csgetblk
+!!$  !! \see psb_s_base_mat_mod::psb_s_base_csgetblk
-  interface 
+!!$  interface 
-    subroutine psb_s_csr_csgetblk(imin,imax,a,b,info,&
+!!$    subroutine psb_s_csr_csgetblk(imin,imax,a,b,info,&
-       & jmin,jmax,iren,append,rscale,cscale)
+!!$       & jmin,jmax,iren,append,rscale,cscale)
-      import :: psb_ipk_, psb_s_csr_sparse_mat, psb_spk_, psb_s_coo_sparse_mat
+!!$      import :: psb_ipk_, psb_s_csr_sparse_mat, psb_spk_, psb_s_coo_sparse_mat
-      class(psb_s_csr_sparse_mat), intent(in) :: a
+!!$      class(psb_s_csr_sparse_mat), intent(in) :: a
-      class(psb_s_coo_sparse_mat), intent(inout) :: b
+!!$      class(psb_s_coo_sparse_mat), intent(inout) :: b
-      integer(psb_ipk_), intent(in)                  :: imin,imax
+!!$      integer(psb_ipk_), intent(in)                  :: imin,imax
-      integer(psb_ipk_),intent(out)                  :: info
+!!$      integer(psb_ipk_),intent(out)                  :: info
-      logical, intent(in), optional        :: append
+!!$      logical, intent(in), optional        :: append
-      integer(psb_ipk_), intent(in), optional        :: iren(:)
+!!$      integer(psb_ipk_), intent(in), optional        :: iren(:)
-      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+!!$      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-      logical, intent(in), optional        :: rscale,cscale
+!!$      logical, intent(in), optional        :: rscale,cscale
-    end subroutine psb_s_csr_csgetblk
+!!$    end subroutine psb_s_csr_csgetblk
-  end interface
+!!$  end interface
     
   !> \memberof psb_s_csr_sparse_mat
   !! \see psb_s_base_mat_mod::psb_s_base_cssv

base/modules/serial/psb_s_mat_mod.f90

@@ -450,7 +450,7 @@ module psb_s_mat_mod
   
   interface 
     subroutine psb_s_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_sspmat_type, psb_spk_
       class(psb_sspmat_type), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -461,7 +461,7 @@ module psb_s_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_s_csgetrow
   end interface
   

base/modules/serial/psb_s_vect_mod.F90

@@ -227,13 +227,14 @@ contains
 
   end subroutine s_vect_bld_n
 
-  function  s_vect_get_vect(x) result(res)
+  function  s_vect_get_vect(x,n) result(res)
     class(psb_s_vect_type), intent(inout)  :: x
     real(psb_spk_), allocatable                 :: res(:)
     integer(psb_ipk_) :: info
+    integer(psb_ipk_), optional :: n
 
     if (allocated(x%v)) then
-      res = x%v%get_vect()
+      res = x%v%get_vect(n)
     end if
   end function s_vect_get_vect
 

base/modules/serial/psb_z_base_mat_mod.f90

@@ -316,7 +316,7 @@ module psb_z_base_mat_mod
   !
   interface 
     subroutine psb_z_base_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_z_base_sparse_mat, psb_dpk_
       class(psb_z_base_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -327,7 +327,7 @@ module psb_z_base_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale, chksz
     end subroutine psb_z_base_csgetrow
   end interface
   
@@ -355,7 +355,7 @@ module psb_z_base_mat_mod
   !
   interface 
     subroutine psb_z_base_csgetblk(imin,imax,a,b,info,&
-         & jmin,jmax,iren,append,rscale,cscale)
+         & jmin,jmax,iren,append,rscale,cscale,chksz)
       import :: psb_ipk_, psb_z_base_sparse_mat, psb_z_coo_sparse_mat, psb_dpk_
       class(psb_z_base_sparse_mat), intent(in) :: a
       class(psb_z_coo_sparse_mat), intent(inout) :: b
@@ -364,7 +364,7 @@ module psb_z_base_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_z_base_csgetblk
   end interface
   
@@ -1550,7 +1550,7 @@ module psb_z_base_mat_mod
   !! \see psb_z_base_mat_mod::psb_z_base_csgetrow
   interface 
     subroutine psb_z_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_z_coo_sparse_mat, psb_dpk_
       class(psb_z_coo_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -1561,7 +1561,7 @@ module psb_z_base_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_z_coo_csgetrow
   end interface
   

base/modules/serial/psb_z_base_vect_mod.f90

@@ -688,19 +688,24 @@ contains
   !! \brief  Extract a copy of the contents
   !!
   !    
-  function  z_base_get_vect(x) result(res)
+  function  z_base_get_vect(x,n) result(res)
     class(psb_z_base_vect_type), intent(inout) :: x
     complex(psb_dpk_), allocatable                 :: res(:)
     integer(psb_ipk_) :: info
+    integer(psb_ipk_), optional :: n
+    ! Local variables
+    integer(psb_ipk_) :: isz
     
     if (.not.allocated(x%v)) return 
     if (.not.x%is_host()) call x%sync()
-    allocate(res(x%get_nrows()),stat=info) 
+    isz = x%get_nrows()
+    if (present(n)) isz = max(0,min(isz,n))
+    allocate(res(isz),stat=info) 
     if (info /= 0) then 
       call psb_errpush(psb_err_alloc_dealloc_,'base_get_vect')
       return
     end if
-    res(:) = x%v(:)
+    res(1:isz) = x%v(1:isz)
   end function z_base_get_vect
     
   !
@@ -744,9 +749,10 @@ contains
         
     integer(psb_ipk_) :: info, first_, last_, nr
 
-    first_=1
+
-    last_=min(psb_size(x%v),size(val))
+    first_                     = 1
     if (present(first)) first_ = max(1,first)
+    last_                      = min(psb_size(x%v),first_+size(val)-1)
     if (present(last))  last_  = min(last,last_)
 
     if (allocated(x%v)) then 

base/modules/serial/psb_z_csc_mat_mod.f90

@@ -313,7 +313,7 @@ module psb_z_csc_mat_mod
   !! \see psb_z_base_mat_mod::psb_z_base_csgetrow
   interface 
     subroutine psb_z_csc_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_z_csc_sparse_mat, psb_dpk_
       class(psb_z_csc_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -324,27 +324,27 @@ module psb_z_csc_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_z_csc_csgetrow
   end interface
 
-  !> \memberof psb_z_csc_sparse_mat
+!!$  !> \memberof psb_z_csc_sparse_mat
-  !! \see psb_z_base_mat_mod::psb_z_base_csgetblk
+!!$  !! \see psb_z_base_mat_mod::psb_z_base_csgetblk
-  interface 
+!!$  interface 
-    subroutine psb_z_csc_csgetblk(imin,imax,a,b,info,&
+!!$    subroutine psb_z_csc_csgetblk(imin,imax,a,b,info,&
-       & jmin,jmax,iren,append,rscale,cscale)
+!!$       & jmin,jmax,iren,append,rscale,cscale,chksz)
-      import :: psb_ipk_, psb_z_csc_sparse_mat, psb_dpk_, psb_z_coo_sparse_mat
+!!$      import :: psb_ipk_, psb_z_csc_sparse_mat, psb_dpk_, psb_z_coo_sparse_mat
-      class(psb_z_csc_sparse_mat), intent(in) :: a
+!!$      class(psb_z_csc_sparse_mat), intent(in) :: a
-      class(psb_z_coo_sparse_mat), intent(inout) :: b
+!!$      class(psb_z_coo_sparse_mat), intent(inout) :: b
-      integer(psb_ipk_), intent(in)                  :: imin,imax
+!!$      integer(psb_ipk_), intent(in)                  :: imin,imax
-      integer(psb_ipk_),intent(out)                  :: info
+!!$      integer(psb_ipk_),intent(out)                  :: info
-      logical, intent(in), optional        :: append
+!!$      logical, intent(in), optional        :: append
-      integer(psb_ipk_), intent(in), optional        :: iren(:)
+!!$      integer(psb_ipk_), intent(in), optional        :: iren(:)
-      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+!!$      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-      logical, intent(in), optional        :: rscale,cscale
+!!$      logical, intent(in), optional        :: rscale,cscale,chksz
-    end subroutine psb_z_csc_csgetblk
+!!$    end subroutine psb_z_csc_csgetblk
-  end interface
+!!$  end interface
-    
+!!$    
   !> \memberof psb_z_csc_sparse_mat
   !! \see psb_z_base_mat_mod::psb_z_base_cssv
   interface 

base/modules/serial/psb_z_csr_mat_mod.f90

@@ -80,6 +80,8 @@ module psb_z_csr_mat_mod
     procedure, pass(a) :: aclsum      => psb_z_csr_aclsum
     procedure, pass(a) :: reallocate_nz => psb_z_csr_reallocate_nz
     procedure, pass(a) :: allocate_mnnz => psb_z_csr_allocate_mnnz
+    procedure, pass(a) :: tril          => psb_z_csr_tril
+    procedure, pass(a) :: triu          => psb_z_csr_triu
     procedure, pass(a) :: cp_to_coo   => psb_z_cp_csr_to_coo
     procedure, pass(a) :: cp_from_coo => psb_z_cp_csr_from_coo
     procedure, pass(a) :: cp_to_fmt   => psb_z_cp_csr_to_fmt
@@ -170,6 +172,93 @@ module psb_z_csr_mat_mod
       integer(psb_ipk_), intent(in), optional     :: ivr(:), ivc(:)
     end subroutine psb_z_csr_print
   end interface
+  !
+  !> Function  tril:
+  !! \memberof  psb_z_base_sparse_mat
+  !! \brief  Copy the lower triangle, i.e. all entries
+  !!         A(I,J) such that J-I <= DIAG
+  !!         default value is DIAG=0, i.e. lower triangle up to
+  !!         the main diagonal.
+  !!         DIAG=-1 means copy the strictly lower triangle
+  !!         DIAG= 1 means copy the lower triangle plus the first diagonal
+  !!                 of the upper triangle.
+  !!         Moreover, apply a clipping by copying entries A(I,J) only if
+  !!         IMIN<=I<=IMAX
+  !!         JMIN<=J<=JMAX
+  !!         
+  !!  \param l     the output (sub)matrix
+  !!  \param info  return code
+  !!  \param diag [0] the last diagonal (J-I) to be considered.
+  !!  \param imin [1] the minimum row index we are interested in 
+  !!  \param imax [a\%get_nrows()] the minimum row index we are interested in 
+  !!  \param jmin [1] minimum col index 
+  !!  \param jmax [a\%get_ncols()] maximum col index 
+  !!  \param iren(:) [none] an array to return renumbered indices (iren(ia(:)),iren(ja(:))
+  !!  \param rscale [false] map [min(ia(:)):max(ia(:))] onto [1:max(ia(:))-min(ia(:))+1]
+  !!  \param cscale [false] map [min(ja(:)):max(ja(:))] onto [1:max(ja(:))-min(ja(:))+1]
+  !!          ( iren cannot be specified with rscale/cscale)
+  !!  \param append [false] append to ia,ja 
+  !!  \param nzin [none]  if append, then first new entry should go in entry nzin+1
+  !!  \param u  [none]  copy of the complementary triangle
+  !!           
+  !
+  interface 
+    subroutine psb_z_csr_tril(a,l,info,diag,imin,imax,&
+         & jmin,jmax,rscale,cscale,u)
+      import :: psb_ipk_, psb_z_csr_sparse_mat, psb_z_coo_sparse_mat, psb_dpk_
+      class(psb_z_csr_sparse_mat), intent(in) :: a
+      class(psb_z_coo_sparse_mat), intent(out) :: l
+      integer(psb_ipk_),intent(out)              :: info
+      integer(psb_ipk_), intent(in), optional    :: diag,imin,imax,jmin,jmax
+      logical, intent(in), optional              :: rscale,cscale
+      class(psb_z_coo_sparse_mat), optional, intent(out) :: u
+    end subroutine psb_z_csr_tril
+  end interface
+  
+  !
+  !> Function  triu:
+  !! \memberof  psb_z_csr_sparse_mat
+  !! \brief  Copy the upper triangle, i.e. all entries
+  !!         A(I,J) such that DIAG <= J-I
+  !!         default value is DIAG=0, i.e. upper triangle from 
+  !!         the main diagonal up.
+  !!         DIAG= 1 means copy the strictly upper triangle
+  !!         DIAG=-1 means copy the upper triangle plus the first diagonal
+  !!                 of the lower triangle.
+  !!         Moreover, apply a clipping by copying entries A(I,J) only if
+  !!         IMIN<=I<=IMAX
+  !!         JMIN<=J<=JMAX
+  !!         Optionally copies the lower triangle at the same time
+  !!         
+  !!  \param u     the output (sub)matrix
+  !!  \param info  return code
+  !!  \param diag [0] the last diagonal (J-I) to be considered.
+  !!  \param imin [1] the minimum row index we are interested in 
+  !!  \param imax [a\%get_nrows()] the minimum row index we are interested in 
+  !!  \param jmin [1] minimum col index 
+  !!  \param jmax [a\%get_ncols()] maximum col index 
+  !!  \param iren(:) [none] an array to return renumbered indices (iren(ia(:)),iren(ja(:))
+  !!  \param rscale [false] map [min(ia(:)):max(ia(:))] onto [1:max(ia(:))-min(ia(:))+1]
+  !!  \param cscale [false] map [min(ja(:)):max(ja(:))] onto [1:max(ja(:))-min(ja(:))+1]
+  !!          ( iren cannot be specified with rscale/cscale)
+  !!  \param append [false] append to ia,ja 
+  !!  \param nzin [none]  if append, then first new entry should go in entry nzin+1
+  !!  \param l  [none]  copy of the complementary triangle
+  !!           
+  !
+  interface 
+    subroutine psb_z_csr_triu(a,u,info,diag,imin,imax,&
+         & jmin,jmax,rscale,cscale,l)
+      import :: psb_ipk_, psb_z_csr_sparse_mat, psb_z_coo_sparse_mat, psb_dpk_
+      class(psb_z_csr_sparse_mat), intent(in) :: a
+      class(psb_z_coo_sparse_mat), intent(out) :: u
+      integer(psb_ipk_),intent(out)              :: info
+      integer(psb_ipk_), intent(in), optional    :: diag,imin,imax,jmin,jmax
+      logical, intent(in), optional              :: rscale,cscale
+      class(psb_z_coo_sparse_mat), optional, intent(out) :: l
+    end subroutine psb_z_csr_triu
+  end interface
+  
   
   !> \memberof psb_z_csr_sparse_mat
   !! \see psb_z_base_mat_mod::psb_z_base_cp_to_coo
@@ -316,7 +405,7 @@ module psb_z_csr_mat_mod
   !! \see psb_z_base_mat_mod::psb_z_base_csgetrow
   interface 
     subroutine psb_z_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_z_csr_sparse_mat, psb_dpk_
       class(psb_z_csr_sparse_mat), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -327,26 +416,26 @@ module psb_z_csr_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_z_csr_csgetrow
   end interface
-
+!!$
-  !> \memberof psb_z_csr_sparse_mat
+!!$  !> \memberof psb_z_csr_sparse_mat
-  !! \see psb_z_base_mat_mod::psb_z_base_csgetblk
+!!$  !! \see psb_z_base_mat_mod::psb_z_base_csgetblk
-  interface 
+!!$  interface 
-    subroutine psb_z_csr_csgetblk(imin,imax,a,b,info,&
+!!$    subroutine psb_z_csr_csgetblk(imin,imax,a,b,info,&
-       & jmin,jmax,iren,append,rscale,cscale)
+!!$       & jmin,jmax,iren,append,rscale,cscale)
-      import :: psb_ipk_, psb_z_csr_sparse_mat, psb_dpk_, psb_z_coo_sparse_mat
+!!$      import :: psb_ipk_, psb_z_csr_sparse_mat, psb_dpk_, psb_z_coo_sparse_mat
-      class(psb_z_csr_sparse_mat), intent(in) :: a
+!!$      class(psb_z_csr_sparse_mat), intent(in) :: a
-      class(psb_z_coo_sparse_mat), intent(inout) :: b
+!!$      class(psb_z_coo_sparse_mat), intent(inout) :: b
-      integer(psb_ipk_), intent(in)                  :: imin,imax
+!!$      integer(psb_ipk_), intent(in)                  :: imin,imax
-      integer(psb_ipk_),intent(out)                  :: info
+!!$      integer(psb_ipk_),intent(out)                  :: info
-      logical, intent(in), optional        :: append
+!!$      logical, intent(in), optional        :: append
-      integer(psb_ipk_), intent(in), optional        :: iren(:)
+!!$      integer(psb_ipk_), intent(in), optional        :: iren(:)
-      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+!!$      integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-      logical, intent(in), optional        :: rscale,cscale
+!!$      logical, intent(in), optional        :: rscale,cscale
-    end subroutine psb_z_csr_csgetblk
+!!$    end subroutine psb_z_csr_csgetblk
-  end interface
+!!$  end interface
     
   !> \memberof psb_z_csr_sparse_mat
   !! \see psb_z_base_mat_mod::psb_z_base_cssv

base/modules/serial/psb_z_mat_mod.f90

@@ -450,7 +450,7 @@ module psb_z_mat_mod
   
   interface 
     subroutine psb_z_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-         & jmin,jmax,iren,append,nzin,rscale,cscale)
+         & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
       import :: psb_ipk_, psb_zspmat_type, psb_dpk_
       class(psb_zspmat_type), intent(in) :: a
       integer(psb_ipk_), intent(in)                  :: imin,imax
@@ -461,7 +461,7 @@ module psb_z_mat_mod
       logical, intent(in), optional        :: append
       integer(psb_ipk_), intent(in), optional        :: iren(:)
       integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-      logical, intent(in), optional        :: rscale,cscale
+      logical, intent(in), optional        :: rscale,cscale,chksz
     end subroutine psb_z_csgetrow
   end interface
   

base/modules/serial/psb_z_vect_mod.F90

@@ -227,13 +227,14 @@ contains
 
   end subroutine z_vect_bld_n
 
-  function  z_vect_get_vect(x) result(res)
+  function  z_vect_get_vect(x,n) result(res)
     class(psb_z_vect_type), intent(inout)  :: x
     complex(psb_dpk_), allocatable                 :: res(:)
     integer(psb_ipk_) :: info
+    integer(psb_ipk_), optional :: n
 
     if (allocated(x%v)) then
-      res = x%v%get_vect()
+      res = x%v%get_vect(n)
     end if
   end function z_vect_get_vect
 

base/serial/impl/psb_c_base_mat_impl.F90

@@ -406,7 +406,7 @@ subroutine psb_c_base_csput_v(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl)
 end subroutine psb_c_base_csput_v
 
 subroutine psb_c_base_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -422,7 +422,7 @@ subroutine psb_c_base_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
   integer(psb_ipk_) :: err_act
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -447,7 +447,7 @@ end subroutine psb_c_base_csgetrow
 ! If performance is critical it can be overridden.
 !
 subroutine psb_c_base_csgetblk(imin,imax,a,b,info,&
-     & jmin,jmax,iren,append,rscale,cscale)
+     & jmin,jmax,iren,append,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -461,7 +461,7 @@ subroutine psb_c_base_csgetblk(imin,imax,a,b,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
   integer(psb_ipk_) :: err_act, nzin, nzout
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -522,7 +522,7 @@ subroutine psb_c_base_csgetblk(imin,imax,a,b,info,&
 
   call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
        & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
-       & nzin=nzin, rscale=rscale, cscale=cscale)
+       & nzin=nzin, rscale=rscale, cscale=cscale, chksz=chksz)
 
   if (info /= psb_success_) goto 9999
 
@@ -645,7 +645,7 @@ subroutine psb_c_base_tril(a,l,info,&
   logical, intent(in), optional            :: rscale,cscale
   class(psb_c_coo_sparse_mat), optional, intent(out) :: u
   
-  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k 
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k, ibk 
   integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
   integer(psb_ipk_), allocatable :: ia(:), ja(:)
   complex(psb_spk_), allocatable    :: val(:)
@@ -653,6 +653,7 @@ subroutine psb_c_base_tril(a,l,info,&
   character(len=20)  :: name='tril'
   logical :: rscale_, cscale_
   logical, parameter :: debug=.false.
+  integer(psb_ipk_), parameter :: nbk=8
 
   call psb_erractionsave(err_act)
   info = psb_success_
@@ -715,12 +716,12 @@ subroutine psb_c_base_tril(a,l,info,&
     call psb_realloc(max(mb,nb),ia,info)
     call psb_realloc(max(mb,nb),ja,info)
     call psb_realloc(max(mb,nb),val,info)
-    do i=imin_,imax_
+    do i=imin_,imax_, nbk
-      call a%csget(i,i,nzout,ia,ja,val,info,&
+      ibk = min(nbk,imax_-i+1)
+      call a%csget(i,i+ibk-1,nzout,ia,ja,val,info,&
            & jmin=jmin_, jmax=jmax_)
       do k=1, nzout
-        j = ja(k) 
+        if ((ja(k)-ia(k))<=diag_) then
-        if (j-i<=diag_) then
           nzlin = nzlin + 1
           l%ia(nzlin)  = ia(k)
           l%ja(nzlin)  = ja(k)
@@ -796,7 +797,7 @@ subroutine psb_c_base_triu(a,u,info,&
   logical, intent(in), optional            :: rscale,cscale
   class(psb_c_coo_sparse_mat), optional, intent(out) :: l
   
-  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k 
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k, ibk
   integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
   integer(psb_ipk_), allocatable :: ia(:), ja(:)
   complex(psb_spk_), allocatable    :: val(:)
@@ -804,6 +805,7 @@ subroutine psb_c_base_triu(a,u,info,&
   character(len=20)  :: name='triu'
   logical :: rscale_, cscale_
   logical, parameter :: debug=.false.
+  integer(psb_ipk_), parameter :: nbk=8
 
   call psb_erractionsave(err_act)
   info = psb_success_
@@ -866,12 +868,12 @@ subroutine psb_c_base_triu(a,u,info,&
     call psb_realloc(max(mb,nb),ia,info)
     call psb_realloc(max(mb,nb),ja,info)
     call psb_realloc(max(mb,nb),val,info)
-    do i=imin_,imax_
+    do i=imin_,imax_, nbk
-      call a%csget(i,i,nzout,ia,ja,val,info,&
+      ibk = min(nbk,imax_-i+1)
+      call a%csget(i,i+ibk-1,nzout,ia,ja,val,info,&
            & jmin=jmin_, jmax=jmax_)
       do k=1, nzout
-        j = ja(k) 
+        if ((ja(k)-ia(k))<diag_) then
-        if (j-i<diag_) then
           nzlin = nzlin + 1
           l%ia(nzlin)  = ia(k)
           l%ja(nzlin)  = ja(k)

base/serial/impl/psb_c_coo_impl.f90

@@ -2259,7 +2259,7 @@ end subroutine psb_c_coo_csgetptn
 ! The output is guaranteed to be sorted
 ! 
 subroutine psb_c_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -2276,9 +2276,9 @@ subroutine psb_c_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
-  logical :: append_, rscale_, cscale_ 
+  logical :: append_, rscale_, cscale_, chksz_
   integer(psb_ipk_) :: nzin_, jmin_, jmax_, err_act, i
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -2321,13 +2321,18 @@ subroutine psb_c_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   else
     cscale_ = .false.
   endif
+  if (present(chksz)) then 
+    chksz_ = chksz
+  else
+    chksz_ = .true.
+  endif
   if ((rscale_.or.cscale_).and.(present(iren))) then 
     info = psb_err_many_optional_arg_
     call psb_errpush(info,name,a_err='iren (rscale.or.cscale)')
     goto 9999
   end if
 
-  call coo_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,info,&
+  call coo_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,chksz_,info,&
        & iren)
 
   if (rscale_) then 
@@ -2352,7 +2357,7 @@ subroutine psb_c_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
 
 contains
 
-  subroutine coo_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,info,&
+  subroutine coo_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,chksz,info,&
        & iren)
 
     use psb_const_mod
@@ -2368,7 +2373,7 @@ contains
     integer(psb_ipk_), allocatable, intent(inout)  :: ia(:), ja(:)
     complex(psb_spk_), allocatable,  intent(inout)    :: val(:)
     integer(psb_ipk_), intent(in)                  :: nzin
-    logical, intent(in)                            :: append
+    logical, intent(in)                            :: append,chksz
     integer(psb_ipk_) :: info
     integer(psb_ipk_), optional                    :: iren(:)
     integer(psb_ipk_) :: nzin_, nza, idx,ip,jp,i,k, nzt, irw, lrw, nra, nca, nrd
@@ -2452,11 +2457,13 @@ contains
         nzt = jp - ip +1 
         nz = 0 
 
-        call psb_ensure_size(nzin_+nzt,ia,info)
+        if (chksz) then 
-        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+          call psb_ensure_size(nzin_+nzt,ia,info)
-        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+          if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-        if (info /= psb_success_) return
+          if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-
+          if (info /= psb_success_) return
+        end if
+        
         if (present(iren)) then 
           do i=ip,jp
             if ((jmin <= a%ja(i)).and.(a%ja(i)<=jmax)) then 
@@ -2488,11 +2495,13 @@ contains
 
       nrd = max(a%get_nrows(),1)
       nzt = ((nza+nrd-1)/nrd)*(lrw-irw+1)
-      call psb_ensure_size(nzin_+nzt,ia,info)
+      if (chksz) then 
-      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+        call psb_ensure_size(nzin_+nzt,ia,info)
-      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-      if (info /= psb_success_) return
+        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-
+        if (info /= psb_success_) return
+      end if
+      
       if (present(iren)) then 
         k = 0 
         do i=1, a%get_nzeros()
@@ -2501,10 +2510,12 @@ contains
             k = k + 1 
             if (k > nzt) then
               nzt = k + nzt
-              call psb_ensure_size(nzin_+nzt,ia,info)
+              if (chksz) then 
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+                call psb_ensure_size(nzin_+nzt,ia,info)
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-              if (info /= psb_success_) return
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+                if (info /= psb_success_) return
+              end if              
             end if
             val(nzin_+k) = a%val(i)
             ia(nzin_+k)  = iren(a%ia(i))
@@ -2519,11 +2530,12 @@ contains
             k = k + 1 
             if (k > nzt) then
               nzt = k + nzt
-              call psb_ensure_size(nzin_+nzt,ia,info)
+              if (chksz) then 
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+                call psb_ensure_size(nzin_+nzt,ia,info)
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-              if (info /= psb_success_) return
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-
+                if (info /= psb_success_) return
+              end if
             end if
             val(nzin_+k) = a%val(i)
             ia(nzin_+k)  = (a%ia(i))

base/serial/impl/psb_c_csc_impl.f90

@@ -1686,7 +1686,7 @@ end subroutine psb_c_csc_csgetptn
 
 
 subroutine psb_c_csc_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -1704,7 +1704,7 @@ subroutine psb_c_csc_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
   logical :: append_, rscale_, cscale_ 
   integer(psb_ipk_) :: nzin_, jmin_, jmax_, err_act, i
@@ -2557,60 +2557,60 @@ end subroutine psb_c_csc_reallocate_nz
 
 
 
-subroutine psb_c_csc_csgetblk(imin,imax,a,b,info,&
+!!$subroutine psb_c_csc_csgetblk(imin,imax,a,b,info,&
-     & jmin,jmax,iren,append,rscale,cscale)
+!!$     & jmin,jmax,iren,append,rscale,cscale)
-  ! Output is always in  COO format 
+!!$  ! Output is always in  COO format 
-  use psb_error_mod
+!!$  use psb_error_mod
-  use psb_const_mod
+!!$  use psb_const_mod
-  use psb_c_csc_mat_mod, psb_protect_name => psb_c_csc_csgetblk
+!!$  use psb_c_csc_mat_mod, psb_protect_name => psb_c_csc_csgetblk
-  implicit none
+!!$  implicit none
-
+!!$
-  class(psb_c_csc_sparse_mat), intent(in) :: a
+!!$  class(psb_c_csc_sparse_mat), intent(in) :: a
-  class(psb_c_coo_sparse_mat), intent(inout) :: b
+!!$  class(psb_c_coo_sparse_mat), intent(inout) :: b
-  integer(psb_ipk_), intent(in)                  :: imin,imax
+!!$  integer(psb_ipk_), intent(in)                  :: imin,imax
-  integer(psb_ipk_),intent(out)                  :: info
+!!$  integer(psb_ipk_),intent(out)                  :: info
-  logical, intent(in), optional        :: append
+!!$  logical, intent(in), optional        :: append
-  integer(psb_ipk_), intent(in), optional        :: iren(:)
+!!$  integer(psb_ipk_), intent(in), optional        :: iren(:)
-  integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+!!$  integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-  logical, intent(in), optional        :: rscale,cscale
+!!$  logical, intent(in), optional        :: rscale,cscale
-  integer(psb_ipk_) :: err_act, nzin, nzout
+!!$  integer(psb_ipk_) :: err_act, nzin, nzout
-  integer(psb_ipk_) :: ierr(5)
+!!$  integer(psb_ipk_) :: ierr(5)
-  character(len=20)  :: name='csget'
+!!$  character(len=20)  :: name='csget'
-  logical :: append_
+!!$  logical :: append_
-  logical, parameter :: debug=.false.
+!!$  logical, parameter :: debug=.false.
-
+!!$
-  call psb_erractionsave(err_act)
+!!$  call psb_erractionsave(err_act)
-  info = psb_success_
+!!$  info = psb_success_
-
+!!$
-  if (present(append)) then 
+!!$  if (present(append)) then 
-    append_ = append
+!!$    append_ = append
-  else
+!!$  else
-    append_ = .false.
+!!$    append_ = .false.
-  endif
+!!$  endif
-  if (append_) then 
+!!$  if (append_) then 
-    nzin = a%get_nzeros()
+!!$    nzin = a%get_nzeros()
-  else
+!!$  else
-    nzin = 0
+!!$    nzin = 0
-  endif
+!!$  endif
-
+!!$
-  call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
+!!$  call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
-       & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
+!!$       & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
-       & nzin=nzin, rscale=rscale, cscale=cscale)
+!!$       & nzin=nzin, rscale=rscale, cscale=cscale)
-
+!!$
-  if (info /= psb_success_) goto 9999
+!!$  if (info /= psb_success_) goto 9999
-
+!!$
-  call b%set_nzeros(nzin+nzout)
+!!$  call b%set_nzeros(nzin+nzout)
-  call b%fix(info)
+!!$  call b%fix(info)
-  if (info /= psb_success_) goto 9999
+!!$  if (info /= psb_success_) goto 9999
-
+!!$
-  call psb_erractionrestore(err_act)
+!!$  call psb_erractionrestore(err_act)
-  return
+!!$  return
-
+!!$
-9999 call psb_error_handler(err_act)
+!!$9999 call psb_error_handler(err_act)
-
+!!$
-  return
+!!$  return
-
+!!$
-end subroutine psb_c_csc_csgetblk
+!!$end subroutine psb_c_csc_csgetblk
 
 subroutine psb_c_csc_reinit(a,clear)
   use psb_error_mod

base/serial/impl/psb_c_csr_impl.f90

@@ -1999,7 +1999,7 @@ end subroutine psb_c_csr_csgetptn
 
 
 subroutine psb_c_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -2017,9 +2017,9 @@ subroutine psb_c_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
-  logical :: append_, rscale_, cscale_ 
+  logical :: append_, rscale_, cscale_, chksz_
   integer(psb_ipk_) :: nzin_, jmin_, jmax_, err_act, i
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -2063,13 +2063,18 @@ subroutine psb_c_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   else
     cscale_ = .false.
   endif
+  if (present(chksz)) then 
+    chksz_ = chksz
+  else
+    chksz_ = .true.
+  endif
   if ((rscale_.or.cscale_).and.(present(iren))) then 
     info = psb_err_many_optional_arg_
     call psb_errpush(info,name,a_err='iren (rscale.or.cscale)')
     goto 9999
   end if
 
-  call csr_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,info,&
+  call csr_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,chksz_,info,&
        & iren)
   
   if (rscale_) then 
@@ -2094,7 +2099,7 @@ subroutine psb_c_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
 
 contains
 
-  subroutine csr_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,info,&
+  subroutine csr_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,chksz,info,&
        & iren)
 
     use psb_const_mod
@@ -2109,7 +2114,7 @@ contains
     integer(psb_ipk_), allocatable, intent(inout)  :: ia(:), ja(:)
     complex(psb_spk_), allocatable,  intent(inout)    :: val(:)
     integer(psb_ipk_), intent(in)                  :: nzin
-    logical, intent(in)                  :: append
+    logical, intent(in)                  :: append, chksz
     integer(psb_ipk_) :: info
     integer(psb_ipk_), optional                    :: iren(:)
     integer(psb_ipk_) :: nzin_, nza, idx,i,j,k, nzt, irw, lrw, icl,lcl, nrd, ncd
@@ -2142,11 +2147,13 @@ contains
     nzt = (a%irp(lrw+1)-a%irp(irw))
     nz = 0 
 
-    call psb_ensure_size(nzin_+nzt,ia,info)
+    if (chksz) then 
-    if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+      call psb_ensure_size(nzin_+nzt,ia,info)
-    if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-
+      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-    if (info /= psb_success_) return
+      
+      if (info /= psb_success_) return
+    end if
     
     if (present(iren)) then 
       do i=irw, lrw
@@ -2178,50 +2185,155 @@ contains
 
 end subroutine psb_c_csr_csgetrow
 
-subroutine psb_c_csr_csgetblk(imin,imax,a,b,info,&
+
-     & jmin,jmax,iren,append,rscale,cscale)
+!
+! CSR implementation of tril/triu
+!
+subroutine psb_c_csr_tril(a,l,info,&
+     & diag,imin,imax,jmin,jmax,rscale,cscale,u)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
-  use psb_c_csr_mat_mod, psb_protect_name => psb_c_csr_csgetblk
+  use psb_c_csr_mat_mod, psb_protect_name => psb_c_csr_tril
   implicit none
 
   class(psb_c_csr_sparse_mat), intent(in) :: a
-  class(psb_c_coo_sparse_mat), intent(inout) :: b
+  class(psb_c_coo_sparse_mat), intent(out) :: l
-  integer(psb_ipk_), intent(in)                  :: imin,imax
+  integer(psb_ipk_),intent(out)            :: info
-  integer(psb_ipk_),intent(out)                  :: info
+  integer(psb_ipk_), intent(in), optional  :: diag,imin,imax,jmin,jmax
-  logical, intent(in), optional        :: append
+  logical, intent(in), optional            :: rscale,cscale
-  integer(psb_ipk_), intent(in), optional        :: iren(:)
+  class(psb_c_coo_sparse_mat), optional, intent(out) :: u
-  integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+  
-  logical, intent(in), optional        :: rscale,cscale
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k
-  integer(psb_ipk_) :: err_act, nzin, nzout
+  integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
+  integer(psb_ipk_), allocatable :: ia(:), ja(:)
+  complex(psb_spk_), allocatable    :: val(:)
   integer(psb_ipk_) :: ierr(5)
-  character(len=20)  :: name='csget'
+  character(len=20)  :: name='tril'
-  logical :: append_
+  logical :: rscale_, cscale_
   logical, parameter :: debug=.false.
 
   call psb_erractionsave(err_act)
   info = psb_success_
 
-  if (present(append)) then 
+  if (present(diag)) then 
-    append_ = append
+    diag_ = diag
   else
-    append_ = .false.
+    diag_ = 0
-  endif
+  end if
-  if (append_) then 
+  if (present(imin)) then 
-    nzin = a%get_nzeros()
+    imin_ = imin
+  else
+    imin_ = 1
+  end if
+  if (present(imax)) then 
+    imax_ = imax
   else
-    nzin = 0
+    imax_ = a%get_nrows()
+  end if
+  if (present(jmin)) then 
+    jmin_ = jmin
+  else
+    jmin_ = 1
+  end if
+  if (present(jmax)) then 
+    jmax_ = jmax
+  else
+    jmax_ = a%get_ncols()
+  end if
+  if (present(rscale)) then 
+    rscale_ = rscale
+  else
+    rscale_ = .true.
+  end if
+  if (present(cscale)) then 
+    cscale_ = cscale
+  else
+    cscale_ = .true.
+  end if
+
+  if (rscale_) then 
+    mb = imax_ - imin_ +1
+  else 
+    mb = imax_ 
+  endif
+  if (cscale_) then 
+    nb = jmax_ - jmin_ +1
+  else 
+    nb = jmax_  
   endif
 
-  call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
-       & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
-       & nzin=nzin, rscale=rscale, cscale=cscale)
 
-  if (info /= psb_success_) goto 9999
+  nz = a%get_nzeros()
+  call l%allocate(mb,nb,nz)
+  
+  if (present(u)) then
+    nzlin = l%get_nzeros() ! At this point it should be 0
+    call u%allocate(mb,nb,nz)
+    nzuin = u%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_
+        do k=irp(i),irp(i+1)-1
+          j = ja(k)
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)<=diag_) then
+              nzlin = nzlin + 1
+              l%ia(nzlin)  = i
+              l%ja(nzlin)  = ja(k)
+              l%val(nzlin) = val(k)
+            else
+              nzuin = nzuin + 1
+              u%ia(nzuin)  = i
+              u%ja(nzuin)  = ja(k)
+              u%val(nzuin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    
+    call l%set_nzeros(nzlin)
+    call u%set_nzeros(nzuin)
+    call u%fix(info)
+    nzout = u%get_nzeros()
+    if (rscale_) &
+         & u%ia(1:nzout) = u%ia(1:nzout) - imin_ + 1
+    if (cscale_) &
+         & u%ja(1:nzout) = u%ja(1:nzout) - jmin_ + 1  
+    if ((diag_ >=-1).and.(imin_ == jmin_)) then 
+      call u%set_triangle(.true.)
+      call u%set_lower(.false.)
+    end if
+  else
+    nzin = l%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_        
+        do k=irp(i),irp(i+1)-1
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)<=diag_) then
+              nzin = nzin + 1
+              l%ia(nzin)  = i
+              l%ja(nzin)  = ja(k)
+              l%val(nzin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    call l%set_nzeros(nzin)
+  end if
+  call l%fix(info)
+  nzout = l%get_nzeros()
+  if (rscale_) &
+       & l%ia(1:nzout) = l%ia(1:nzout) - imin_ + 1
+  if (cscale_) &
+       & l%ja(1:nzout) = l%ja(1:nzout) - jmin_ + 1
+  
+  if ((diag_ <= 0).and.(imin_ == jmin_)) then 
+    call l%set_triangle(.true.)
+    call l%set_lower(.true.)
+  end if
 
-  call b%set_nzeros(nzin+nzout)
-  call b%fix(info)
   if (info /= psb_success_) goto 9999
 
   call psb_erractionrestore(err_act)
@@ -2231,8 +2343,162 @@ subroutine psb_c_csr_csgetblk(imin,imax,a,b,info,&
 
   return
 
-end subroutine psb_c_csr_csgetblk
+end subroutine psb_c_csr_tril
+
+subroutine psb_c_csr_triu(a,u,info,&
+     & diag,imin,imax,jmin,jmax,rscale,cscale,l)
+  ! Output is always in  COO format 
+  use psb_error_mod
+  use psb_const_mod
+  use psb_c_csr_mat_mod, psb_protect_name => psb_c_csr_triu
+  implicit none
+
+  class(psb_c_csr_sparse_mat), intent(in) :: a
+  class(psb_c_coo_sparse_mat), intent(out) :: u
+  integer(psb_ipk_),intent(out)            :: info
+  integer(psb_ipk_), intent(in), optional  :: diag,imin,imax,jmin,jmax
+  logical, intent(in), optional            :: rscale,cscale
+  class(psb_c_coo_sparse_mat), optional, intent(out) :: l
+  
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k
+  integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
+  integer(psb_ipk_), allocatable :: ia(:), ja(:)
+  complex(psb_spk_), allocatable    :: val(:)
+  integer(psb_ipk_) :: ierr(5)
+  character(len=20)  :: name='triu'
+  logical :: rscale_, cscale_
+  logical, parameter :: debug=.false.
+
+  call psb_erractionsave(err_act)
+  info = psb_success_
+
+  if (present(diag)) then 
+    diag_ = diag
+  else
+    diag_ = 0
+  end if
+  if (present(imin)) then 
+    imin_ = imin
+  else
+    imin_ = 1
+  end if
+  if (present(imax)) then 
+    imax_ = imax
+  else
+    imax_ = a%get_nrows()
+  end if
+  if (present(jmin)) then 
+    jmin_ = jmin
+  else
+    jmin_ = 1
+  end if
+  if (present(jmax)) then 
+    jmax_ = jmax
+  else
+    jmax_ = a%get_ncols()
+  end if
+  if (present(rscale)) then 
+    rscale_ = rscale
+  else
+    rscale_ = .true.
+  end if
+  if (present(cscale)) then 
+    cscale_ = cscale
+  else
+    cscale_ = .true.
+  end if
+
+  if (rscale_) then 
+    mb = imax_ - imin_ +1
+  else 
+    mb = imax_ 
+  endif
+  if (cscale_) then 
+    nb = jmax_ - jmin_ +1
+  else 
+    nb = jmax_  
+  endif
+
+
+  nz = a%get_nzeros()
+  call u%allocate(mb,nb,nz)
+  
+  if (present(l)) then
+    nzuin = u%get_nzeros() ! At this point it should be 0
+    call l%allocate(mb,nb,nz)
+    nzlin = l%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_
+        do k=irp(i),irp(i+1)-1
+          j = ja(k)
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)<diag_) then
+              nzlin = nzlin + 1
+              l%ia(nzlin)  = i
+              l%ja(nzlin)  = ja(k)
+              l%val(nzlin) = val(k)
+            else
+              nzuin = nzuin + 1
+              u%ia(nzuin)  = i
+              u%ja(nzuin)  = ja(k)
+              u%val(nzuin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    call u%set_nzeros(nzuin)
+    call l%set_nzeros(nzlin)
+    call l%fix(info)
+    nzout = l%get_nzeros()
+    if (rscale_) &
+         & l%ia(1:nzout) = l%ia(1:nzout) - imin_ + 1
+    if (cscale_) &
+         & l%ja(1:nzout) = l%ja(1:nzout) - jmin_ + 1  
+    if ((diag_ <=1).and.(imin_ == jmin_)) then 
+      call l%set_triangle(.true.)
+      call l%set_lower(.true.)
+    end if
+  else
+    nzin = u%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_        
+        do k=irp(i),irp(i+1)-1
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)>=diag_) then
+              nzin = nzin + 1
+              u%ia(nzin)  = i
+              u%ja(nzin)  = ja(k)
+              u%val(nzin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    call u%set_nzeros(nzin)
+  end if
+  call u%fix(info)
+  nzout = u%get_nzeros()
+  if (rscale_) &
+       & u%ia(1:nzout) = u%ia(1:nzout) - imin_ + 1
+  if (cscale_) &
+       & u%ja(1:nzout) = u%ja(1:nzout) - jmin_ + 1
+  
+  if ((diag_ >= 0).and.(imin_ == jmin_)) then 
+    call u%set_triangle(.true.)
+    call u%set_upper(.true.)
+  end if
+
+  if (info /= psb_success_) goto 9999
+
+  call psb_erractionrestore(err_act)
+  return
+
+9999 call psb_error_handler(err_act)
+
+  return
 
+end subroutine psb_c_csr_triu
 
 
 subroutine psb_c_csr_csput_a(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl) 

base/serial/impl/psb_c_mat_impl.F90

@@ -793,7 +793,7 @@ end subroutine psb_c_csgetptn
 
 
 subroutine psb_c_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -810,7 +810,7 @@ subroutine psb_c_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
   integer(psb_ipk_) :: err_act
   character(len=20)  :: name='csget'
@@ -826,7 +826,7 @@ subroutine psb_c_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
 
 
   call a%a%csget(imin,imax,nz,ia,ja,val,info,&
-       & jmin,jmax,iren,append,nzin,rscale,cscale)
+       & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   if (info /= psb_success_) goto 9999 
 
   call psb_erractionrestore(err_act)

base/serial/impl/psb_d_base_mat_impl.F90

@@ -406,7 +406,7 @@ subroutine psb_d_base_csput_v(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl)
 end subroutine psb_d_base_csput_v
 
 subroutine psb_d_base_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -422,7 +422,7 @@ subroutine psb_d_base_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
   integer(psb_ipk_) :: err_act
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -447,7 +447,7 @@ end subroutine psb_d_base_csgetrow
 ! If performance is critical it can be overridden.
 !
 subroutine psb_d_base_csgetblk(imin,imax,a,b,info,&
-     & jmin,jmax,iren,append,rscale,cscale)
+     & jmin,jmax,iren,append,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -461,7 +461,7 @@ subroutine psb_d_base_csgetblk(imin,imax,a,b,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
   integer(psb_ipk_) :: err_act, nzin, nzout
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -522,7 +522,7 @@ subroutine psb_d_base_csgetblk(imin,imax,a,b,info,&
 
   call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
        & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
-       & nzin=nzin, rscale=rscale, cscale=cscale)
+       & nzin=nzin, rscale=rscale, cscale=cscale, chksz=chksz)
 
   if (info /= psb_success_) goto 9999
 
@@ -645,7 +645,7 @@ subroutine psb_d_base_tril(a,l,info,&
   logical, intent(in), optional            :: rscale,cscale
   class(psb_d_coo_sparse_mat), optional, intent(out) :: u
   
-  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k 
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k, ibk 
   integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
   integer(psb_ipk_), allocatable :: ia(:), ja(:)
   real(psb_dpk_), allocatable    :: val(:)
@@ -653,6 +653,7 @@ subroutine psb_d_base_tril(a,l,info,&
   character(len=20)  :: name='tril'
   logical :: rscale_, cscale_
   logical, parameter :: debug=.false.
+  integer(psb_ipk_), parameter :: nbk=8
 
   call psb_erractionsave(err_act)
   info = psb_success_
@@ -715,12 +716,12 @@ subroutine psb_d_base_tril(a,l,info,&
     call psb_realloc(max(mb,nb),ia,info)
     call psb_realloc(max(mb,nb),ja,info)
     call psb_realloc(max(mb,nb),val,info)
-    do i=imin_,imax_
+    do i=imin_,imax_, nbk
-      call a%csget(i,i,nzout,ia,ja,val,info,&
+      ibk = min(nbk,imax_-i+1)
+      call a%csget(i,i+ibk-1,nzout,ia,ja,val,info,&
            & jmin=jmin_, jmax=jmax_)
       do k=1, nzout
-        j = ja(k) 
+        if ((ja(k)-ia(k))<=diag_) then
-        if (j-i<=diag_) then
           nzlin = nzlin + 1
           l%ia(nzlin)  = ia(k)
           l%ja(nzlin)  = ja(k)
@@ -796,7 +797,7 @@ subroutine psb_d_base_triu(a,u,info,&
   logical, intent(in), optional            :: rscale,cscale
   class(psb_d_coo_sparse_mat), optional, intent(out) :: l
   
-  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k 
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k, ibk
   integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
   integer(psb_ipk_), allocatable :: ia(:), ja(:)
   real(psb_dpk_), allocatable    :: val(:)
@@ -804,6 +805,7 @@ subroutine psb_d_base_triu(a,u,info,&
   character(len=20)  :: name='triu'
   logical :: rscale_, cscale_
   logical, parameter :: debug=.false.
+  integer(psb_ipk_), parameter :: nbk=8
 
   call psb_erractionsave(err_act)
   info = psb_success_
@@ -866,12 +868,12 @@ subroutine psb_d_base_triu(a,u,info,&
     call psb_realloc(max(mb,nb),ia,info)
     call psb_realloc(max(mb,nb),ja,info)
     call psb_realloc(max(mb,nb),val,info)
-    do i=imin_,imax_
+    do i=imin_,imax_, nbk
-      call a%csget(i,i,nzout,ia,ja,val,info,&
+      ibk = min(nbk,imax_-i+1)
+      call a%csget(i,i+ibk-1,nzout,ia,ja,val,info,&
            & jmin=jmin_, jmax=jmax_)
       do k=1, nzout
-        j = ja(k) 
+        if ((ja(k)-ia(k))<diag_) then
-        if (j-i<diag_) then
           nzlin = nzlin + 1
           l%ia(nzlin)  = ia(k)
           l%ja(nzlin)  = ja(k)

base/serial/impl/psb_d_coo_impl.f90

@@ -2259,7 +2259,7 @@ end subroutine psb_d_coo_csgetptn
 ! The output is guaranteed to be sorted
 ! 
 subroutine psb_d_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -2276,9 +2276,9 @@ subroutine psb_d_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
-  logical :: append_, rscale_, cscale_ 
+  logical :: append_, rscale_, cscale_, chksz_
   integer(psb_ipk_) :: nzin_, jmin_, jmax_, err_act, i
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -2321,13 +2321,18 @@ subroutine psb_d_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   else
     cscale_ = .false.
   endif
+  if (present(chksz)) then 
+    chksz_ = chksz
+  else
+    chksz_ = .true.
+  endif
   if ((rscale_.or.cscale_).and.(present(iren))) then 
     info = psb_err_many_optional_arg_
     call psb_errpush(info,name,a_err='iren (rscale.or.cscale)')
     goto 9999
   end if
 
-  call coo_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,info,&
+  call coo_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,chksz_,info,&
        & iren)
 
   if (rscale_) then 
@@ -2352,7 +2357,7 @@ subroutine psb_d_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
 
 contains
 
-  subroutine coo_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,info,&
+  subroutine coo_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,chksz,info,&
        & iren)
 
     use psb_const_mod
@@ -2368,7 +2373,7 @@ contains
     integer(psb_ipk_), allocatable, intent(inout)  :: ia(:), ja(:)
     real(psb_dpk_), allocatable,  intent(inout)    :: val(:)
     integer(psb_ipk_), intent(in)                  :: nzin
-    logical, intent(in)                            :: append
+    logical, intent(in)                            :: append,chksz
     integer(psb_ipk_) :: info
     integer(psb_ipk_), optional                    :: iren(:)
     integer(psb_ipk_) :: nzin_, nza, idx,ip,jp,i,k, nzt, irw, lrw, nra, nca, nrd
@@ -2452,11 +2457,13 @@ contains
         nzt = jp - ip +1 
         nz = 0 
 
-        call psb_ensure_size(nzin_+nzt,ia,info)
+        if (chksz) then 
-        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+          call psb_ensure_size(nzin_+nzt,ia,info)
-        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+          if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-        if (info /= psb_success_) return
+          if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-
+          if (info /= psb_success_) return
+        end if
+        
         if (present(iren)) then 
           do i=ip,jp
             if ((jmin <= a%ja(i)).and.(a%ja(i)<=jmax)) then 
@@ -2488,11 +2495,13 @@ contains
 
       nrd = max(a%get_nrows(),1)
       nzt = ((nza+nrd-1)/nrd)*(lrw-irw+1)
-      call psb_ensure_size(nzin_+nzt,ia,info)
+      if (chksz) then 
-      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+        call psb_ensure_size(nzin_+nzt,ia,info)
-      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-      if (info /= psb_success_) return
+        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-
+        if (info /= psb_success_) return
+      end if
+      
       if (present(iren)) then 
         k = 0 
         do i=1, a%get_nzeros()
@@ -2501,10 +2510,12 @@ contains
             k = k + 1 
             if (k > nzt) then
               nzt = k + nzt
-              call psb_ensure_size(nzin_+nzt,ia,info)
+              if (chksz) then 
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+                call psb_ensure_size(nzin_+nzt,ia,info)
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-              if (info /= psb_success_) return
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+                if (info /= psb_success_) return
+              end if              
             end if
             val(nzin_+k) = a%val(i)
             ia(nzin_+k)  = iren(a%ia(i))
@@ -2519,11 +2530,12 @@ contains
             k = k + 1 
             if (k > nzt) then
               nzt = k + nzt
-              call psb_ensure_size(nzin_+nzt,ia,info)
+              if (chksz) then 
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+                call psb_ensure_size(nzin_+nzt,ia,info)
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-              if (info /= psb_success_) return
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-
+                if (info /= psb_success_) return
+              end if
             end if
             val(nzin_+k) = a%val(i)
             ia(nzin_+k)  = (a%ia(i))

base/serial/impl/psb_d_csc_impl.f90

@@ -1686,7 +1686,7 @@ end subroutine psb_d_csc_csgetptn
 
 
 subroutine psb_d_csc_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -1704,7 +1704,7 @@ subroutine psb_d_csc_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
   logical :: append_, rscale_, cscale_ 
   integer(psb_ipk_) :: nzin_, jmin_, jmax_, err_act, i
@@ -2557,60 +2557,60 @@ end subroutine psb_d_csc_reallocate_nz
 
 
 
-subroutine psb_d_csc_csgetblk(imin,imax,a,b,info,&
+!!$subroutine psb_d_csc_csgetblk(imin,imax,a,b,info,&
-     & jmin,jmax,iren,append,rscale,cscale)
+!!$     & jmin,jmax,iren,append,rscale,cscale)
-  ! Output is always in  COO format 
+!!$  ! Output is always in  COO format 
-  use psb_error_mod
+!!$  use psb_error_mod
-  use psb_const_mod
+!!$  use psb_const_mod
-  use psb_d_csc_mat_mod, psb_protect_name => psb_d_csc_csgetblk
+!!$  use psb_d_csc_mat_mod, psb_protect_name => psb_d_csc_csgetblk
-  implicit none
+!!$  implicit none
-
+!!$
-  class(psb_d_csc_sparse_mat), intent(in) :: a
+!!$  class(psb_d_csc_sparse_mat), intent(in) :: a
-  class(psb_d_coo_sparse_mat), intent(inout) :: b
+!!$  class(psb_d_coo_sparse_mat), intent(inout) :: b
-  integer(psb_ipk_), intent(in)                  :: imin,imax
+!!$  integer(psb_ipk_), intent(in)                  :: imin,imax
-  integer(psb_ipk_),intent(out)                  :: info
+!!$  integer(psb_ipk_),intent(out)                  :: info
-  logical, intent(in), optional        :: append
+!!$  logical, intent(in), optional        :: append
-  integer(psb_ipk_), intent(in), optional        :: iren(:)
+!!$  integer(psb_ipk_), intent(in), optional        :: iren(:)
-  integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+!!$  integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-  logical, intent(in), optional        :: rscale,cscale
+!!$  logical, intent(in), optional        :: rscale,cscale
-  integer(psb_ipk_) :: err_act, nzin, nzout
+!!$  integer(psb_ipk_) :: err_act, nzin, nzout
-  integer(psb_ipk_) :: ierr(5)
+!!$  integer(psb_ipk_) :: ierr(5)
-  character(len=20)  :: name='csget'
+!!$  character(len=20)  :: name='csget'
-  logical :: append_
+!!$  logical :: append_
-  logical, parameter :: debug=.false.
+!!$  logical, parameter :: debug=.false.
-
+!!$
-  call psb_erractionsave(err_act)
+!!$  call psb_erractionsave(err_act)
-  info = psb_success_
+!!$  info = psb_success_
-
+!!$
-  if (present(append)) then 
+!!$  if (present(append)) then 
-    append_ = append
+!!$    append_ = append
-  else
+!!$  else
-    append_ = .false.
+!!$    append_ = .false.
-  endif
+!!$  endif
-  if (append_) then 
+!!$  if (append_) then 
-    nzin = a%get_nzeros()
+!!$    nzin = a%get_nzeros()
-  else
+!!$  else
-    nzin = 0
+!!$    nzin = 0
-  endif
+!!$  endif
-
+!!$
-  call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
+!!$  call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
-       & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
+!!$       & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
-       & nzin=nzin, rscale=rscale, cscale=cscale)
+!!$       & nzin=nzin, rscale=rscale, cscale=cscale)
-
+!!$
-  if (info /= psb_success_) goto 9999
+!!$  if (info /= psb_success_) goto 9999
-
+!!$
-  call b%set_nzeros(nzin+nzout)
+!!$  call b%set_nzeros(nzin+nzout)
-  call b%fix(info)
+!!$  call b%fix(info)
-  if (info /= psb_success_) goto 9999
+!!$  if (info /= psb_success_) goto 9999
-
+!!$
-  call psb_erractionrestore(err_act)
+!!$  call psb_erractionrestore(err_act)
-  return
+!!$  return
-
+!!$
-9999 call psb_error_handler(err_act)
+!!$9999 call psb_error_handler(err_act)
-
+!!$
-  return
+!!$  return
-
+!!$
-end subroutine psb_d_csc_csgetblk
+!!$end subroutine psb_d_csc_csgetblk
 
 subroutine psb_d_csc_reinit(a,clear)
   use psb_error_mod

base/serial/impl/psb_d_csr_impl.f90

@@ -1999,7 +1999,7 @@ end subroutine psb_d_csr_csgetptn
 
 
 subroutine psb_d_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -2017,9 +2017,9 @@ subroutine psb_d_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
-  logical :: append_, rscale_, cscale_ 
+  logical :: append_, rscale_, cscale_, chksz_
   integer(psb_ipk_) :: nzin_, jmin_, jmax_, err_act, i
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -2063,13 +2063,18 @@ subroutine psb_d_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   else
     cscale_ = .false.
   endif
+  if (present(chksz)) then 
+    chksz_ = chksz
+  else
+    chksz_ = .true.
+  endif
   if ((rscale_.or.cscale_).and.(present(iren))) then 
     info = psb_err_many_optional_arg_
     call psb_errpush(info,name,a_err='iren (rscale.or.cscale)')
     goto 9999
   end if
 
-  call csr_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,info,&
+  call csr_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,chksz_,info,&
        & iren)
   
   if (rscale_) then 
@@ -2094,7 +2099,7 @@ subroutine psb_d_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
 
 contains
 
-  subroutine csr_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,info,&
+  subroutine csr_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,chksz,info,&
        & iren)
 
     use psb_const_mod
@@ -2109,7 +2114,7 @@ contains
     integer(psb_ipk_), allocatable, intent(inout)  :: ia(:), ja(:)
     real(psb_dpk_), allocatable,  intent(inout)    :: val(:)
     integer(psb_ipk_), intent(in)                  :: nzin
-    logical, intent(in)                  :: append
+    logical, intent(in)                  :: append, chksz
     integer(psb_ipk_) :: info
     integer(psb_ipk_), optional                    :: iren(:)
     integer(psb_ipk_) :: nzin_, nza, idx,i,j,k, nzt, irw, lrw, icl,lcl, nrd, ncd
@@ -2142,11 +2147,13 @@ contains
     nzt = (a%irp(lrw+1)-a%irp(irw))
     nz = 0 
 
-    call psb_ensure_size(nzin_+nzt,ia,info)
+    if (chksz) then 
-    if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+      call psb_ensure_size(nzin_+nzt,ia,info)
-    if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-
+      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-    if (info /= psb_success_) return
+      
+      if (info /= psb_success_) return
+    end if
     
     if (present(iren)) then 
       do i=irw, lrw
@@ -2178,50 +2185,155 @@ contains
 
 end subroutine psb_d_csr_csgetrow
 
-subroutine psb_d_csr_csgetblk(imin,imax,a,b,info,&
+
-     & jmin,jmax,iren,append,rscale,cscale)
+!
+! CSR implementation of tril/triu
+!
+subroutine psb_d_csr_tril(a,l,info,&
+     & diag,imin,imax,jmin,jmax,rscale,cscale,u)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
-  use psb_d_csr_mat_mod, psb_protect_name => psb_d_csr_csgetblk
+  use psb_d_csr_mat_mod, psb_protect_name => psb_d_csr_tril
   implicit none
 
   class(psb_d_csr_sparse_mat), intent(in) :: a
-  class(psb_d_coo_sparse_mat), intent(inout) :: b
+  class(psb_d_coo_sparse_mat), intent(out) :: l
-  integer(psb_ipk_), intent(in)                  :: imin,imax
+  integer(psb_ipk_),intent(out)            :: info
-  integer(psb_ipk_),intent(out)                  :: info
+  integer(psb_ipk_), intent(in), optional  :: diag,imin,imax,jmin,jmax
-  logical, intent(in), optional        :: append
+  logical, intent(in), optional            :: rscale,cscale
-  integer(psb_ipk_), intent(in), optional        :: iren(:)
+  class(psb_d_coo_sparse_mat), optional, intent(out) :: u
-  integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+  
-  logical, intent(in), optional        :: rscale,cscale
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k
-  integer(psb_ipk_) :: err_act, nzin, nzout
+  integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
+  integer(psb_ipk_), allocatable :: ia(:), ja(:)
+  real(psb_dpk_), allocatable    :: val(:)
   integer(psb_ipk_) :: ierr(5)
-  character(len=20)  :: name='csget'
+  character(len=20)  :: name='tril'
-  logical :: append_
+  logical :: rscale_, cscale_
   logical, parameter :: debug=.false.
 
   call psb_erractionsave(err_act)
   info = psb_success_
 
-  if (present(append)) then 
+  if (present(diag)) then 
-    append_ = append
+    diag_ = diag
   else
-    append_ = .false.
+    diag_ = 0
-  endif
+  end if
-  if (append_) then 
+  if (present(imin)) then 
-    nzin = a%get_nzeros()
+    imin_ = imin
+  else
+    imin_ = 1
+  end if
+  if (present(imax)) then 
+    imax_ = imax
   else
-    nzin = 0
+    imax_ = a%get_nrows()
+  end if
+  if (present(jmin)) then 
+    jmin_ = jmin
+  else
+    jmin_ = 1
+  end if
+  if (present(jmax)) then 
+    jmax_ = jmax
+  else
+    jmax_ = a%get_ncols()
+  end if
+  if (present(rscale)) then 
+    rscale_ = rscale
+  else
+    rscale_ = .true.
+  end if
+  if (present(cscale)) then 
+    cscale_ = cscale
+  else
+    cscale_ = .true.
+  end if
+
+  if (rscale_) then 
+    mb = imax_ - imin_ +1
+  else 
+    mb = imax_ 
+  endif
+  if (cscale_) then 
+    nb = jmax_ - jmin_ +1
+  else 
+    nb = jmax_  
   endif
 
-  call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
-       & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
-       & nzin=nzin, rscale=rscale, cscale=cscale)
 
-  if (info /= psb_success_) goto 9999
+  nz = a%get_nzeros()
+  call l%allocate(mb,nb,nz)
+  
+  if (present(u)) then
+    nzlin = l%get_nzeros() ! At this point it should be 0
+    call u%allocate(mb,nb,nz)
+    nzuin = u%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_
+        do k=irp(i),irp(i+1)-1
+          j = ja(k)
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)<=diag_) then
+              nzlin = nzlin + 1
+              l%ia(nzlin)  = i
+              l%ja(nzlin)  = ja(k)
+              l%val(nzlin) = val(k)
+            else
+              nzuin = nzuin + 1
+              u%ia(nzuin)  = i
+              u%ja(nzuin)  = ja(k)
+              u%val(nzuin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    
+    call l%set_nzeros(nzlin)
+    call u%set_nzeros(nzuin)
+    call u%fix(info)
+    nzout = u%get_nzeros()
+    if (rscale_) &
+         & u%ia(1:nzout) = u%ia(1:nzout) - imin_ + 1
+    if (cscale_) &
+         & u%ja(1:nzout) = u%ja(1:nzout) - jmin_ + 1  
+    if ((diag_ >=-1).and.(imin_ == jmin_)) then 
+      call u%set_triangle(.true.)
+      call u%set_lower(.false.)
+    end if
+  else
+    nzin = l%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_        
+        do k=irp(i),irp(i+1)-1
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)<=diag_) then
+              nzin = nzin + 1
+              l%ia(nzin)  = i
+              l%ja(nzin)  = ja(k)
+              l%val(nzin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    call l%set_nzeros(nzin)
+  end if
+  call l%fix(info)
+  nzout = l%get_nzeros()
+  if (rscale_) &
+       & l%ia(1:nzout) = l%ia(1:nzout) - imin_ + 1
+  if (cscale_) &
+       & l%ja(1:nzout) = l%ja(1:nzout) - jmin_ + 1
+  
+  if ((diag_ <= 0).and.(imin_ == jmin_)) then 
+    call l%set_triangle(.true.)
+    call l%set_lower(.true.)
+  end if
 
-  call b%set_nzeros(nzin+nzout)
-  call b%fix(info)
   if (info /= psb_success_) goto 9999
 
   call psb_erractionrestore(err_act)
@@ -2231,8 +2343,162 @@ subroutine psb_d_csr_csgetblk(imin,imax,a,b,info,&
 
   return
 
-end subroutine psb_d_csr_csgetblk
+end subroutine psb_d_csr_tril
+
+subroutine psb_d_csr_triu(a,u,info,&
+     & diag,imin,imax,jmin,jmax,rscale,cscale,l)
+  ! Output is always in  COO format 
+  use psb_error_mod
+  use psb_const_mod
+  use psb_d_csr_mat_mod, psb_protect_name => psb_d_csr_triu
+  implicit none
+
+  class(psb_d_csr_sparse_mat), intent(in) :: a
+  class(psb_d_coo_sparse_mat), intent(out) :: u
+  integer(psb_ipk_),intent(out)            :: info
+  integer(psb_ipk_), intent(in), optional  :: diag,imin,imax,jmin,jmax
+  logical, intent(in), optional            :: rscale,cscale
+  class(psb_d_coo_sparse_mat), optional, intent(out) :: l
+  
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k
+  integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
+  integer(psb_ipk_), allocatable :: ia(:), ja(:)
+  real(psb_dpk_), allocatable    :: val(:)
+  integer(psb_ipk_) :: ierr(5)
+  character(len=20)  :: name='triu'
+  logical :: rscale_, cscale_
+  logical, parameter :: debug=.false.
+
+  call psb_erractionsave(err_act)
+  info = psb_success_
+
+  if (present(diag)) then 
+    diag_ = diag
+  else
+    diag_ = 0
+  end if
+  if (present(imin)) then 
+    imin_ = imin
+  else
+    imin_ = 1
+  end if
+  if (present(imax)) then 
+    imax_ = imax
+  else
+    imax_ = a%get_nrows()
+  end if
+  if (present(jmin)) then 
+    jmin_ = jmin
+  else
+    jmin_ = 1
+  end if
+  if (present(jmax)) then 
+    jmax_ = jmax
+  else
+    jmax_ = a%get_ncols()
+  end if
+  if (present(rscale)) then 
+    rscale_ = rscale
+  else
+    rscale_ = .true.
+  end if
+  if (present(cscale)) then 
+    cscale_ = cscale
+  else
+    cscale_ = .true.
+  end if
+
+  if (rscale_) then 
+    mb = imax_ - imin_ +1
+  else 
+    mb = imax_ 
+  endif
+  if (cscale_) then 
+    nb = jmax_ - jmin_ +1
+  else 
+    nb = jmax_  
+  endif
+
+
+  nz = a%get_nzeros()
+  call u%allocate(mb,nb,nz)
+  
+  if (present(l)) then
+    nzuin = u%get_nzeros() ! At this point it should be 0
+    call l%allocate(mb,nb,nz)
+    nzlin = l%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_
+        do k=irp(i),irp(i+1)-1
+          j = ja(k)
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)<diag_) then
+              nzlin = nzlin + 1
+              l%ia(nzlin)  = i
+              l%ja(nzlin)  = ja(k)
+              l%val(nzlin) = val(k)
+            else
+              nzuin = nzuin + 1
+              u%ia(nzuin)  = i
+              u%ja(nzuin)  = ja(k)
+              u%val(nzuin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    call u%set_nzeros(nzuin)
+    call l%set_nzeros(nzlin)
+    call l%fix(info)
+    nzout = l%get_nzeros()
+    if (rscale_) &
+         & l%ia(1:nzout) = l%ia(1:nzout) - imin_ + 1
+    if (cscale_) &
+         & l%ja(1:nzout) = l%ja(1:nzout) - jmin_ + 1  
+    if ((diag_ <=1).and.(imin_ == jmin_)) then 
+      call l%set_triangle(.true.)
+      call l%set_lower(.true.)
+    end if
+  else
+    nzin = u%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_        
+        do k=irp(i),irp(i+1)-1
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)>=diag_) then
+              nzin = nzin + 1
+              u%ia(nzin)  = i
+              u%ja(nzin)  = ja(k)
+              u%val(nzin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    call u%set_nzeros(nzin)
+  end if
+  call u%fix(info)
+  nzout = u%get_nzeros()
+  if (rscale_) &
+       & u%ia(1:nzout) = u%ia(1:nzout) - imin_ + 1
+  if (cscale_) &
+       & u%ja(1:nzout) = u%ja(1:nzout) - jmin_ + 1
+  
+  if ((diag_ >= 0).and.(imin_ == jmin_)) then 
+    call u%set_triangle(.true.)
+    call u%set_upper(.true.)
+  end if
+
+  if (info /= psb_success_) goto 9999
+
+  call psb_erractionrestore(err_act)
+  return
+
+9999 call psb_error_handler(err_act)
+
+  return
 
+end subroutine psb_d_csr_triu
 
 
 subroutine psb_d_csr_csput_a(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl) 

base/serial/impl/psb_d_mat_impl.F90

@@ -793,7 +793,7 @@ end subroutine psb_d_csgetptn
 
 
 subroutine psb_d_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -810,7 +810,7 @@ subroutine psb_d_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
   integer(psb_ipk_) :: err_act
   character(len=20)  :: name='csget'
@@ -826,7 +826,7 @@ subroutine psb_d_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
 
 
   call a%a%csget(imin,imax,nz,ia,ja,val,info,&
-       & jmin,jmax,iren,append,nzin,rscale,cscale)
+       & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   if (info /= psb_success_) goto 9999 
 
   call psb_erractionrestore(err_act)

base/serial/impl/psb_s_base_mat_impl.F90

@@ -406,7 +406,7 @@ subroutine psb_s_base_csput_v(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl)
 end subroutine psb_s_base_csput_v
 
 subroutine psb_s_base_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -422,7 +422,7 @@ subroutine psb_s_base_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
   integer(psb_ipk_) :: err_act
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -447,7 +447,7 @@ end subroutine psb_s_base_csgetrow
 ! If performance is critical it can be overridden.
 !
 subroutine psb_s_base_csgetblk(imin,imax,a,b,info,&
-     & jmin,jmax,iren,append,rscale,cscale)
+     & jmin,jmax,iren,append,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -461,7 +461,7 @@ subroutine psb_s_base_csgetblk(imin,imax,a,b,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
   integer(psb_ipk_) :: err_act, nzin, nzout
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -522,7 +522,7 @@ subroutine psb_s_base_csgetblk(imin,imax,a,b,info,&
 
   call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
        & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
-       & nzin=nzin, rscale=rscale, cscale=cscale)
+       & nzin=nzin, rscale=rscale, cscale=cscale, chksz=chksz)
 
   if (info /= psb_success_) goto 9999
 
@@ -645,7 +645,7 @@ subroutine psb_s_base_tril(a,l,info,&
   logical, intent(in), optional            :: rscale,cscale
   class(psb_s_coo_sparse_mat), optional, intent(out) :: u
   
-  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k 
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k, ibk 
   integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
   integer(psb_ipk_), allocatable :: ia(:), ja(:)
   real(psb_spk_), allocatable    :: val(:)
@@ -653,6 +653,7 @@ subroutine psb_s_base_tril(a,l,info,&
   character(len=20)  :: name='tril'
   logical :: rscale_, cscale_
   logical, parameter :: debug=.false.
+  integer(psb_ipk_), parameter :: nbk=8
 
   call psb_erractionsave(err_act)
   info = psb_success_
@@ -715,12 +716,12 @@ subroutine psb_s_base_tril(a,l,info,&
     call psb_realloc(max(mb,nb),ia,info)
     call psb_realloc(max(mb,nb),ja,info)
     call psb_realloc(max(mb,nb),val,info)
-    do i=imin_,imax_
+    do i=imin_,imax_, nbk
-      call a%csget(i,i,nzout,ia,ja,val,info,&
+      ibk = min(nbk,imax_-i+1)
+      call a%csget(i,i+ibk-1,nzout,ia,ja,val,info,&
            & jmin=jmin_, jmax=jmax_)
       do k=1, nzout
-        j = ja(k) 
+        if ((ja(k)-ia(k))<=diag_) then
-        if (j-i<=diag_) then
           nzlin = nzlin + 1
           l%ia(nzlin)  = ia(k)
           l%ja(nzlin)  = ja(k)
@@ -796,7 +797,7 @@ subroutine psb_s_base_triu(a,u,info,&
   logical, intent(in), optional            :: rscale,cscale
   class(psb_s_coo_sparse_mat), optional, intent(out) :: l
   
-  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k 
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k, ibk
   integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
   integer(psb_ipk_), allocatable :: ia(:), ja(:)
   real(psb_spk_), allocatable    :: val(:)
@@ -804,6 +805,7 @@ subroutine psb_s_base_triu(a,u,info,&
   character(len=20)  :: name='triu'
   logical :: rscale_, cscale_
   logical, parameter :: debug=.false.
+  integer(psb_ipk_), parameter :: nbk=8
 
   call psb_erractionsave(err_act)
   info = psb_success_
@@ -866,12 +868,12 @@ subroutine psb_s_base_triu(a,u,info,&
     call psb_realloc(max(mb,nb),ia,info)
     call psb_realloc(max(mb,nb),ja,info)
     call psb_realloc(max(mb,nb),val,info)
-    do i=imin_,imax_
+    do i=imin_,imax_, nbk
-      call a%csget(i,i,nzout,ia,ja,val,info,&
+      ibk = min(nbk,imax_-i+1)
+      call a%csget(i,i+ibk-1,nzout,ia,ja,val,info,&
            & jmin=jmin_, jmax=jmax_)
       do k=1, nzout
-        j = ja(k) 
+        if ((ja(k)-ia(k))<diag_) then
-        if (j-i<diag_) then
           nzlin = nzlin + 1
           l%ia(nzlin)  = ia(k)
           l%ja(nzlin)  = ja(k)

base/serial/impl/psb_s_coo_impl.f90

@@ -2259,7 +2259,7 @@ end subroutine psb_s_coo_csgetptn
 ! The output is guaranteed to be sorted
 ! 
 subroutine psb_s_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -2276,9 +2276,9 @@ subroutine psb_s_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
-  logical :: append_, rscale_, cscale_ 
+  logical :: append_, rscale_, cscale_, chksz_
   integer(psb_ipk_) :: nzin_, jmin_, jmax_, err_act, i
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -2321,13 +2321,18 @@ subroutine psb_s_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   else
     cscale_ = .false.
   endif
+  if (present(chksz)) then 
+    chksz_ = chksz
+  else
+    chksz_ = .true.
+  endif
   if ((rscale_.or.cscale_).and.(present(iren))) then 
     info = psb_err_many_optional_arg_
     call psb_errpush(info,name,a_err='iren (rscale.or.cscale)')
     goto 9999
   end if
 
-  call coo_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,info,&
+  call coo_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,chksz_,info,&
        & iren)
 
   if (rscale_) then 
@@ -2352,7 +2357,7 @@ subroutine psb_s_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
 
 contains
 
-  subroutine coo_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,info,&
+  subroutine coo_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,chksz,info,&
        & iren)
 
     use psb_const_mod
@@ -2368,7 +2373,7 @@ contains
     integer(psb_ipk_), allocatable, intent(inout)  :: ia(:), ja(:)
     real(psb_spk_), allocatable,  intent(inout)    :: val(:)
     integer(psb_ipk_), intent(in)                  :: nzin
-    logical, intent(in)                            :: append
+    logical, intent(in)                            :: append,chksz
     integer(psb_ipk_) :: info
     integer(psb_ipk_), optional                    :: iren(:)
     integer(psb_ipk_) :: nzin_, nza, idx,ip,jp,i,k, nzt, irw, lrw, nra, nca, nrd
@@ -2452,11 +2457,13 @@ contains
         nzt = jp - ip +1 
         nz = 0 
 
-        call psb_ensure_size(nzin_+nzt,ia,info)
+        if (chksz) then 
-        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+          call psb_ensure_size(nzin_+nzt,ia,info)
-        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+          if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-        if (info /= psb_success_) return
+          if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-
+          if (info /= psb_success_) return
+        end if
+        
         if (present(iren)) then 
           do i=ip,jp
             if ((jmin <= a%ja(i)).and.(a%ja(i)<=jmax)) then 
@@ -2488,11 +2495,13 @@ contains
 
       nrd = max(a%get_nrows(),1)
       nzt = ((nza+nrd-1)/nrd)*(lrw-irw+1)
-      call psb_ensure_size(nzin_+nzt,ia,info)
+      if (chksz) then 
-      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+        call psb_ensure_size(nzin_+nzt,ia,info)
-      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-      if (info /= psb_success_) return
+        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-
+        if (info /= psb_success_) return
+      end if
+      
       if (present(iren)) then 
         k = 0 
         do i=1, a%get_nzeros()
@@ -2501,10 +2510,12 @@ contains
             k = k + 1 
             if (k > nzt) then
               nzt = k + nzt
-              call psb_ensure_size(nzin_+nzt,ia,info)
+              if (chksz) then 
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+                call psb_ensure_size(nzin_+nzt,ia,info)
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-              if (info /= psb_success_) return
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+                if (info /= psb_success_) return
+              end if              
             end if
             val(nzin_+k) = a%val(i)
             ia(nzin_+k)  = iren(a%ia(i))
@@ -2519,11 +2530,12 @@ contains
             k = k + 1 
             if (k > nzt) then
               nzt = k + nzt
-              call psb_ensure_size(nzin_+nzt,ia,info)
+              if (chksz) then 
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+                call psb_ensure_size(nzin_+nzt,ia,info)
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-              if (info /= psb_success_) return
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-
+                if (info /= psb_success_) return
+              end if
             end if
             val(nzin_+k) = a%val(i)
             ia(nzin_+k)  = (a%ia(i))

base/serial/impl/psb_s_csc_impl.f90

@@ -1686,7 +1686,7 @@ end subroutine psb_s_csc_csgetptn
 
 
 subroutine psb_s_csc_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -1704,7 +1704,7 @@ subroutine psb_s_csc_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
   logical :: append_, rscale_, cscale_ 
   integer(psb_ipk_) :: nzin_, jmin_, jmax_, err_act, i
@@ -2557,60 +2557,60 @@ end subroutine psb_s_csc_reallocate_nz
 
 
 
-subroutine psb_s_csc_csgetblk(imin,imax,a,b,info,&
+!!$subroutine psb_s_csc_csgetblk(imin,imax,a,b,info,&
-     & jmin,jmax,iren,append,rscale,cscale)
+!!$     & jmin,jmax,iren,append,rscale,cscale)
-  ! Output is always in  COO format 
+!!$  ! Output is always in  COO format 
-  use psb_error_mod
+!!$  use psb_error_mod
-  use psb_const_mod
+!!$  use psb_const_mod
-  use psb_s_csc_mat_mod, psb_protect_name => psb_s_csc_csgetblk
+!!$  use psb_s_csc_mat_mod, psb_protect_name => psb_s_csc_csgetblk
-  implicit none
+!!$  implicit none
-
+!!$
-  class(psb_s_csc_sparse_mat), intent(in) :: a
+!!$  class(psb_s_csc_sparse_mat), intent(in) :: a
-  class(psb_s_coo_sparse_mat), intent(inout) :: b
+!!$  class(psb_s_coo_sparse_mat), intent(inout) :: b
-  integer(psb_ipk_), intent(in)                  :: imin,imax
+!!$  integer(psb_ipk_), intent(in)                  :: imin,imax
-  integer(psb_ipk_),intent(out)                  :: info
+!!$  integer(psb_ipk_),intent(out)                  :: info
-  logical, intent(in), optional        :: append
+!!$  logical, intent(in), optional        :: append
-  integer(psb_ipk_), intent(in), optional        :: iren(:)
+!!$  integer(psb_ipk_), intent(in), optional        :: iren(:)
-  integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+!!$  integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-  logical, intent(in), optional        :: rscale,cscale
+!!$  logical, intent(in), optional        :: rscale,cscale
-  integer(psb_ipk_) :: err_act, nzin, nzout
+!!$  integer(psb_ipk_) :: err_act, nzin, nzout
-  integer(psb_ipk_) :: ierr(5)
+!!$  integer(psb_ipk_) :: ierr(5)
-  character(len=20)  :: name='csget'
+!!$  character(len=20)  :: name='csget'
-  logical :: append_
+!!$  logical :: append_
-  logical, parameter :: debug=.false.
+!!$  logical, parameter :: debug=.false.
-
+!!$
-  call psb_erractionsave(err_act)
+!!$  call psb_erractionsave(err_act)
-  info = psb_success_
+!!$  info = psb_success_
-
+!!$
-  if (present(append)) then 
+!!$  if (present(append)) then 
-    append_ = append
+!!$    append_ = append
-  else
+!!$  else
-    append_ = .false.
+!!$    append_ = .false.
-  endif
+!!$  endif
-  if (append_) then 
+!!$  if (append_) then 
-    nzin = a%get_nzeros()
+!!$    nzin = a%get_nzeros()
-  else
+!!$  else
-    nzin = 0
+!!$    nzin = 0
-  endif
+!!$  endif
-
+!!$
-  call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
+!!$  call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
-       & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
+!!$       & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
-       & nzin=nzin, rscale=rscale, cscale=cscale)
+!!$       & nzin=nzin, rscale=rscale, cscale=cscale)
-
+!!$
-  if (info /= psb_success_) goto 9999
+!!$  if (info /= psb_success_) goto 9999
-
+!!$
-  call b%set_nzeros(nzin+nzout)
+!!$  call b%set_nzeros(nzin+nzout)
-  call b%fix(info)
+!!$  call b%fix(info)
-  if (info /= psb_success_) goto 9999
+!!$  if (info /= psb_success_) goto 9999
-
+!!$
-  call psb_erractionrestore(err_act)
+!!$  call psb_erractionrestore(err_act)
-  return
+!!$  return
-
+!!$
-9999 call psb_error_handler(err_act)
+!!$9999 call psb_error_handler(err_act)
-
+!!$
-  return
+!!$  return
-
+!!$
-end subroutine psb_s_csc_csgetblk
+!!$end subroutine psb_s_csc_csgetblk
 
 subroutine psb_s_csc_reinit(a,clear)
   use psb_error_mod

base/serial/impl/psb_s_csr_impl.f90

@@ -1999,7 +1999,7 @@ end subroutine psb_s_csr_csgetptn
 
 
 subroutine psb_s_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -2017,9 +2017,9 @@ subroutine psb_s_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
-  logical :: append_, rscale_, cscale_ 
+  logical :: append_, rscale_, cscale_, chksz_
   integer(psb_ipk_) :: nzin_, jmin_, jmax_, err_act, i
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -2063,13 +2063,18 @@ subroutine psb_s_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   else
     cscale_ = .false.
   endif
+  if (present(chksz)) then 
+    chksz_ = chksz
+  else
+    chksz_ = .true.
+  endif
   if ((rscale_.or.cscale_).and.(present(iren))) then 
     info = psb_err_many_optional_arg_
     call psb_errpush(info,name,a_err='iren (rscale.or.cscale)')
     goto 9999
   end if
 
-  call csr_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,info,&
+  call csr_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,chksz_,info,&
        & iren)
   
   if (rscale_) then 
@@ -2094,7 +2099,7 @@ subroutine psb_s_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
 
 contains
 
-  subroutine csr_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,info,&
+  subroutine csr_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,chksz,info,&
        & iren)
 
     use psb_const_mod
@@ -2109,7 +2114,7 @@ contains
     integer(psb_ipk_), allocatable, intent(inout)  :: ia(:), ja(:)
     real(psb_spk_), allocatable,  intent(inout)    :: val(:)
     integer(psb_ipk_), intent(in)                  :: nzin
-    logical, intent(in)                  :: append
+    logical, intent(in)                  :: append, chksz
     integer(psb_ipk_) :: info
     integer(psb_ipk_), optional                    :: iren(:)
     integer(psb_ipk_) :: nzin_, nza, idx,i,j,k, nzt, irw, lrw, icl,lcl, nrd, ncd
@@ -2142,11 +2147,13 @@ contains
     nzt = (a%irp(lrw+1)-a%irp(irw))
     nz = 0 
 
-    call psb_ensure_size(nzin_+nzt,ia,info)
+    if (chksz) then 
-    if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+      call psb_ensure_size(nzin_+nzt,ia,info)
-    if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-
+      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-    if (info /= psb_success_) return
+      
+      if (info /= psb_success_) return
+    end if
     
     if (present(iren)) then 
       do i=irw, lrw
@@ -2178,50 +2185,155 @@ contains
 
 end subroutine psb_s_csr_csgetrow
 
-subroutine psb_s_csr_csgetblk(imin,imax,a,b,info,&
+
-     & jmin,jmax,iren,append,rscale,cscale)
+!
+! CSR implementation of tril/triu
+!
+subroutine psb_s_csr_tril(a,l,info,&
+     & diag,imin,imax,jmin,jmax,rscale,cscale,u)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
-  use psb_s_csr_mat_mod, psb_protect_name => psb_s_csr_csgetblk
+  use psb_s_csr_mat_mod, psb_protect_name => psb_s_csr_tril
   implicit none
 
   class(psb_s_csr_sparse_mat), intent(in) :: a
-  class(psb_s_coo_sparse_mat), intent(inout) :: b
+  class(psb_s_coo_sparse_mat), intent(out) :: l
-  integer(psb_ipk_), intent(in)                  :: imin,imax
+  integer(psb_ipk_),intent(out)            :: info
-  integer(psb_ipk_),intent(out)                  :: info
+  integer(psb_ipk_), intent(in), optional  :: diag,imin,imax,jmin,jmax
-  logical, intent(in), optional        :: append
+  logical, intent(in), optional            :: rscale,cscale
-  integer(psb_ipk_), intent(in), optional        :: iren(:)
+  class(psb_s_coo_sparse_mat), optional, intent(out) :: u
-  integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+  
-  logical, intent(in), optional        :: rscale,cscale
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k
-  integer(psb_ipk_) :: err_act, nzin, nzout
+  integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
+  integer(psb_ipk_), allocatable :: ia(:), ja(:)
+  real(psb_spk_), allocatable    :: val(:)
   integer(psb_ipk_) :: ierr(5)
-  character(len=20)  :: name='csget'
+  character(len=20)  :: name='tril'
-  logical :: append_
+  logical :: rscale_, cscale_
   logical, parameter :: debug=.false.
 
   call psb_erractionsave(err_act)
   info = psb_success_
 
-  if (present(append)) then 
+  if (present(diag)) then 
-    append_ = append
+    diag_ = diag
   else
-    append_ = .false.
+    diag_ = 0
-  endif
+  end if
-  if (append_) then 
+  if (present(imin)) then 
-    nzin = a%get_nzeros()
+    imin_ = imin
+  else
+    imin_ = 1
+  end if
+  if (present(imax)) then 
+    imax_ = imax
   else
-    nzin = 0
+    imax_ = a%get_nrows()
+  end if
+  if (present(jmin)) then 
+    jmin_ = jmin
+  else
+    jmin_ = 1
+  end if
+  if (present(jmax)) then 
+    jmax_ = jmax
+  else
+    jmax_ = a%get_ncols()
+  end if
+  if (present(rscale)) then 
+    rscale_ = rscale
+  else
+    rscale_ = .true.
+  end if
+  if (present(cscale)) then 
+    cscale_ = cscale
+  else
+    cscale_ = .true.
+  end if
+
+  if (rscale_) then 
+    mb = imax_ - imin_ +1
+  else 
+    mb = imax_ 
+  endif
+  if (cscale_) then 
+    nb = jmax_ - jmin_ +1
+  else 
+    nb = jmax_  
   endif
 
-  call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
-       & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
-       & nzin=nzin, rscale=rscale, cscale=cscale)
 
-  if (info /= psb_success_) goto 9999
+  nz = a%get_nzeros()
+  call l%allocate(mb,nb,nz)
+  
+  if (present(u)) then
+    nzlin = l%get_nzeros() ! At this point it should be 0
+    call u%allocate(mb,nb,nz)
+    nzuin = u%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_
+        do k=irp(i),irp(i+1)-1
+          j = ja(k)
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)<=diag_) then
+              nzlin = nzlin + 1
+              l%ia(nzlin)  = i
+              l%ja(nzlin)  = ja(k)
+              l%val(nzlin) = val(k)
+            else
+              nzuin = nzuin + 1
+              u%ia(nzuin)  = i
+              u%ja(nzuin)  = ja(k)
+              u%val(nzuin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    
+    call l%set_nzeros(nzlin)
+    call u%set_nzeros(nzuin)
+    call u%fix(info)
+    nzout = u%get_nzeros()
+    if (rscale_) &
+         & u%ia(1:nzout) = u%ia(1:nzout) - imin_ + 1
+    if (cscale_) &
+         & u%ja(1:nzout) = u%ja(1:nzout) - jmin_ + 1  
+    if ((diag_ >=-1).and.(imin_ == jmin_)) then 
+      call u%set_triangle(.true.)
+      call u%set_lower(.false.)
+    end if
+  else
+    nzin = l%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_        
+        do k=irp(i),irp(i+1)-1
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)<=diag_) then
+              nzin = nzin + 1
+              l%ia(nzin)  = i
+              l%ja(nzin)  = ja(k)
+              l%val(nzin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    call l%set_nzeros(nzin)
+  end if
+  call l%fix(info)
+  nzout = l%get_nzeros()
+  if (rscale_) &
+       & l%ia(1:nzout) = l%ia(1:nzout) - imin_ + 1
+  if (cscale_) &
+       & l%ja(1:nzout) = l%ja(1:nzout) - jmin_ + 1
+  
+  if ((diag_ <= 0).and.(imin_ == jmin_)) then 
+    call l%set_triangle(.true.)
+    call l%set_lower(.true.)
+  end if
 
-  call b%set_nzeros(nzin+nzout)
-  call b%fix(info)
   if (info /= psb_success_) goto 9999
 
   call psb_erractionrestore(err_act)
@@ -2231,8 +2343,162 @@ subroutine psb_s_csr_csgetblk(imin,imax,a,b,info,&
 
   return
 
-end subroutine psb_s_csr_csgetblk
+end subroutine psb_s_csr_tril
+
+subroutine psb_s_csr_triu(a,u,info,&
+     & diag,imin,imax,jmin,jmax,rscale,cscale,l)
+  ! Output is always in  COO format 
+  use psb_error_mod
+  use psb_const_mod
+  use psb_s_csr_mat_mod, psb_protect_name => psb_s_csr_triu
+  implicit none
+
+  class(psb_s_csr_sparse_mat), intent(in) :: a
+  class(psb_s_coo_sparse_mat), intent(out) :: u
+  integer(psb_ipk_),intent(out)            :: info
+  integer(psb_ipk_), intent(in), optional  :: diag,imin,imax,jmin,jmax
+  logical, intent(in), optional            :: rscale,cscale
+  class(psb_s_coo_sparse_mat), optional, intent(out) :: l
+  
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k
+  integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
+  integer(psb_ipk_), allocatable :: ia(:), ja(:)
+  real(psb_spk_), allocatable    :: val(:)
+  integer(psb_ipk_) :: ierr(5)
+  character(len=20)  :: name='triu'
+  logical :: rscale_, cscale_
+  logical, parameter :: debug=.false.
+
+  call psb_erractionsave(err_act)
+  info = psb_success_
+
+  if (present(diag)) then 
+    diag_ = diag
+  else
+    diag_ = 0
+  end if
+  if (present(imin)) then 
+    imin_ = imin
+  else
+    imin_ = 1
+  end if
+  if (present(imax)) then 
+    imax_ = imax
+  else
+    imax_ = a%get_nrows()
+  end if
+  if (present(jmin)) then 
+    jmin_ = jmin
+  else
+    jmin_ = 1
+  end if
+  if (present(jmax)) then 
+    jmax_ = jmax
+  else
+    jmax_ = a%get_ncols()
+  end if
+  if (present(rscale)) then 
+    rscale_ = rscale
+  else
+    rscale_ = .true.
+  end if
+  if (present(cscale)) then 
+    cscale_ = cscale
+  else
+    cscale_ = .true.
+  end if
+
+  if (rscale_) then 
+    mb = imax_ - imin_ +1
+  else 
+    mb = imax_ 
+  endif
+  if (cscale_) then 
+    nb = jmax_ - jmin_ +1
+  else 
+    nb = jmax_  
+  endif
+
+
+  nz = a%get_nzeros()
+  call u%allocate(mb,nb,nz)
+  
+  if (present(l)) then
+    nzuin = u%get_nzeros() ! At this point it should be 0
+    call l%allocate(mb,nb,nz)
+    nzlin = l%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_
+        do k=irp(i),irp(i+1)-1
+          j = ja(k)
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)<diag_) then
+              nzlin = nzlin + 1
+              l%ia(nzlin)  = i
+              l%ja(nzlin)  = ja(k)
+              l%val(nzlin) = val(k)
+            else
+              nzuin = nzuin + 1
+              u%ia(nzuin)  = i
+              u%ja(nzuin)  = ja(k)
+              u%val(nzuin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    call u%set_nzeros(nzuin)
+    call l%set_nzeros(nzlin)
+    call l%fix(info)
+    nzout = l%get_nzeros()
+    if (rscale_) &
+         & l%ia(1:nzout) = l%ia(1:nzout) - imin_ + 1
+    if (cscale_) &
+         & l%ja(1:nzout) = l%ja(1:nzout) - jmin_ + 1  
+    if ((diag_ <=1).and.(imin_ == jmin_)) then 
+      call l%set_triangle(.true.)
+      call l%set_lower(.true.)
+    end if
+  else
+    nzin = u%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_        
+        do k=irp(i),irp(i+1)-1
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)>=diag_) then
+              nzin = nzin + 1
+              u%ia(nzin)  = i
+              u%ja(nzin)  = ja(k)
+              u%val(nzin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    call u%set_nzeros(nzin)
+  end if
+  call u%fix(info)
+  nzout = u%get_nzeros()
+  if (rscale_) &
+       & u%ia(1:nzout) = u%ia(1:nzout) - imin_ + 1
+  if (cscale_) &
+       & u%ja(1:nzout) = u%ja(1:nzout) - jmin_ + 1
+  
+  if ((diag_ >= 0).and.(imin_ == jmin_)) then 
+    call u%set_triangle(.true.)
+    call u%set_upper(.true.)
+  end if
+
+  if (info /= psb_success_) goto 9999
+
+  call psb_erractionrestore(err_act)
+  return
+
+9999 call psb_error_handler(err_act)
+
+  return
 
+end subroutine psb_s_csr_triu
 
 
 subroutine psb_s_csr_csput_a(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl) 

base/serial/impl/psb_s_mat_impl.F90

@@ -793,7 +793,7 @@ end subroutine psb_s_csgetptn
 
 
 subroutine psb_s_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -810,7 +810,7 @@ subroutine psb_s_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
   integer(psb_ipk_) :: err_act
   character(len=20)  :: name='csget'
@@ -826,7 +826,7 @@ subroutine psb_s_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
 
 
   call a%a%csget(imin,imax,nz,ia,ja,val,info,&
-       & jmin,jmax,iren,append,nzin,rscale,cscale)
+       & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   if (info /= psb_success_) goto 9999 
 
   call psb_erractionrestore(err_act)

base/serial/impl/psb_z_base_mat_impl.F90

@@ -406,7 +406,7 @@ subroutine psb_z_base_csput_v(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl)
 end subroutine psb_z_base_csput_v
 
 subroutine psb_z_base_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -422,7 +422,7 @@ subroutine psb_z_base_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
   integer(psb_ipk_) :: err_act
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -447,7 +447,7 @@ end subroutine psb_z_base_csgetrow
 ! If performance is critical it can be overridden.
 !
 subroutine psb_z_base_csgetblk(imin,imax,a,b,info,&
-     & jmin,jmax,iren,append,rscale,cscale)
+     & jmin,jmax,iren,append,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -461,7 +461,7 @@ subroutine psb_z_base_csgetblk(imin,imax,a,b,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
   integer(psb_ipk_) :: err_act, nzin, nzout
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -522,7 +522,7 @@ subroutine psb_z_base_csgetblk(imin,imax,a,b,info,&
 
   call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
        & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
-       & nzin=nzin, rscale=rscale, cscale=cscale)
+       & nzin=nzin, rscale=rscale, cscale=cscale, chksz=chksz)
 
   if (info /= psb_success_) goto 9999
 
@@ -645,7 +645,7 @@ subroutine psb_z_base_tril(a,l,info,&
   logical, intent(in), optional            :: rscale,cscale
   class(psb_z_coo_sparse_mat), optional, intent(out) :: u
   
-  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k 
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k, ibk 
   integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
   integer(psb_ipk_), allocatable :: ia(:), ja(:)
   complex(psb_dpk_), allocatable    :: val(:)
@@ -653,6 +653,7 @@ subroutine psb_z_base_tril(a,l,info,&
   character(len=20)  :: name='tril'
   logical :: rscale_, cscale_
   logical, parameter :: debug=.false.
+  integer(psb_ipk_), parameter :: nbk=8
 
   call psb_erractionsave(err_act)
   info = psb_success_
@@ -715,12 +716,12 @@ subroutine psb_z_base_tril(a,l,info,&
     call psb_realloc(max(mb,nb),ia,info)
     call psb_realloc(max(mb,nb),ja,info)
     call psb_realloc(max(mb,nb),val,info)
-    do i=imin_,imax_
+    do i=imin_,imax_, nbk
-      call a%csget(i,i,nzout,ia,ja,val,info,&
+      ibk = min(nbk,imax_-i+1)
+      call a%csget(i,i+ibk-1,nzout,ia,ja,val,info,&
            & jmin=jmin_, jmax=jmax_)
       do k=1, nzout
-        j = ja(k) 
+        if ((ja(k)-ia(k))<=diag_) then
-        if (j-i<=diag_) then
           nzlin = nzlin + 1
           l%ia(nzlin)  = ia(k)
           l%ja(nzlin)  = ja(k)
@@ -796,7 +797,7 @@ subroutine psb_z_base_triu(a,u,info,&
   logical, intent(in), optional            :: rscale,cscale
   class(psb_z_coo_sparse_mat), optional, intent(out) :: l
   
-  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k 
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k, ibk
   integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
   integer(psb_ipk_), allocatable :: ia(:), ja(:)
   complex(psb_dpk_), allocatable    :: val(:)
@@ -804,6 +805,7 @@ subroutine psb_z_base_triu(a,u,info,&
   character(len=20)  :: name='triu'
   logical :: rscale_, cscale_
   logical, parameter :: debug=.false.
+  integer(psb_ipk_), parameter :: nbk=8
 
   call psb_erractionsave(err_act)
   info = psb_success_
@@ -866,12 +868,12 @@ subroutine psb_z_base_triu(a,u,info,&
     call psb_realloc(max(mb,nb),ia,info)
     call psb_realloc(max(mb,nb),ja,info)
     call psb_realloc(max(mb,nb),val,info)
-    do i=imin_,imax_
+    do i=imin_,imax_, nbk
-      call a%csget(i,i,nzout,ia,ja,val,info,&
+      ibk = min(nbk,imax_-i+1)
+      call a%csget(i,i+ibk-1,nzout,ia,ja,val,info,&
            & jmin=jmin_, jmax=jmax_)
       do k=1, nzout
-        j = ja(k) 
+        if ((ja(k)-ia(k))<diag_) then
-        if (j-i<diag_) then
           nzlin = nzlin + 1
           l%ia(nzlin)  = ia(k)
           l%ja(nzlin)  = ja(k)

base/serial/impl/psb_z_coo_impl.f90

@@ -2259,7 +2259,7 @@ end subroutine psb_z_coo_csgetptn
 ! The output is guaranteed to be sorted
 ! 
 subroutine psb_z_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -2276,9 +2276,9 @@ subroutine psb_z_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
-  logical :: append_, rscale_, cscale_ 
+  logical :: append_, rscale_, cscale_, chksz_
   integer(psb_ipk_) :: nzin_, jmin_, jmax_, err_act, i
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -2321,13 +2321,18 @@ subroutine psb_z_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   else
     cscale_ = .false.
   endif
+  if (present(chksz)) then 
+    chksz_ = chksz
+  else
+    chksz_ = .true.
+  endif
   if ((rscale_.or.cscale_).and.(present(iren))) then 
     info = psb_err_many_optional_arg_
     call psb_errpush(info,name,a_err='iren (rscale.or.cscale)')
     goto 9999
   end if
 
-  call coo_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,info,&
+  call coo_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,chksz_,info,&
        & iren)
 
   if (rscale_) then 
@@ -2352,7 +2357,7 @@ subroutine psb_z_coo_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
 
 contains
 
-  subroutine coo_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,info,&
+  subroutine coo_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,chksz,info,&
        & iren)
 
     use psb_const_mod
@@ -2368,7 +2373,7 @@ contains
     integer(psb_ipk_), allocatable, intent(inout)  :: ia(:), ja(:)
     complex(psb_dpk_), allocatable,  intent(inout)    :: val(:)
     integer(psb_ipk_), intent(in)                  :: nzin
-    logical, intent(in)                            :: append
+    logical, intent(in)                            :: append,chksz
     integer(psb_ipk_) :: info
     integer(psb_ipk_), optional                    :: iren(:)
     integer(psb_ipk_) :: nzin_, nza, idx,ip,jp,i,k, nzt, irw, lrw, nra, nca, nrd
@@ -2452,11 +2457,13 @@ contains
         nzt = jp - ip +1 
         nz = 0 
 
-        call psb_ensure_size(nzin_+nzt,ia,info)
+        if (chksz) then 
-        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+          call psb_ensure_size(nzin_+nzt,ia,info)
-        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+          if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-        if (info /= psb_success_) return
+          if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-
+          if (info /= psb_success_) return
+        end if
+        
         if (present(iren)) then 
           do i=ip,jp
             if ((jmin <= a%ja(i)).and.(a%ja(i)<=jmax)) then 
@@ -2488,11 +2495,13 @@ contains
 
       nrd = max(a%get_nrows(),1)
       nzt = ((nza+nrd-1)/nrd)*(lrw-irw+1)
-      call psb_ensure_size(nzin_+nzt,ia,info)
+      if (chksz) then 
-      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+        call psb_ensure_size(nzin_+nzt,ia,info)
-      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-      if (info /= psb_success_) return
+        if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-
+        if (info /= psb_success_) return
+      end if
+      
       if (present(iren)) then 
         k = 0 
         do i=1, a%get_nzeros()
@@ -2501,10 +2510,12 @@ contains
             k = k + 1 
             if (k > nzt) then
               nzt = k + nzt
-              call psb_ensure_size(nzin_+nzt,ia,info)
+              if (chksz) then 
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+                call psb_ensure_size(nzin_+nzt,ia,info)
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-              if (info /= psb_success_) return
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+                if (info /= psb_success_) return
+              end if              
             end if
             val(nzin_+k) = a%val(i)
             ia(nzin_+k)  = iren(a%ia(i))
@@ -2519,11 +2530,12 @@ contains
             k = k + 1 
             if (k > nzt) then
               nzt = k + nzt
-              call psb_ensure_size(nzin_+nzt,ia,info)
+              if (chksz) then 
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+                call psb_ensure_size(nzin_+nzt,ia,info)
-              if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-              if (info /= psb_success_) return
+                if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-
+                if (info /= psb_success_) return
+              end if
             end if
             val(nzin_+k) = a%val(i)
             ia(nzin_+k)  = (a%ia(i))

base/serial/impl/psb_z_csc_impl.f90

@@ -1686,7 +1686,7 @@ end subroutine psb_z_csc_csgetptn
 
 
 subroutine psb_z_csc_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -1704,7 +1704,7 @@ subroutine psb_z_csc_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
   logical :: append_, rscale_, cscale_ 
   integer(psb_ipk_) :: nzin_, jmin_, jmax_, err_act, i
@@ -2557,60 +2557,60 @@ end subroutine psb_z_csc_reallocate_nz
 
 
 
-subroutine psb_z_csc_csgetblk(imin,imax,a,b,info,&
+!!$subroutine psb_z_csc_csgetblk(imin,imax,a,b,info,&
-     & jmin,jmax,iren,append,rscale,cscale)
+!!$     & jmin,jmax,iren,append,rscale,cscale)
-  ! Output is always in  COO format 
+!!$  ! Output is always in  COO format 
-  use psb_error_mod
+!!$  use psb_error_mod
-  use psb_const_mod
+!!$  use psb_const_mod
-  use psb_z_csc_mat_mod, psb_protect_name => psb_z_csc_csgetblk
+!!$  use psb_z_csc_mat_mod, psb_protect_name => psb_z_csc_csgetblk
-  implicit none
+!!$  implicit none
-
+!!$
-  class(psb_z_csc_sparse_mat), intent(in) :: a
+!!$  class(psb_z_csc_sparse_mat), intent(in) :: a
-  class(psb_z_coo_sparse_mat), intent(inout) :: b
+!!$  class(psb_z_coo_sparse_mat), intent(inout) :: b
-  integer(psb_ipk_), intent(in)                  :: imin,imax
+!!$  integer(psb_ipk_), intent(in)                  :: imin,imax
-  integer(psb_ipk_),intent(out)                  :: info
+!!$  integer(psb_ipk_),intent(out)                  :: info
-  logical, intent(in), optional        :: append
+!!$  logical, intent(in), optional        :: append
-  integer(psb_ipk_), intent(in), optional        :: iren(:)
+!!$  integer(psb_ipk_), intent(in), optional        :: iren(:)
-  integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+!!$  integer(psb_ipk_), intent(in), optional        :: jmin,jmax
-  logical, intent(in), optional        :: rscale,cscale
+!!$  logical, intent(in), optional        :: rscale,cscale
-  integer(psb_ipk_) :: err_act, nzin, nzout
+!!$  integer(psb_ipk_) :: err_act, nzin, nzout
-  integer(psb_ipk_) :: ierr(5)
+!!$  integer(psb_ipk_) :: ierr(5)
-  character(len=20)  :: name='csget'
+!!$  character(len=20)  :: name='csget'
-  logical :: append_
+!!$  logical :: append_
-  logical, parameter :: debug=.false.
+!!$  logical, parameter :: debug=.false.
-
+!!$
-  call psb_erractionsave(err_act)
+!!$  call psb_erractionsave(err_act)
-  info = psb_success_
+!!$  info = psb_success_
-
+!!$
-  if (present(append)) then 
+!!$  if (present(append)) then 
-    append_ = append
+!!$    append_ = append
-  else
+!!$  else
-    append_ = .false.
+!!$    append_ = .false.
-  endif
+!!$  endif
-  if (append_) then 
+!!$  if (append_) then 
-    nzin = a%get_nzeros()
+!!$    nzin = a%get_nzeros()
-  else
+!!$  else
-    nzin = 0
+!!$    nzin = 0
-  endif
+!!$  endif
-
+!!$
-  call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
+!!$  call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
-       & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
+!!$       & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
-       & nzin=nzin, rscale=rscale, cscale=cscale)
+!!$       & nzin=nzin, rscale=rscale, cscale=cscale)
-
+!!$
-  if (info /= psb_success_) goto 9999
+!!$  if (info /= psb_success_) goto 9999
-
+!!$
-  call b%set_nzeros(nzin+nzout)
+!!$  call b%set_nzeros(nzin+nzout)
-  call b%fix(info)
+!!$  call b%fix(info)
-  if (info /= psb_success_) goto 9999
+!!$  if (info /= psb_success_) goto 9999
-
+!!$
-  call psb_erractionrestore(err_act)
+!!$  call psb_erractionrestore(err_act)
-  return
+!!$  return
-
+!!$
-9999 call psb_error_handler(err_act)
+!!$9999 call psb_error_handler(err_act)
-
+!!$
-  return
+!!$  return
-
+!!$
-end subroutine psb_z_csc_csgetblk
+!!$end subroutine psb_z_csc_csgetblk
 
 subroutine psb_z_csc_reinit(a,clear)
   use psb_error_mod

base/serial/impl/psb_z_csr_impl.f90

@@ -1999,7 +1999,7 @@ end subroutine psb_z_csr_csgetptn
 
 
 subroutine psb_z_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -2017,9 +2017,9 @@ subroutine psb_z_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
-  logical :: append_, rscale_, cscale_ 
+  logical :: append_, rscale_, cscale_, chksz_
   integer(psb_ipk_) :: nzin_, jmin_, jmax_, err_act, i
   integer(psb_ipk_) :: ierr(5)
   character(len=20)  :: name='csget'
@@ -2063,13 +2063,18 @@ subroutine psb_z_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   else
     cscale_ = .false.
   endif
+  if (present(chksz)) then 
+    chksz_ = chksz
+  else
+    chksz_ = .true.
+  endif
   if ((rscale_.or.cscale_).and.(present(iren))) then 
     info = psb_err_many_optional_arg_
     call psb_errpush(info,name,a_err='iren (rscale.or.cscale)')
     goto 9999
   end if
 
-  call csr_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,info,&
+  call csr_getrow(imin,imax,jmin_,jmax_,a,nz,ia,ja,val,nzin_,append_,chksz_,info,&
        & iren)
   
   if (rscale_) then 
@@ -2094,7 +2099,7 @@ subroutine psb_z_csr_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
 
 contains
 
-  subroutine csr_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,info,&
+  subroutine csr_getrow(imin,imax,jmin,jmax,a,nz,ia,ja,val,nzin,append,chksz,info,&
        & iren)
 
     use psb_const_mod
@@ -2109,7 +2114,7 @@ contains
     integer(psb_ipk_), allocatable, intent(inout)  :: ia(:), ja(:)
     complex(psb_dpk_), allocatable,  intent(inout)    :: val(:)
     integer(psb_ipk_), intent(in)                  :: nzin
-    logical, intent(in)                  :: append
+    logical, intent(in)                  :: append, chksz
     integer(psb_ipk_) :: info
     integer(psb_ipk_), optional                    :: iren(:)
     integer(psb_ipk_) :: nzin_, nza, idx,i,j,k, nzt, irw, lrw, icl,lcl, nrd, ncd
@@ -2142,11 +2147,13 @@ contains
     nzt = (a%irp(lrw+1)-a%irp(irw))
     nz = 0 
 
-    call psb_ensure_size(nzin_+nzt,ia,info)
+    if (chksz) then 
-    if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
+      call psb_ensure_size(nzin_+nzt,ia,info)
-    if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
+      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,ja,info)
-
+      if (info == psb_success_) call psb_ensure_size(nzin_+nzt,val,info)
-    if (info /= psb_success_) return
+      
+      if (info /= psb_success_) return
+    end if
     
     if (present(iren)) then 
       do i=irw, lrw
@@ -2178,50 +2185,155 @@ contains
 
 end subroutine psb_z_csr_csgetrow
 
-subroutine psb_z_csr_csgetblk(imin,imax,a,b,info,&
+
-     & jmin,jmax,iren,append,rscale,cscale)
+!
+! CSR implementation of tril/triu
+!
+subroutine psb_z_csr_tril(a,l,info,&
+     & diag,imin,imax,jmin,jmax,rscale,cscale,u)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
-  use psb_z_csr_mat_mod, psb_protect_name => psb_z_csr_csgetblk
+  use psb_z_csr_mat_mod, psb_protect_name => psb_z_csr_tril
   implicit none
 
   class(psb_z_csr_sparse_mat), intent(in) :: a
-  class(psb_z_coo_sparse_mat), intent(inout) :: b
+  class(psb_z_coo_sparse_mat), intent(out) :: l
-  integer(psb_ipk_), intent(in)                  :: imin,imax
+  integer(psb_ipk_),intent(out)            :: info
-  integer(psb_ipk_),intent(out)                  :: info
+  integer(psb_ipk_), intent(in), optional  :: diag,imin,imax,jmin,jmax
-  logical, intent(in), optional        :: append
+  logical, intent(in), optional            :: rscale,cscale
-  integer(psb_ipk_), intent(in), optional        :: iren(:)
+  class(psb_z_coo_sparse_mat), optional, intent(out) :: u
-  integer(psb_ipk_), intent(in), optional        :: jmin,jmax
+  
-  logical, intent(in), optional        :: rscale,cscale
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k
-  integer(psb_ipk_) :: err_act, nzin, nzout
+  integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
+  integer(psb_ipk_), allocatable :: ia(:), ja(:)
+  complex(psb_dpk_), allocatable    :: val(:)
   integer(psb_ipk_) :: ierr(5)
-  character(len=20)  :: name='csget'
+  character(len=20)  :: name='tril'
-  logical :: append_
+  logical :: rscale_, cscale_
   logical, parameter :: debug=.false.
 
   call psb_erractionsave(err_act)
   info = psb_success_
 
-  if (present(append)) then 
+  if (present(diag)) then 
-    append_ = append
+    diag_ = diag
   else
-    append_ = .false.
+    diag_ = 0
-  endif
+  end if
-  if (append_) then 
+  if (present(imin)) then 
-    nzin = a%get_nzeros()
+    imin_ = imin
+  else
+    imin_ = 1
+  end if
+  if (present(imax)) then 
+    imax_ = imax
   else
-    nzin = 0
+    imax_ = a%get_nrows()
+  end if
+  if (present(jmin)) then 
+    jmin_ = jmin
+  else
+    jmin_ = 1
+  end if
+  if (present(jmax)) then 
+    jmax_ = jmax
+  else
+    jmax_ = a%get_ncols()
+  end if
+  if (present(rscale)) then 
+    rscale_ = rscale
+  else
+    rscale_ = .true.
+  end if
+  if (present(cscale)) then 
+    cscale_ = cscale
+  else
+    cscale_ = .true.
+  end if
+
+  if (rscale_) then 
+    mb = imax_ - imin_ +1
+  else 
+    mb = imax_ 
+  endif
+  if (cscale_) then 
+    nb = jmax_ - jmin_ +1
+  else 
+    nb = jmax_  
   endif
 
-  call a%csget(imin,imax,nzout,b%ia,b%ja,b%val,info,&
-       & jmin=jmin, jmax=jmax, iren=iren, append=append_, &
-       & nzin=nzin, rscale=rscale, cscale=cscale)
 
-  if (info /= psb_success_) goto 9999
+  nz = a%get_nzeros()
+  call l%allocate(mb,nb,nz)
+  
+  if (present(u)) then
+    nzlin = l%get_nzeros() ! At this point it should be 0
+    call u%allocate(mb,nb,nz)
+    nzuin = u%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_
+        do k=irp(i),irp(i+1)-1
+          j = ja(k)
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)<=diag_) then
+              nzlin = nzlin + 1
+              l%ia(nzlin)  = i
+              l%ja(nzlin)  = ja(k)
+              l%val(nzlin) = val(k)
+            else
+              nzuin = nzuin + 1
+              u%ia(nzuin)  = i
+              u%ja(nzuin)  = ja(k)
+              u%val(nzuin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    
+    call l%set_nzeros(nzlin)
+    call u%set_nzeros(nzuin)
+    call u%fix(info)
+    nzout = u%get_nzeros()
+    if (rscale_) &
+         & u%ia(1:nzout) = u%ia(1:nzout) - imin_ + 1
+    if (cscale_) &
+         & u%ja(1:nzout) = u%ja(1:nzout) - jmin_ + 1  
+    if ((diag_ >=-1).and.(imin_ == jmin_)) then 
+      call u%set_triangle(.true.)
+      call u%set_lower(.false.)
+    end if
+  else
+    nzin = l%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_        
+        do k=irp(i),irp(i+1)-1
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)<=diag_) then
+              nzin = nzin + 1
+              l%ia(nzin)  = i
+              l%ja(nzin)  = ja(k)
+              l%val(nzin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    call l%set_nzeros(nzin)
+  end if
+  call l%fix(info)
+  nzout = l%get_nzeros()
+  if (rscale_) &
+       & l%ia(1:nzout) = l%ia(1:nzout) - imin_ + 1
+  if (cscale_) &
+       & l%ja(1:nzout) = l%ja(1:nzout) - jmin_ + 1
+  
+  if ((diag_ <= 0).and.(imin_ == jmin_)) then 
+    call l%set_triangle(.true.)
+    call l%set_lower(.true.)
+  end if
 
-  call b%set_nzeros(nzin+nzout)
-  call b%fix(info)
   if (info /= psb_success_) goto 9999
 
   call psb_erractionrestore(err_act)
@@ -2231,8 +2343,162 @@ subroutine psb_z_csr_csgetblk(imin,imax,a,b,info,&
 
   return
 
-end subroutine psb_z_csr_csgetblk
+end subroutine psb_z_csr_tril
+
+subroutine psb_z_csr_triu(a,u,info,&
+     & diag,imin,imax,jmin,jmax,rscale,cscale,l)
+  ! Output is always in  COO format 
+  use psb_error_mod
+  use psb_const_mod
+  use psb_z_csr_mat_mod, psb_protect_name => psb_z_csr_triu
+  implicit none
+
+  class(psb_z_csr_sparse_mat), intent(in) :: a
+  class(psb_z_coo_sparse_mat), intent(out) :: u
+  integer(psb_ipk_),intent(out)            :: info
+  integer(psb_ipk_), intent(in), optional  :: diag,imin,imax,jmin,jmax
+  logical, intent(in), optional            :: rscale,cscale
+  class(psb_z_coo_sparse_mat), optional, intent(out) :: l
+  
+  integer(psb_ipk_) :: err_act, nzin, nzout, i, j, k
+  integer(psb_ipk_) :: imin_, imax_, jmin_, jmax_, mb,nb, diag_, nzlin, nzuin, nz
+  integer(psb_ipk_), allocatable :: ia(:), ja(:)
+  complex(psb_dpk_), allocatable    :: val(:)
+  integer(psb_ipk_) :: ierr(5)
+  character(len=20)  :: name='triu'
+  logical :: rscale_, cscale_
+  logical, parameter :: debug=.false.
+
+  call psb_erractionsave(err_act)
+  info = psb_success_
+
+  if (present(diag)) then 
+    diag_ = diag
+  else
+    diag_ = 0
+  end if
+  if (present(imin)) then 
+    imin_ = imin
+  else
+    imin_ = 1
+  end if
+  if (present(imax)) then 
+    imax_ = imax
+  else
+    imax_ = a%get_nrows()
+  end if
+  if (present(jmin)) then 
+    jmin_ = jmin
+  else
+    jmin_ = 1
+  end if
+  if (present(jmax)) then 
+    jmax_ = jmax
+  else
+    jmax_ = a%get_ncols()
+  end if
+  if (present(rscale)) then 
+    rscale_ = rscale
+  else
+    rscale_ = .true.
+  end if
+  if (present(cscale)) then 
+    cscale_ = cscale
+  else
+    cscale_ = .true.
+  end if
+
+  if (rscale_) then 
+    mb = imax_ - imin_ +1
+  else 
+    mb = imax_ 
+  endif
+  if (cscale_) then 
+    nb = jmax_ - jmin_ +1
+  else 
+    nb = jmax_  
+  endif
+
+
+  nz = a%get_nzeros()
+  call u%allocate(mb,nb,nz)
+  
+  if (present(l)) then
+    nzuin = u%get_nzeros() ! At this point it should be 0
+    call l%allocate(mb,nb,nz)
+    nzlin = l%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_
+        do k=irp(i),irp(i+1)-1
+          j = ja(k)
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)<diag_) then
+              nzlin = nzlin + 1
+              l%ia(nzlin)  = i
+              l%ja(nzlin)  = ja(k)
+              l%val(nzlin) = val(k)
+            else
+              nzuin = nzuin + 1
+              u%ia(nzuin)  = i
+              u%ja(nzuin)  = ja(k)
+              u%val(nzuin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    call u%set_nzeros(nzuin)
+    call l%set_nzeros(nzlin)
+    call l%fix(info)
+    nzout = l%get_nzeros()
+    if (rscale_) &
+         & l%ia(1:nzout) = l%ia(1:nzout) - imin_ + 1
+    if (cscale_) &
+         & l%ja(1:nzout) = l%ja(1:nzout) - jmin_ + 1  
+    if ((diag_ <=1).and.(imin_ == jmin_)) then 
+      call l%set_triangle(.true.)
+      call l%set_lower(.true.)
+    end if
+  else
+    nzin = u%get_nzeros() ! At this point it should be 0
+    associate(val =>a%val, ja => a%ja, irp=>a%irp)
+      do i=imin_,imax_        
+        do k=irp(i),irp(i+1)-1
+          if ((jmin_<=j).and.(j<=jmax_)) then 
+            if ((ja(k)-i)>=diag_) then
+              nzin = nzin + 1
+              u%ia(nzin)  = i
+              u%ja(nzin)  = ja(k)
+              u%val(nzin) = val(k)
+            end if
+          end if
+        end do
+      end do
+    end associate
+    call u%set_nzeros(nzin)
+  end if
+  call u%fix(info)
+  nzout = u%get_nzeros()
+  if (rscale_) &
+       & u%ia(1:nzout) = u%ia(1:nzout) - imin_ + 1
+  if (cscale_) &
+       & u%ja(1:nzout) = u%ja(1:nzout) - jmin_ + 1
+  
+  if ((diag_ >= 0).and.(imin_ == jmin_)) then 
+    call u%set_triangle(.true.)
+    call u%set_upper(.true.)
+  end if
+
+  if (info /= psb_success_) goto 9999
+
+  call psb_erractionrestore(err_act)
+  return
+
+9999 call psb_error_handler(err_act)
+
+  return
 
+end subroutine psb_z_csr_triu
 
 
 subroutine psb_z_csr_csput_a(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl) 

base/serial/impl/psb_z_mat_impl.F90

@@ -793,7 +793,7 @@ end subroutine psb_z_csgetptn
 
 
 subroutine psb_z_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
-     & jmin,jmax,iren,append,nzin,rscale,cscale)
+     & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   ! Output is always in  COO format 
   use psb_error_mod
   use psb_const_mod
@@ -810,7 +810,7 @@ subroutine psb_z_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
   logical, intent(in), optional        :: append
   integer(psb_ipk_), intent(in), optional        :: iren(:)
   integer(psb_ipk_), intent(in), optional        :: jmin,jmax, nzin
-  logical, intent(in), optional        :: rscale,cscale
+  logical, intent(in), optional        :: rscale,cscale,chksz
 
   integer(psb_ipk_) :: err_act
   character(len=20)  :: name='csget'
@@ -826,7 +826,7 @@ subroutine psb_z_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
 
 
   call a%a%csget(imin,imax,nz,ia,ja,val,info,&
-       & jmin,jmax,iren,append,nzin,rscale,cscale)
+       & jmin,jmax,iren,append,nzin,rscale,cscale,chksz)
   if (info /= psb_success_) goto 9999 
 
   call psb_erractionrestore(err_act)

base/tools/psb_cdprt.f90

@@ -93,7 +93,11 @@ subroutine psb_cdprt(iout,desc_p,glob,short, verbosity)
     if (me == i)  then
       write(iout,*) me,': Local descriptor data: points:',local_points,&
            & ' halo:',local_halo
-      write(iout,*) me,': Volume to surface ratio:',real(local_points,psb_dpk_)/real(local_halo,psb_dpk_)
+      if (local_halo>0) then 
+        write(iout,*) me,': Volume to surface ratio:',real(local_points,psb_dpk_)/real(local_halo,psb_dpk_)
+      else
+        write(iout,*) me,': Volume to surface ratio:',0.0_psb_dpk_
+      end if
     end if
     call psb_barrier(ictxt)
   end do

base/tools/psb_csphalo.F90

@@ -199,8 +199,7 @@ Subroutine psb_csphalo(a,desc_a,blk,info,rowcnv,colcnv,&
        & rvsz,1,psb_mpi_def_integer,icomm,minfo)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='mpi_alltoall'
+    call psb_errpush(info,name,a_err='mpi_alltoall')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -228,8 +227,7 @@ Subroutine psb_csphalo(a,desc_a,blk,info,rowcnv,colcnv,&
        & ' Send:',sdsz(:),' Receive:',rvsz(:)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psb_sp_reall'
+    call psb_errpush(info,name,a_err='psb_sp_reall')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
   mat_recv = iszr
@@ -237,6 +235,11 @@ Subroutine psb_csphalo(a,desc_a,blk,info,rowcnv,colcnv,&
   call psb_ensure_size(max(iszs,1),iasnd,info)
   if (info == psb_success_) call psb_ensure_size(max(iszs,1),jasnd,info)
   if (info == psb_success_) call psb_ensure_size(max(iszs,1),valsnd,info)
+  if (info /= psb_success_) then
+    info=psb_err_from_subroutine_
+    call psb_errpush(info,name,a_err='ensure_size')
+    goto 9999
+  end if
 
   l1  = 0
   ipx = 1
@@ -258,8 +261,7 @@ Subroutine psb_csphalo(a,desc_a,blk,info,rowcnv,colcnv,&
            &  append=.true.,nzin=tot_elem)
       if (info /= psb_success_) then
         info=psb_err_from_subroutine_
-        ch_err='psb_sp_getrow'
+        call psb_errpush(info,name,a_err='psb_sp_getrow')
-        call psb_errpush(info,name,a_err=ch_err)
         goto 9999
       end if
       tot_elem=tot_elem+n_elem
@@ -273,8 +275,7 @@ Subroutine psb_csphalo(a,desc_a,blk,info,rowcnv,colcnv,&
   if (colcnv_) call psb_loc_to_glob(jasnd(1:nz),desc_a,info,iact='I')
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psb_loc_to_glob'
+    call psb_errpush(info,name,a_err='psb_loc_to_glob')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -287,8 +288,7 @@ Subroutine psb_csphalo(a,desc_a,blk,info,rowcnv,colcnv,&
        & acoo%ja,rvsz,brvindx,psb_mpi_ipk_integer,icomm,minfo)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='mpi_alltoallv'
+    call psb_errpush(info,name,a_err='mpi_alltoallv')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -300,8 +300,7 @@ Subroutine psb_csphalo(a,desc_a,blk,info,rowcnv,colcnv,&
 
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psbglob_to_loc'
+    call psb_errpush(info,name,a_err='psbglob_to_loc')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -353,8 +352,7 @@ Subroutine psb_csphalo(a,desc_a,blk,info,rowcnv,colcnv,&
   call blk%cscnv(info,type=outfmt_,dupl=psb_dupl_add_)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psb_spcnv'
+    call psb_errpush(info,name,a_err='psb_spcnv')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 

base/tools/psb_dsphalo.F90

@@ -199,8 +199,7 @@ Subroutine psb_dsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
        & rvsz,1,psb_mpi_def_integer,icomm,minfo)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='mpi_alltoall'
+    call psb_errpush(info,name,a_err='mpi_alltoall')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -228,8 +227,7 @@ Subroutine psb_dsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
        & ' Send:',sdsz(:),' Receive:',rvsz(:)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psb_sp_reall'
+    call psb_errpush(info,name,a_err='psb_sp_reall')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
   mat_recv = iszr
@@ -237,6 +235,11 @@ Subroutine psb_dsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
   call psb_ensure_size(max(iszs,1),iasnd,info)
   if (info == psb_success_) call psb_ensure_size(max(iszs,1),jasnd,info)
   if (info == psb_success_) call psb_ensure_size(max(iszs,1),valsnd,info)
+  if (info /= psb_success_) then
+    info=psb_err_from_subroutine_
+    call psb_errpush(info,name,a_err='ensure_size')
+    goto 9999
+  end if
 
   l1  = 0
   ipx = 1
@@ -258,8 +261,7 @@ Subroutine psb_dsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
            &  append=.true.,nzin=tot_elem)
       if (info /= psb_success_) then
         info=psb_err_from_subroutine_
-        ch_err='psb_sp_getrow'
+        call psb_errpush(info,name,a_err='psb_sp_getrow')
-        call psb_errpush(info,name,a_err=ch_err)
         goto 9999
       end if
       tot_elem=tot_elem+n_elem
@@ -273,8 +275,7 @@ Subroutine psb_dsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
   if (colcnv_) call psb_loc_to_glob(jasnd(1:nz),desc_a,info,iact='I')
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psb_loc_to_glob'
+    call psb_errpush(info,name,a_err='psb_loc_to_glob')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -287,8 +288,7 @@ Subroutine psb_dsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
        & acoo%ja,rvsz,brvindx,psb_mpi_ipk_integer,icomm,minfo)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='mpi_alltoallv'
+    call psb_errpush(info,name,a_err='mpi_alltoallv')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -300,8 +300,7 @@ Subroutine psb_dsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
 
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psbglob_to_loc'
+    call psb_errpush(info,name,a_err='psbglob_to_loc')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -353,8 +352,7 @@ Subroutine psb_dsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
   call blk%cscnv(info,type=outfmt_,dupl=psb_dupl_add_)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psb_spcnv'
+    call psb_errpush(info,name,a_err='psb_spcnv')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 

base/tools/psb_ssphalo.F90

@@ -199,8 +199,7 @@ Subroutine psb_ssphalo(a,desc_a,blk,info,rowcnv,colcnv,&
        & rvsz,1,psb_mpi_def_integer,icomm,minfo)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='mpi_alltoall'
+    call psb_errpush(info,name,a_err='mpi_alltoall')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -228,8 +227,7 @@ Subroutine psb_ssphalo(a,desc_a,blk,info,rowcnv,colcnv,&
        & ' Send:',sdsz(:),' Receive:',rvsz(:)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psb_sp_reall'
+    call psb_errpush(info,name,a_err='psb_sp_reall')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
   mat_recv = iszr
@@ -237,6 +235,11 @@ Subroutine psb_ssphalo(a,desc_a,blk,info,rowcnv,colcnv,&
   call psb_ensure_size(max(iszs,1),iasnd,info)
   if (info == psb_success_) call psb_ensure_size(max(iszs,1),jasnd,info)
   if (info == psb_success_) call psb_ensure_size(max(iszs,1),valsnd,info)
+  if (info /= psb_success_) then
+    info=psb_err_from_subroutine_
+    call psb_errpush(info,name,a_err='ensure_size')
+    goto 9999
+  end if
 
   l1  = 0
   ipx = 1
@@ -258,8 +261,7 @@ Subroutine psb_ssphalo(a,desc_a,blk,info,rowcnv,colcnv,&
            &  append=.true.,nzin=tot_elem)
       if (info /= psb_success_) then
         info=psb_err_from_subroutine_
-        ch_err='psb_sp_getrow'
+        call psb_errpush(info,name,a_err='psb_sp_getrow')
-        call psb_errpush(info,name,a_err=ch_err)
         goto 9999
       end if
       tot_elem=tot_elem+n_elem
@@ -273,8 +275,7 @@ Subroutine psb_ssphalo(a,desc_a,blk,info,rowcnv,colcnv,&
   if (colcnv_) call psb_loc_to_glob(jasnd(1:nz),desc_a,info,iact='I')
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psb_loc_to_glob'
+    call psb_errpush(info,name,a_err='psb_loc_to_glob')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -287,8 +288,7 @@ Subroutine psb_ssphalo(a,desc_a,blk,info,rowcnv,colcnv,&
        & acoo%ja,rvsz,brvindx,psb_mpi_ipk_integer,icomm,minfo)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='mpi_alltoallv'
+    call psb_errpush(info,name,a_err='mpi_alltoallv')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -300,8 +300,7 @@ Subroutine psb_ssphalo(a,desc_a,blk,info,rowcnv,colcnv,&
 
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psbglob_to_loc'
+    call psb_errpush(info,name,a_err='psbglob_to_loc')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -353,8 +352,7 @@ Subroutine psb_ssphalo(a,desc_a,blk,info,rowcnv,colcnv,&
   call blk%cscnv(info,type=outfmt_,dupl=psb_dupl_add_)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psb_spcnv'
+    call psb_errpush(info,name,a_err='psb_spcnv')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 

base/tools/psb_zsphalo.F90

@@ -199,8 +199,7 @@ Subroutine psb_zsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
        & rvsz,1,psb_mpi_def_integer,icomm,minfo)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='mpi_alltoall'
+    call psb_errpush(info,name,a_err='mpi_alltoall')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -228,8 +227,7 @@ Subroutine psb_zsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
        & ' Send:',sdsz(:),' Receive:',rvsz(:)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psb_sp_reall'
+    call psb_errpush(info,name,a_err='psb_sp_reall')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
   mat_recv = iszr
@@ -237,6 +235,11 @@ Subroutine psb_zsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
   call psb_ensure_size(max(iszs,1),iasnd,info)
   if (info == psb_success_) call psb_ensure_size(max(iszs,1),jasnd,info)
   if (info == psb_success_) call psb_ensure_size(max(iszs,1),valsnd,info)
+  if (info /= psb_success_) then
+    info=psb_err_from_subroutine_
+    call psb_errpush(info,name,a_err='ensure_size')
+    goto 9999
+  end if
 
   l1  = 0
   ipx = 1
@@ -258,8 +261,7 @@ Subroutine psb_zsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
            &  append=.true.,nzin=tot_elem)
       if (info /= psb_success_) then
         info=psb_err_from_subroutine_
-        ch_err='psb_sp_getrow'
+        call psb_errpush(info,name,a_err='psb_sp_getrow')
-        call psb_errpush(info,name,a_err=ch_err)
         goto 9999
       end if
       tot_elem=tot_elem+n_elem
@@ -273,8 +275,7 @@ Subroutine psb_zsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
   if (colcnv_) call psb_loc_to_glob(jasnd(1:nz),desc_a,info,iact='I')
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psb_loc_to_glob'
+    call psb_errpush(info,name,a_err='psb_loc_to_glob')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -287,8 +288,7 @@ Subroutine psb_zsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
        & acoo%ja,rvsz,brvindx,psb_mpi_ipk_integer,icomm,minfo)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='mpi_alltoallv'
+    call psb_errpush(info,name,a_err='mpi_alltoallv')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -300,8 +300,7 @@ Subroutine psb_zsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
 
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psbglob_to_loc'
+    call psb_errpush(info,name,a_err='psbglob_to_loc')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 
@@ -353,8 +352,7 @@ Subroutine psb_zsphalo(a,desc_a,blk,info,rowcnv,colcnv,&
   call blk%cscnv(info,type=outfmt_,dupl=psb_dupl_add_)
   if (info /= psb_success_) then
     info=psb_err_from_subroutine_
-    ch_err='psb_spcnv'
+    call psb_errpush(info,name,a_err='psb_spcnv')
-    call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
 

cbind/base/psb_c_base.h

@@ -1,3 +1,4 @@
+
 #ifndef PSB_C_BASE__
 #define PSB_C_BASE__
 #ifdef __cplusplus

cbind/krylov/psb_ckrylov_cbind_mod.f90

@@ -13,7 +13,7 @@ contains
     use psb_prec_cbind_mod
     use psb_base_string_cbind_mod
     implicit none 
-    integer(c_int)          :: res
+    integer(psb_c_int)          :: res
     type(psb_c_cspmat)    :: ah
     type(psb_c_descriptor)  :: cdh
     type(psb_c_cprec)       :: ph
@@ -38,14 +38,14 @@ contains
     use psb_prec_cbind_mod
     use psb_base_string_cbind_mod
     implicit none 
-    integer(c_int)          :: res
+    integer(psb_c_int)          :: res
     type(psb_c_cspmat)    :: ah
     type(psb_c_descriptor)  :: cdh
     type(psb_c_cprec)       :: ph
     type(psb_c_cvector)     :: bh,xh
-    integer(c_int), value :: itmax,itrace,irst,istop
+    integer(psb_c_int), value :: itmax,itrace,irst,istop
     real(c_double), value :: eps
-    integer(c_int)        :: iter
+    integer(psb_c_int)        :: iter
     real(c_double)        :: err
     character(c_char)       :: methd(*)
     type(solveroptions)     :: options

cbind/krylov/psb_dkrylov_cbind_mod.f90

@@ -13,7 +13,7 @@ contains
     use psb_prec_cbind_mod
     use psb_base_string_cbind_mod
     implicit none 
-    integer(c_int)          :: res
+    integer(psb_c_int)          :: res
     type(psb_c_dspmat)    :: ah
     type(psb_c_descriptor)  :: cdh
     type(psb_c_dprec)       :: ph
@@ -21,7 +21,6 @@ contains
     character(c_char)       :: methd(*)
     type(solveroptions)     :: options
 
-    write(0,*) 'psb_c_dkrylov options ', options%eps
     res= psb_c_dkrylov_opt(methd, ah, ph, bh, xh, options%eps,cdh,  &
          & itmax=options%itmax, iter=options%iter,&
          & itrace=options%itrace, istop=options%istop,&
@@ -39,14 +38,14 @@ contains
     use psb_prec_cbind_mod
     use psb_base_string_cbind_mod
     implicit none 
-    integer(c_int)          :: res
+    integer(psb_c_int)          :: res
     type(psb_c_dspmat)    :: ah
     type(psb_c_descriptor)  :: cdh
     type(psb_c_dprec)       :: ph
     type(psb_c_dvector)     :: bh,xh
-    integer(c_int), value :: itmax,itrace,irst,istop
+    integer(psb_c_int), value :: itmax,itrace,irst,istop
     real(c_double), value :: eps
-    integer(c_int)        :: iter
+    integer(psb_c_int)        :: iter
     real(c_double)        :: err
     character(c_char)       :: methd(*)
     type(solveroptions)     :: options
@@ -59,7 +58,6 @@ contains
     character(len=20)     :: fmethd
     real(psb_dpk_)      :: feps,ferr
 
-    write(0,*) 'psb_c_dkrylov_opt options ', eps
     res = -1
     if (c_associated(cdh%item)) then 
       call c_f_pointer(cdh%item,descp)

cbind/krylov/psb_skrylov_cbind_mod.f90

@@ -13,7 +13,7 @@ contains
     use psb_prec_cbind_mod
     use psb_base_string_cbind_mod
     implicit none 
-    integer(c_int)          :: res
+    integer(psb_c_int)          :: res
     type(psb_c_sspmat)    :: ah
     type(psb_c_descriptor)  :: cdh
     type(psb_c_sprec)       :: ph
@@ -38,14 +38,14 @@ contains
     use psb_prec_cbind_mod
     use psb_base_string_cbind_mod
     implicit none 
-    integer(c_int)          :: res
+    integer(psb_c_int)          :: res
     type(psb_c_sspmat)    :: ah
     type(psb_c_descriptor)  :: cdh
     type(psb_c_sprec)       :: ph
     type(psb_c_svector)     :: bh,xh
-    integer(c_int), value :: itmax,itrace,irst,istop
+    integer(psb_c_int), value :: itmax,itrace,irst,istop
     real(c_double), value :: eps
-    integer(c_int)        :: iter
+    integer(psb_c_int)        :: iter
     real(c_double)        :: err
     character(c_char)       :: methd(*)
     type(solveroptions)     :: options

cbind/krylov/psb_zkrylov_cbind_mod.f90

@@ -13,7 +13,7 @@ contains
     use psb_prec_cbind_mod
     use psb_base_string_cbind_mod
     implicit none 
-    integer(c_int)          :: res
+    integer(psb_c_int)          :: res
     type(psb_c_zspmat)    :: ah
     type(psb_c_descriptor)  :: cdh
     type(psb_c_zprec)       :: ph
@@ -38,14 +38,14 @@ contains
     use psb_prec_cbind_mod
     use psb_base_string_cbind_mod
     implicit none 
-    integer(c_int)          :: res
+    integer(psb_c_int)          :: res
     type(psb_c_zspmat)    :: ah
     type(psb_c_descriptor)  :: cdh
     type(psb_c_zprec)       :: ph
     type(psb_c_zvector)     :: bh,xh
-    integer(c_int), value :: itmax,itrace,irst,istop
+    integer(psb_c_int), value :: itmax,itrace,irst,istop
     real(c_double), value :: eps
-    integer(c_int)        :: iter
+    integer(psb_c_int)        :: iter
     real(c_double)        :: err
     character(c_char)       :: methd(*)
     type(solveroptions)     :: options

cbind/prec/psb_c_cprec.h

@@ -8,15 +8,15 @@
 extern "C" {
 #endif
 
-typedef struct PSB_C_CPREC {
+  typedef struct PSB_C_CPREC {
-  void *cprec;
+    void *cprec;
-} psb_c_cprec; 
+  } psb_c_cprec; 
-
+  
-psb_c_cprec* psb_c_new_cprec(); 
+  psb_c_cprec* psb_c_new_cprec(); 
-
+  
-psb_i_t  psb_c_cprecinit(psb_c_cprec *ph, const char *ptype);
+  psb_i_t  psb_c_cprecinit(psb_i_t ictxt,psb_c_cprec *ph, const char *ptype);
-psb_i_t  psb_c_cprecbld(psb_c_cspmat *ah, psb_c_descriptor *cdh, psb_c_cprec *ph);
+  psb_i_t  psb_c_cprecbld(psb_c_cspmat *ah, psb_c_descriptor *cdh, psb_c_cprec *ph);
-psb_i_t  psb_c_cprecfree(psb_c_cprec *ph);
+  psb_i_t  psb_c_cprecfree(psb_c_cprec *ph);
 #ifdef __cplusplus
 }
 #endif

cbind/prec/psb_c_dprec.h

@@ -8,15 +8,15 @@
 extern "C" {
 #endif
 
-typedef struct PSB_C_DPREC {
+  typedef struct PSB_C_DPREC {
-  void *dprec;
+    void *dprec;
-} psb_c_dprec; 
+  } psb_c_dprec; 
-
+  
-psb_c_dprec* psb_c_new_dprec();
+  psb_c_dprec* psb_c_new_dprec();
-
+  
-psb_i_t  psb_c_dprecinit(psb_c_dprec *ph, const char *ptype);
+  psb_i_t  psb_c_dprecinit(psb_i_t ictxt, psb_c_dprec *ph, const char *ptype);
-psb_i_t  psb_c_dprecbld(psb_c_dspmat *ah, psb_c_descriptor *cdh, psb_c_dprec *ph);
+  psb_i_t  psb_c_dprecbld(psb_c_dspmat *ah, psb_c_descriptor *cdh, psb_c_dprec *ph);
-psb_i_t  psb_c_dprecfree(psb_c_dprec *ph);
+  psb_i_t  psb_c_dprecfree(psb_c_dprec *ph);
 #ifdef __cplusplus
 }
 #endif

cbind/prec/psb_c_sprec.h

@@ -8,15 +8,15 @@
 extern "C" {
 #endif
 
-typedef struct PSB_C_SPREC {
+  typedef struct PSB_C_SPREC {
-  void *sprec;
+    void *sprec;
-} psb_c_sprec; 
+  } psb_c_sprec; 
-
+  
-psb_c_sprec* psb_c_new_sprec();
+  psb_c_sprec* psb_c_new_sprec();
-
+  
-psb_i_t  psb_c_sprecinit(psb_c_sprec *ph, const char *ptype);
+  psb_i_t  psb_c_sprecinit(psb_i_t ictxt, psb_c_sprec *ph, const char *ptype);
-psb_i_t  psb_c_sprecbld(psb_c_sspmat *ah, psb_c_descriptor *cdh, psb_c_sprec *ph);
+  psb_i_t  psb_c_sprecbld(psb_c_sspmat *ah, psb_c_descriptor *cdh, psb_c_sprec *ph);
-psb_i_t  psb_c_sprecfree(psb_c_sprec *ph);
+  psb_i_t  psb_c_sprecfree(psb_c_sprec *ph);
 #ifdef __cplusplus
 }
 #endif

cbind/prec/psb_c_zprec.h

@@ -8,15 +8,15 @@
 extern "C" {
 #endif
 
-typedef struct PSB_C_ZPREC {
+  typedef struct PSB_C_ZPREC {
-  void *zprec;
+    void *zprec;
-} psb_c_zprec; 
+  } psb_c_zprec; 
-
+  
-psb_c_zprec* psb_c_new_zprec();
+  psb_c_zprec* psb_c_new_zprec();
-
+  
-psb_i_t  psb_c_zprecinit(psb_c_zprec *ph, const char *ptype);
+  psb_i_t  psb_c_zprecinit(psb_i_t ictxt, psb_c_zprec *ph, const char *ptype);
-psb_i_t  psb_c_zprecbld(psb_c_zspmat *ah, psb_c_descriptor *cdh, psb_c_zprec *ph);
+  psb_i_t  psb_c_zprecbld(psb_c_zspmat *ah, psb_c_descriptor *cdh, psb_c_zprec *ph);
-psb_i_t  psb_c_zprecfree(psb_c_zprec *ph);
+  psb_i_t  psb_c_zprecfree(psb_c_zprec *ph);
 #ifdef __cplusplus
 }
 #endif

cbind/prec/psb_cprec_cbind_mod.f90

@@ -13,12 +13,13 @@ module psb_cprec_cbind_mod
 contains 
 
 
-  function  psb_c_cprecinit(ph,ptype) bind(c) result(res)
+  function  psb_c_cprecinit(ictxt,ph,ptype) bind(c) result(res)
     use psb_base_mod
     use psb_prec_mod
     use psb_base_string_cbind_mod
     implicit none 
-    integer(c_int)          :: res
+    integer(psb_c_int), value :: ictxt   
+    integer(psb_c_int)          :: res
     type(psb_c_cprec) :: ph
     character(c_char)       :: ptype(*)
     type(psb_cprec_type), pointer :: precp
@@ -36,7 +37,7 @@ contains
 
     call stringc2f(ptype,fptype)
     
-    call psb_precinit(precp,fptype,info) 
+    call psb_precinit(ictxt,precp,fptype,info) 
     
     res = min(0,info)
     return
@@ -51,7 +52,7 @@ contains
     use psb_base_string_cbind_mod
     implicit none 
     
-    integer(c_int) :: res
+    integer(psb_c_int) :: res
     type(psb_c_cspmat) :: ah
     type(psb_c_cprec) :: ph
     type(psb_c_descriptor) :: cdh
@@ -94,7 +95,7 @@ contains
     use psb_base_string_cbind_mod
     implicit none 
     
-    integer(c_int) :: res
+    integer(psb_c_int) :: res
     type(psb_c_cprec) :: ph
 
     type(psb_cprec_type), pointer :: precp

cbind/prec/psb_dprec_cbind_mod.f90

@@ -13,12 +13,13 @@ module psb_dprec_cbind_mod
 contains 
 
 
-  function  psb_c_dprecinit(ph,ptype) bind(c) result(res)
+  function  psb_c_dprecinit(ictxt,ph,ptype) bind(c) result(res)
     use psb_base_mod
     use psb_prec_mod
     use psb_base_string_cbind_mod
     implicit none 
-    integer(c_int)          :: res
+    integer(psb_c_int), value :: ictxt   
+    integer(psb_c_int)          :: res
     type(psb_c_dprec) :: ph
     character(c_char)       :: ptype(*)
     type(psb_dprec_type), pointer :: precp
@@ -36,7 +37,7 @@ contains
 
     call stringc2f(ptype,fptype)
     
-    call psb_precinit(precp,fptype,info) 
+    call psb_precinit(ictxt,precp,fptype,info) 
     
     res = min(0,info)
     return
@@ -51,7 +52,7 @@ contains
     use psb_base_string_cbind_mod
     implicit none 
     
-    integer(c_int) :: res
+    integer(psb_c_int) :: res
     type(psb_c_dspmat) :: ah
     type(psb_c_dprec) :: ph
     type(psb_c_descriptor) :: cdh
@@ -94,7 +95,7 @@ contains
     use psb_base_string_cbind_mod
     implicit none 
     
-    integer(c_int) :: res
+    integer(psb_c_int) :: res
     type(psb_c_dprec) :: ph
 
     type(psb_dprec_type), pointer :: precp

cbind/prec/psb_sprec_cbind_mod.f90

@@ -13,12 +13,13 @@ module psb_sprec_cbind_mod
 contains 
 
 
-  function  psb_c_sprecinit(ph,ptype) bind(c) result(res)
+  function  psb_c_sprecinit(ictxt,ph,ptype) bind(c) result(res)
     use psb_base_mod
     use psb_prec_mod
     use psb_base_string_cbind_mod
     implicit none 
-    integer(c_int)          :: res
+    integer(psb_c_int), value :: ictxt   
+    integer(psb_c_int)          :: res
     type(psb_c_sprec) :: ph
     character(c_char)       :: ptype(*)
     type(psb_sprec_type), pointer :: precp
@@ -36,7 +37,7 @@ contains
 
     call stringc2f(ptype,fptype)
     
-    call psb_precinit(precp,fptype,info) 
+    call psb_precinit(ictxt,precp,fptype,info) 
     
     res = min(0,info)
     return
@@ -51,7 +52,7 @@ contains
     use psb_base_string_cbind_mod
     implicit none 
     
-    integer(c_int) :: res
+    integer(psb_c_int) :: res
     type(psb_c_sspmat) :: ah
     type(psb_c_sprec) :: ph
     type(psb_c_descriptor) :: cdh
@@ -94,7 +95,7 @@ contains
     use psb_base_string_cbind_mod
     implicit none 
     
-    integer(c_int) :: res
+    integer(psb_c_int) :: res
     type(psb_c_sprec) :: ph
 
     type(psb_sprec_type), pointer :: precp

cbind/prec/psb_zprec_cbind_mod.f90

@@ -13,12 +13,13 @@ module psb_zprec_cbind_mod
 contains 
 
 
-  function  psb_c_zprecinit(ph,ptype) bind(c) result(res)
+  function  psb_c_zprecinit(ictxt,ph,ptype) bind(c) result(res)
     use psb_base_mod
     use psb_prec_mod
     use psb_base_string_cbind_mod
     implicit none 
-    integer(c_int)          :: res
+    integer(psb_c_int), value :: ictxt   
+    integer(psb_c_int)          :: res
     type(psb_c_zprec) :: ph
     character(c_char)       :: ptype(*)
     type(psb_zprec_type), pointer :: precp
@@ -36,7 +37,7 @@ contains
 
     call stringc2f(ptype,fptype)
     
-    call psb_precinit(precp,fptype,info) 
+    call psb_precinit(ictxt,precp,fptype,info) 
     
     res = min(0,info)
     return
@@ -51,7 +52,7 @@ contains
     use psb_base_string_cbind_mod
     implicit none 
     
-    integer(c_int) :: res
+    integer(psb_c_int) :: res
     type(psb_c_zspmat) :: ah
     type(psb_c_zprec) :: ph
     type(psb_c_descriptor) :: cdh
@@ -94,7 +95,7 @@ contains
     use psb_base_string_cbind_mod
     implicit none 
     
-    integer(c_int) :: res
+    integer(psb_c_int) :: res
     type(psb_c_zprec) :: ph
 
     type(psb_zprec_type), pointer :: precp

cbind/test/pargen/ppdec.c

@@ -320,7 +320,7 @@ int main(int argc, char *argv[])
   psb_c_barrier(ictxt);
   /* Set up the preconditioner */ 
   ph  = psb_c_new_dprec();
-  psb_c_dprecinit(ph,ptype);
+  psb_c_dprecinit(ictxt,ph,ptype);
   ret=psb_c_dprecbld(ah,cdh,ph);
   fprintf(stderr,"From psb_c_dprecbld: %d\n",ret); 
 

configure.ac

@@ -144,7 +144,7 @@ else
 AC_LANG([C])
 if test "X$MPICC" = "X" ; then
     # This is our MPICC compiler preference: it will override ACX_MPI's first try.
-    AC_CHECK_PROGS([MPICC],[mpxlc mpcc mpiicc mpicc cc])
+    AC_CHECK_PROGS([MPICC],[mpxlc mpiicc mpcc mpicc cc])
 fi
 ACX_MPI([], [AC_MSG_ERROR([[Cannot find any suitable MPI implementation for C]])])
 AC_PROG_CC_STDC
@@ -380,8 +380,8 @@ if test "X$FCOPT" == "X" ; then
    FCOPT="-g -O3  $FCOPT"
  elif test "X$psblas_cv_fc" == X"xlf" ; then
    # XL compiler : consider using -qarch=auto
-     FCOPT="-O3 -qarch=auto  -qfixed -qsuffix=f=f:cpp=F -qlanglvl=extended $FCOPT"
+     FCOPT="-O3 -qarch=auto  -qlanglvl=extended -qxlf2003=polymorphic:autorealloc $FCOPT"
-     FCFLAGS="-qhalt=e $FCFLAGS"
+     FCFLAGS="-qhalt=e -qlanglvl=extended -qxlf2003=polymorphic:autorealloc $FCFLAGS"
  elif test "X$psblas_cv_fc" == X"ifc" ; then
    # other compilers ..
      FCOPT="-O3 $FCOPT"

docs/html/footnode.html

@@ -1,6 +1,6 @@
 <!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.0 Transitional//EN">
 
-<!--Converted with LaTeX2HTML 2017.2 (Released Jan 23, 2017) -->
+<!--Converted with LaTeX2HTML 2018 (Released Feb 1, 2018) -->
 <HTML>
 <HEAD>
 <TITLE>Footnotes</TITLE>
@@ -9,7 +9,7 @@
 <META NAME="resource-type" CONTENT="document">
 <META NAME="distribution" CONTENT="global">
 
-<META NAME="Generator" CONTENT="LaTeX2HTML v2017.2">
+<META NAME="Generator" CONTENT="LaTeX2HTML v2018">
 <META HTTP-EQUIV="Content-Style-Type" CONTENT="text/css">
 
 <LINK REL="STYLESHEET" HREF="userhtml.css">
@@ -137,7 +137,7 @@ sample scatter/gather routines.
 .
 </PRE>
 </DD>
-<DT><A NAME="foot7809">... follows</A><A
+<DT><A NAME="foot7826">... follows</A><A
  HREF="node126.html#tex2html31"><SUP><SPAN CLASS="arabic">4</SPAN></SUP></A></DT>
 <DD>The string is case-insensitive
 
@@ -173,7 +173,7 @@ sample scatter/gather routines.
 .
 </PRE>
 </DD>
-<DT><A NAME="foot8258">... method</A><A
+<DT><A NAME="foot8281">... method</A><A
  HREF="node133.html#tex2html32"><SUP><SPAN CLASS="arabic">5</SPAN></SUP></A></DT>
 <DD>Note:
       the implementation is for <SPAN CLASS="MATH"><IMG