diff --git a/opt/Makefile b/opt/Makefile
index 2663bf58..a416e459 100644
--- a/opt/Makefile
+++ b/opt/Makefile
@@ -16,9 +16,10 @@ FINCLUDES=$(FMFLAG)$(LIBDIR) $(FMFLAG).
 EXEDIR=./runs
 
 all: libopt.a
-libopt.a: psb_d_ell_impl.o psb_d_ell_mat_mod.o
-	ar cur libopt.a psb_d_ell_impl.o psb_d_ell_mat_mod.o
+libopt.a: psb_d_ell_impl.o psb_d_ell_mat_mod.o rsb_mod.o psb_d_rsb_mat_mod.o
+	ar cur libopt.a psb_d_ell_impl.o psb_d_ell_mat_mod.o rsb_mod.o psb_d_rsb_mat_mod.o
 psb_d_ell_impl.o: psb_d_ell_mat_mod.o
+psb_d_rsb_mat_mod.o:  rsb_mod.o
 
 clean: 
 	/bin/rm -f  psb_d_ell_impl.o psb_d_ell_mat_mod.o *$(.mod)
diff --git a/opt/psb_d_rsb_mat_mod.F03 b/opt/psb_d_rsb_mat_mod.F03
new file mode 100644
index 00000000..f30c4ce9
--- /dev/null
+++ b/opt/psb_d_rsb_mat_mod.F03
@@ -0,0 +1,768 @@
+!
+!
+! FIXME:
+! * some RSB constants are used in their value form, and with no explanation
+! * error handling
+! * PSBLAS interface adherence
+! * should test and fix all the problems that for sure will occur
+! * ..
+! 
+module psb_d_rsb_mat_mod
+  use psb_d_base_mat_mod
+  use rsb_mod
+#ifdef HAVE_LIBRSB
+  use iso_c_binding  
+#endif
+#if 1
+#define PSBRSB_DEBUG(MSG) write(*,*) __FILE__,':',__LINE__,':',MSG
+#else
+#define PSBRSB_DEBUG(MSG)
+#endif
+  integer :: c_f_order=2        ! FIXME: here should use RSB_FLAG_WANT_COLUMN_MAJOR_ORDER
+  integer,parameter :: c_f_index=1 ! 0x001000 ! FIXME: here should use RSB_FLAG_FORTRAN_INDICES_INTERFACE
+  integer,parameter  :: c_d_typecode=68 ! FIXME: here should use ..
+  integer,parameter  :: c_def_flags =c_f_index ! FIXME: here should use ..
+  integer :: c_srt_flags =4 ! flags if rsb input is row major sorted ..
+  integer :: c_own_flags =2 ! flags if rsb input shall not be freed by rsb
+  integer :: c_upd_flags =c_f_index ! flags for when updating the assembled rsb matrix
+  type, extends(psb_d_base_sparse_mat) :: psb_d_rsb_sparse_mat
+#ifdef HAVE_LIBRSB
+    type(c_ptr) :: rsbmptr=c_null_ptr
+  contains 
+    procedure, pass(a) :: get_size     => d_rsb_get_size
+    procedure, pass(a) :: get_nzeros   => d_rsb_get_nzeros
+    procedure, pass(a) :: get_ncols   => d_rsb_get_ncols
+    procedure, pass(a) :: get_nrows   => d_rsb_get_nrows
+    procedure, pass(a) :: get_fmt      => d_rsb_get_fmt
+    procedure, pass(a) :: sizeof       => d_rsb_sizeof
+    procedure, pass(a) :: d_csmm       => psb_d_rsb_csmm
+    procedure, pass(a) :: d_csmv       => psb_d_rsb_csmv
+    procedure, pass(a) :: d_inner_cssm => psb_d_rsb_cssm
+    procedure, pass(a) :: d_inner_cssv => psb_d_rsb_cssv
+    procedure, pass(a) :: d_scals      => psb_d_rsb_scals
+    procedure, pass(a) :: d_scal       => psb_d_rsb_scal
+    procedure, pass(a) :: csnmi        => psb_d_rsb_csnmi
+    procedure, pass(a) :: csnm1        => psb_d_rsb_csnm1
+    procedure, pass(a) :: rowsum       => psb_d_rsb_rowsum
+    procedure, pass(a) :: arwsum       => psb_d_rsb_arwsum
+    procedure, pass(a) :: colsum       => psb_d_rsb_colsum
+    procedure, pass(a) :: aclsum       => psb_d_rsb_aclsum
+!    procedure, pass(a) :: reallocate_nz => psb_d_rsb_reallocate_nz ! FIXME
+!    procedure, pass(a) :: allocate_mnnz => psb_d_rsb_allocate_mnnz ! FIXME
+    procedure, pass(a) :: cp_to_coo    => psb_d_cp_rsb_to_coo
+    procedure, pass(a) :: cp_from_coo  => psb_d_cp_rsb_from_coo
+    procedure, pass(a) :: cp_to_fmt    => psb_d_cp_rsb_to_fmt
+    procedure, pass(a) :: cp_from_fmt  => psb_d_cp_rsb_from_fmt
+    procedure, pass(a) :: mv_to_coo    => psb_d_mv_rsb_to_coo
+    procedure, pass(a) :: mv_from_coo  => psb_d_mv_rsb_from_coo
+    procedure, pass(a) :: mv_to_fmt    => psb_d_mv_rsb_to_fmt
+    procedure, pass(a) :: mv_from_fmt  => psb_d_mv_rsb_from_fmt
+    procedure, pass(a) :: csput        => psb_d_rsb_csput
+    procedure, pass(a) :: get_diag     => psb_d_rsb_get_diag
+    procedure, pass(a) :: csgetptn     => psb_d_rsb_csgetptn
+    procedure, pass(a) :: d_csgetrow   => psb_d_rsb_csgetrow
+    procedure, pass(a) :: get_nz_row   => d_rsb_get_nz_row
+    procedure, pass(a) :: reinit       => psb_d_rsb_reinit
+    procedure, pass(a) :: trim         => psb_d_rsb_trim ! evil
+    procedure, pass(a) :: print        => psb_d_rsb_print
+    procedure, pass(a) :: free         => d_rsb_free
+    procedure, pass(a) :: mold         => psb_d_rsb_mold
+    procedure, pass(a) :: psb_d_rsb_cp_from
+    generic, public    :: cp_from => psb_d_rsb_cp_from
+    procedure, pass(a) :: psb_d_rsb_mv_from
+    generic, public    :: mv_from => psb_d_rsb_mv_from
+
+#endif
+  end type psb_d_rsb_sparse_mat
+  ! FIXME: complete the following
+  !private :: d_rsb_get_nzeros, d_rsb_get_fmt
+  private :: d_rsb_to_psb_info
+#ifdef HAVE_LIBRSB
+  contains 
+
+  function psv_rsb_mat_init() result(res)
+    implicit none 
+    integer :: res
+    !PSBRSB_DEBUG('')
+    res=-1 ! FIXME
+#ifdef HAVE_LIBRSB
+    res=d_rsb_to_psb_info(rsb_init(c_null_ptr))
+#endif
+  end function psv_rsb_mat_init
+
+ function d_rsb_to_psb_info(info) result(res)
+    implicit none 
+    integer :: res,info
+    !PSBRSB_DEBUG('')
+    res=info
+  end function d_rsb_to_psb_info
+
+  function d_rsb_get_nzeros(a) result(res)
+    implicit none 
+    class(psb_d_rsb_sparse_mat), intent(in) :: a
+    integer :: res
+    !PSBRSB_DEBUG('')
+    res=rsb_get_matrix_nnz(a%rsbmptr)
+  end function d_rsb_get_nzeros
+
+  function d_rsb_get_nrows(a) result(res)
+    implicit none 
+    class(psb_d_rsb_sparse_mat), intent(in) :: a
+    integer :: res
+    !PSBRSB_DEBUG('')
+    res=rsb_get_matrix_n_rows(a%rsbmptr)
+  end function d_rsb_get_nrows
+
+  function d_rsb_get_ncols(a) result(res)
+    implicit none 
+    class(psb_d_rsb_sparse_mat), intent(in) :: a
+    integer :: res
+    !PSBRSB_DEBUG('')
+    res=rsb_get_matrix_n_columns(a%rsbmptr)
+  end function d_rsb_get_ncols
+
+  function d_rsb_get_fmt(a) result(res)
+    implicit none 
+    class(psb_d_rsb_sparse_mat), intent(in) :: a
+    character(len=5) :: res
+    !the following printout is harmful, here, if happening during a write :) (causes a deadlock)
+    !PSBRSB_DEBUG('')
+    res = 'RSB'
+  end function d_rsb_get_fmt
+  
+  function d_rsb_get_size(a) result(res)
+    implicit none 
+    class(psb_d_rsb_sparse_mat), intent(in) :: a
+    integer :: res
+    !PSBRSB_DEBUG('')
+    res = d_rsb_get_nzeros(a)
+  end function d_rsb_get_size
+
+  function d_rsb_sizeof(a) result(res)
+    implicit none 
+    class(psb_d_rsb_sparse_mat), intent(in) :: a
+    integer(psb_long_int_k_) :: res
+    !PSBRSB_DEBUG('')
+    res=rsb_sizeof(a%rsbmptr)
+  end function d_rsb_sizeof
+
+subroutine psb_d_rsb_csmv(alpha,a,x,beta,y,info,trans) 
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(in) :: a
+  real(psb_dpk_), intent(in) :: alpha, beta, x(:)
+  real(psb_dpk_), intent(inout)       :: y(:)
+  integer, intent(out)                :: info
+  character, optional, intent(in)     :: trans
+  character :: trans_
+    PSBRSB_DEBUG('')
+  info = psb_success_
+
+  if (present(trans)) then
+    trans_ = trans
+  else
+    trans_ = 'N'
+  end if
+  info=d_rsb_to_psb_info(rsb_spmv(a%rsbmptr,x,y,alpha,beta,1,1,rsb_psblas_trans_to_rsb_trans(trans_)))
+end subroutine psb_d_rsb_csmv
+
+subroutine psb_d_rsb_cssv(alpha,a,x,beta,y,info,trans) 
+  ! FIXME: and what when x is an alias of y ?
+  ! FIXME: ignoring beta
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(in) :: a
+  real(psb_dpk_), intent(in) :: alpha, beta, x(:)
+  real(psb_dpk_), intent(inout)       :: y(:)
+  integer, intent(out)                :: info
+  character, optional, intent(in)     :: trans
+  character :: trans_
+    PSBRSB_DEBUG('')
+  info = psb_success_
+
+  if (present(trans)) then
+    trans_ = trans
+  else
+    trans_ = 'N'
+  end if
+  info=d_rsb_to_psb_info(rsb_spsv(a%rsbmptr,x,y,alpha,1,1,rsb_psblas_trans_to_rsb_trans(trans_)))
+end subroutine psb_d_rsb_cssv
+
+subroutine psb_d_rsb_scals(d,a,info) 
+  use psb_sparse_mod
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(inout) :: a
+  real(psb_dpk_), intent(in)      :: d
+  integer, intent(out)            :: info
+    PSBRSB_DEBUG('')
+  info=d_rsb_to_psb_info(rsb_elemental_scale(a%rsbmptr,d))
+end subroutine psb_d_rsb_scals
+
+subroutine psb_d_rsb_scal(d,a,info) 
+  use psb_sparse_mod
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(inout) :: a
+  real(psb_dpk_), intent(in)      :: d(:)
+  integer, intent(out)            :: info
+    PSBRSB_DEBUG('')
+  info=d_rsb_to_psb_info(rsb_scale_rows(a%rsbmptr,d))
+end subroutine psb_d_rsb_scal
+
+  subroutine  d_rsb_free(a) 
+    implicit none 
+    class(psb_d_rsb_sparse_mat), intent(inout) :: a
+    type(c_ptr) :: dummy
+    !PSBRSB_DEBUG('freeing RSB matrix')
+    dummy=rsb_free_sparse_matrix(a%rsbmptr)
+  end subroutine d_rsb_free
+
+subroutine  psb_d_rsb_trim(a)
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(inout) :: a
+    !PSBRSB_DEBUG('')
+  ! FIXME: this is supposed to remain empty for RSB
+end subroutine  psb_d_rsb_trim
+
+    subroutine psb_d_rsb_print(iout,a,iv,eirs,eics,head,ivr,ivc)
+      integer, intent(in)               :: iout
+      class(psb_d_rsb_sparse_mat), intent(in) :: a   
+      integer, intent(in), optional     :: iv(:)
+      integer, intent(in), optional     :: eirs,eics
+      character(len=*), optional        :: head
+      integer, intent(in), optional     :: ivr(:), ivc(:)
+      integer             :: info
+    PSBRSB_DEBUG('')
+      ! FIXME: UNFINISHED
+      info=rsb_print_matrix_t(a%rsbmptr)
+    end subroutine psb_d_rsb_print
+
+    subroutine psb_d_rsb_get_diag(a,d,info) 
+      class(psb_d_rsb_sparse_mat), intent(in) :: a
+      real(psb_dpk_), intent(out)     :: d(:)
+      integer, intent(out)            :: info
+      !PSBRSB_DEBUG('')
+      info=rsb_getdiag(a%rsbmptr,d)
+    end subroutine psb_d_rsb_get_diag
+  
+function psb_d_rsb_csnmi(a) result(res)
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(in) :: a
+  real(psb_dpk_),target         :: res
+  real(psb_dpk_)         :: resa(1)
+  integer             :: info
+    !PSBRSB_DEBUG('')
+  info=rsb_infinity_norm(a%rsbmptr,resa,rsb_psblas_trans_to_rsb_trans('N'))
+  !info=rsb_infinity_norm(a%rsbmptr,c_loc(res),rsb_psblas_trans_to_rsb_trans('N'))
+  res=resa(1)
+end function psb_d_rsb_csnmi
+
+function psb_d_rsb_csnm1(a) result(res)
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(in) :: a
+  real(psb_dpk_)         :: res
+  real(psb_dpk_)         :: resa(1)
+  integer             :: info
+    PSBRSB_DEBUG('')
+  info=rsb_one_norm(a%rsbmptr,resa,rsb_psblas_trans_to_rsb_trans('N'))
+  !info=rsb_one_norm(a%rsbmptr,res,rsb_psblas_trans_to_rsb_trans('N'))
+end function psb_d_rsb_csnm1
+
+subroutine psb_d_rsb_aclsum(d,a) 
+  use psb_sparse_mod
+  class(psb_d_rsb_sparse_mat), intent(in) :: a
+  real(psb_dpk_), intent(out)              :: d(:)
+    PSBRSB_DEBUG('')
+  info=rsb_absolute_columns_sums(a%rsbmptr,d)
+end subroutine psb_d_rsb_aclsum
+
+subroutine psb_d_rsb_arwsum(d,a) 
+  use psb_sparse_mod
+  class(psb_d_rsb_sparse_mat), intent(in) :: a
+  real(psb_dpk_), intent(out)              :: d(:)
+    PSBRSB_DEBUG('')
+  info=rsb_absolute_rows_sums(a%rsbmptr,d)
+end subroutine psb_d_rsb_arwsum
+
+subroutine psb_d_rsb_csmm(alpha,a,x,beta,y,info,trans) 
+  use psb_sparse_mod
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(in) :: a
+  real(psb_dpk_), intent(in)          :: alpha, beta, x(:,:)
+  real(psb_dpk_), intent(inout)       :: y(:,:)
+  integer, intent(out)                :: info
+  character, optional, intent(in)     :: trans
+
+  character :: trans_
+  integer :: ldy,ldx,nc
+    PSBRSB_DEBUG('')
+
+  if (present(trans)) then
+    trans_ = trans
+  else
+    trans_ = 'N'
+  end if
+  ldx=size(x,1); ldy=size(y,1)
+  nc=min(size(x,2),size(y,2) )
+!!$  info=d_rsb_to_psb_info(rsb_spmm(a%rsbmptr,x,y,ldx,ldy,nc,rsb_psblas_trans_to_rsb_trans(trans_),alpha,beta,c_f_order))
+end subroutine psb_d_rsb_csmm
+
+subroutine psb_d_rsb_cssm(alpha,a,x,beta,y,info,trans) 
+  use psb_sparse_mod
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(in) :: a
+  real(psb_dpk_), intent(in)          :: alpha, beta, x(:,:)
+  real(psb_dpk_), intent(inout)       :: y(:,:)
+  integer, intent(out)                :: info
+  character, optional, intent(in)     :: trans
+  integer :: ldy,ldx,nc
+  character :: trans_
+    PSBRSB_DEBUG('')
+  if (present(trans)) then
+    trans_ = trans
+  else
+    trans_ = 'N'
+  end if
+  ldx=size(x,1); ldy=size(y,1)
+  nc=min(size(x,2),size(y,2) )
+!!$  info=d_rsb_to_psb_info(rsb_spsm(a%rsbmptr,y,ldy,nc,rsb_psblas_trans_to_rsb_trans(trans_),alpha,beta,c_f_order))
+end subroutine
+
+subroutine psb_d_rsb_rowsum(d,a) 
+  use psb_sparse_mod
+  class(psb_d_rsb_sparse_mat), intent(in) :: a
+  real(psb_dpk_), intent(out)             :: d(:)
+  integer :: info
+    PSBRSB_DEBUG('')
+  info=d_rsb_to_psb_info(rsb_rows_sums(a%rsbmptr,d))
+end subroutine psb_d_rsb_rowsum
+
+subroutine psb_d_rsb_colsum(d,a) 
+  use psb_sparse_mod
+  class(psb_d_rsb_sparse_mat), intent(in) :: a
+  real(psb_dpk_), intent(out)             :: d(:)
+  integer :: info
+    PSBRSB_DEBUG('')
+  info=d_rsb_to_psb_info(rsb_columns_sums(a%rsbmptr,d))
+end subroutine psb_d_rsb_colsum
+
+subroutine psb_d_rsb_mold(a,b,info) 
+  use psb_sparse_mod
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(in)  :: a
+  class(psb_d_base_sparse_mat), intent(out), allocatable  :: b
+  integer, intent(out)                    :: info
+  Integer :: err_act
+  character(len=20)  :: name='reallocate_nz'
+  logical, parameter :: debug=.false.
+    PSBRSB_DEBUG('')
+
+  call psb_get_erraction(err_act)
+  
+  allocate(psb_d_rsb_sparse_mat :: b, stat=info)
+
+  if (info /= psb_success_) then 
+    info = psb_err_alloc_dealloc_ 
+    call psb_errpush(info, name)
+    goto 9999
+  end if
+  return
+9999 continue
+  if (err_act /= psb_act_ret_) then
+    call psb_error()
+  end if
+  return
+end subroutine psb_d_rsb_mold
+
+subroutine psb_d_rsb_reinit(a,clear)
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(inout) :: a   
+  logical, intent(in), optional :: clear
+  Integer :: info
+    PSBRSB_DEBUG('')
+  info=d_rsb_to_psb_info(rsb_reinit_matrix(a%rsbmptr))
+end subroutine psb_d_rsb_reinit
+
+
+  function  d_rsb_get_nz_row(idx,a) result(res)
+    implicit none
+    class(psb_d_rsb_sparse_mat), intent(in) :: a
+    integer, intent(in)                  :: idx
+    integer                              :: res
+    integer                              :: info
+    PSBRSB_DEBUG('')
+    res=0
+    res=rsb_get_rows_nnz(a%rsbmptr,idx,idx,c_f_index,info)
+    info=d_rsb_to_psb_info(info)
+    if(info.ne.0)res=0
+  end function d_rsb_get_nz_row
+
+subroutine psb_d_cp_rsb_to_coo(a,b,info) 
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(in)  :: a
+  class(psb_d_coo_sparse_mat), intent(inout) :: b
+  integer, intent(out)                      :: info
+
+  integer, allocatable :: itemp(:)
+  !locals
+  logical             :: rwshr_
+  Integer             :: nza, nr, nc,i,j,irw, idl,err_act
+  integer             :: debug_level, debug_unit
+  character(len=20)   :: name
+    PSBRSB_DEBUG('')
+  info = psb_success_
+  nr  = a%get_nrows()
+  nc  = a%get_ncols()
+  nza = a%get_nzeros()
+  call b%allocate(nr,nc,nza)
+  call b%psb_d_base_sparse_mat%cp_from(a%psb_d_base_sparse_mat)
+  info=d_rsb_to_psb_info(rsb_get_coo(a%rsbmptr,b%val,b%ia,b%ja,c_f_index))
+  call b%set_nzeros(a%get_nzeros())
+  call b%set_nrows(a%get_nrows())
+  call b%set_ncols(a%get_ncols())
+  call b%fix(info)
+  !write(*,*)b%val
+  !write(*,*)b%ia
+  !write(*,*)b%ja
+  !write(*,*)b%get_nrows()
+  !write(*,*)b%get_ncols()
+  !write(*,*)b%get_nzeros()
+  !write(*,*)a%get_nrows()
+  !write(*,*)a%get_ncols()
+  !write(*,*)a%get_nzeros()
+end subroutine psb_d_cp_rsb_to_coo
+
+subroutine psb_d_cp_rsb_to_fmt(a,b,info) 
+  use psb_sparse_mod
+  implicit none 
+
+  class(psb_d_rsb_sparse_mat), intent(in)   :: a
+  class(psb_d_base_sparse_mat), intent(inout) :: b
+  integer, intent(out)                       :: info
+
+  !locals
+  type(psb_d_coo_sparse_mat) :: tmp
+  logical             :: rwshr_
+  Integer             :: nza, nr, i,j,irw, idl,err_act, nc
+  integer              :: debug_level, debug_unit
+  character(len=20)   :: name
+    PSBRSB_DEBUG('')
+
+  info = psb_success_
+
+  select type (b)
+  type is (psb_d_coo_sparse_mat) 
+    call a%cp_to_coo(b,info)
+
+  type is (psb_d_rsb_sparse_mat) 
+    call b%psb_d_base_sparse_mat%cp_from(a%psb_d_base_sparse_mat)! FIXME: ?
+    b%rsbmptr=rsb_clone(a%rsbmptr) ! FIXME is thi enough ?
+    ! FIXME: error handling needed here
+
+  class default
+    call a%cp_to_coo(tmp,info)
+    if (info == psb_success_) call b%mv_from_coo(tmp,info)
+  end select
+end subroutine psb_d_cp_rsb_to_fmt
+
+subroutine psb_d_cp_rsb_from_coo(a,b,info) 
+  use psb_sparse_mod
+  implicit none 
+
+  class(psb_d_rsb_sparse_mat), intent(inout) :: a
+  class(psb_d_coo_sparse_mat), intent(in)    :: b
+  integer, intent(out)                        :: info
+
+  integer, allocatable :: itemp(:)
+  !locals
+  logical             :: rwshr_
+  Integer             :: nza, nr, i,j,irw, idl,err_act, nc
+  integer             :: debug_level, debug_unit
+  integer             :: flags=c_def_flags
+  character(len=20)   :: name
+    PSBRSB_DEBUG('')
+
+   info = psb_success_
+   ! This is to have fix_coo called behind the scenes
+   if(b%is_sorted()) flags=flags+c_srt_flags 
+   !write (*,*) b%val
+   ! FIXME: and if sorted ? the process could be speeded up !
+   a%rsbmptr=rsb_allocate_rsb_sparse_matrix_const&
+  &(b%val,b%ia,b%ja,b%get_nzeros(),c_d_typecode,b%get_nrows(),b%get_ncols(),1,1,flags,info)
+    info=d_rsb_to_psb_info(info)
+   ! FIXME: should destroy tmp ?
+end subroutine psb_d_cp_rsb_from_coo
+
+subroutine psb_d_cp_rsb_from_fmt(a,b,info) 
+  use psb_sparse_mod
+  implicit none 
+
+  class(psb_d_rsb_sparse_mat), intent(inout) :: a
+  class(psb_d_base_sparse_mat), intent(in)   :: b
+  integer, intent(out)                        :: info
+
+  !locals
+  type(psb_d_coo_sparse_mat) :: tmp
+  logical             :: rwshr_
+  Integer             :: nz, nr, i,j,irw, idl,err_act, nc
+  integer              :: debug_level, debug_unit
+  character(len=20)   :: name
+    PSBRSB_DEBUG('')
+
+  info = psb_success_
+
+  select type (b)
+  type is (psb_d_coo_sparse_mat) 
+    call a%cp_from_coo(b,info)
+
+  type is (psb_d_csr_sparse_mat) 
+  a%rsbmptr=rsb_allocate_rsb_sparse_matrix_from_csr_const&
+  &(b%val,b%irp,b%ja,b%get_nzeros(),c_d_typecode,b%get_nrows(),b%get_ncols(),1,1,c_def_flags+c_srt_flags,info)
+    info=d_rsb_to_psb_info(info)
+
+  type is (psb_d_rsb_sparse_mat) 
+    call b%cp_to_fmt(a,info)  ! FIXME
+    ! FIXME: missing error handling 
+
+  class default
+    call b%cp_to_coo(tmp,info)
+    if (info == psb_success_) call a%mv_from_coo(tmp,info)
+  end select
+end subroutine psb_d_cp_rsb_from_fmt
+
+
+subroutine psb_d_rsb_csgetrow(imin,imax,a,nz,ia,ja,val,info,&
+     & jmin,jmax,iren,append,nzin,rscale,cscale)
+  use psb_sparse_mod
+  implicit none
+
+  class(psb_d_rsb_sparse_mat), intent(in) :: a
+  integer, intent(in)                  :: imin,imax
+  integer, intent(out)                 :: nz
+  integer, allocatable, intent(inout)  :: ia(:), ja(:)
+  real(psb_dpk_), allocatable,  intent(inout)    :: val(:)
+  integer,intent(out)                  :: info
+  logical, intent(in), optional        :: append
+  integer, intent(in), optional        :: iren(:)
+  integer, intent(in), optional        :: jmin,jmax, nzin
+  logical, intent(in), optional        :: rscale,cscale
+
+  logical :: append_, rscale_, cscale_ 
+  integer :: nzin_, jmin_, jmax_, err_act, i
+  character(len=20)  :: name='csget'
+  logical, parameter :: debug=.false.
+  ! FIXME: MISSING THE HANDLING OF OPTIONS, HERE
+    PSBRSB_DEBUG('')
+
+  call psb_erractionsave(err_act)
+  info = psb_success_
+  
+  if (present(jmin)) then
+    jmin_ = jmin
+  else
+    jmin_ = 1
+  endif
+  if (present(jmax)) then
+    jmax_ = jmax
+  else
+    jmax_ = a%get_ncols()
+  endif
+
+  if ((imax<imin).or.(jmax_<jmin_)) then 
+    nz = 0
+    return
+  end if
+
+  if (present(append)) then
+    append_=append
+  else
+    append_=.false.
+  endif
+  if ((append_).and.(present(nzin))) then 
+    nzin_ = nzin
+  else
+    nzin_ = 0
+  endif
+  if (present(rscale)) then 
+    rscale_ = rscale
+  else
+    rscale_ = .false.
+  endif
+  if (present(cscale)) then 
+    cscale_ = cscale
+  else
+    cscale_ = .false.
+  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
+
+    nzin_=rsb_get_block_nnz(a%rsbmptr,imin,imax,imin,imax,c_f_index,info)
+    ! FIXME: unfinished; missing error handling ..
+
+    call psb_ensure_size(nzin_,ia,info)
+    if (info == psb_success_) call psb_ensure_size(nzin_,ja,info)
+    if (info == psb_success_) call psb_ensure_size(nzin_,val,info)
+    if (info /= psb_success_) return
+    
+    info=d_rsb_to_psb_info(rsb_get_block_sparse&
+    &(a%rsbmptr,val,imin,imax,jmin,jmax,ia,ja,c_null_ptr,c_null_ptr,nzin_,c_f_index))
+    ! FIXME: unfinished; missing error handling ..
+
+9999 continue
+  call psb_erractionrestore(err_act)
+
+  if (err_act == psb_act_abort_) then
+    call psb_error()
+    return
+  end if
+
+  end subroutine psb_d_rsb_csgetrow
+
+subroutine psb_d_rsb_csgetptn(imin,imax,a,nz,ia,ja,info,&
+     & jmin,jmax,iren,append,nzin,rscale,cscale)
+  use psb_sparse_mod
+  implicit none
+
+  class(psb_d_rsb_sparse_mat), intent(in) :: a
+  integer, intent(in)                  :: imin,imax
+  integer, intent(out)                 :: nz
+  integer, allocatable, intent(inout)  :: ia(:), ja(:)
+  integer,intent(out)                  :: info
+  logical, intent(in), optional        :: append
+  integer, intent(in), optional        :: iren(:)
+  integer, intent(in), optional        :: jmin,jmax, nzin
+  logical, intent(in), optional        :: rscale,cscale
+
+  logical :: append_, rscale_, cscale_ 
+  integer :: nzin_, jmin_, jmax_, err_act, i
+  character(len=20)  :: name='csget'
+  logical, parameter :: debug=.false.
+    PSBRSB_DEBUG('')
+
+    if (append) then 
+      nzin_ = nzin
+    else
+      nzin_ = 0
+    endif
+
+    !nzt = ..
+    nz = 0 
+
+    call psb_ensure_size(nzin_,ia,info)
+    if (info == psb_success_) call psb_ensure_size(nzin_,ja,info)
+
+    if (info /= psb_success_) return
+    nzin_=rsb_get_block_nnz(a%rsbmptr,imin,imax,imin,imax,c_f_index,info)
+    ! FIXME: unfinished; missing error handling ..
+
+    call psb_ensure_size(nzin_,ia,info)
+    if (info == psb_success_) call psb_ensure_size(nzin_,ja,info)
+    if (info /= psb_success_) return
+    
+    info=d_rsb_to_psb_info(rsb_get_block_sparse_pattern&
+    &(a%rsbmptr,imin,imax,jmin,jmax,ia,ja,c_null_ptr,c_null_ptr,nzin_,c_f_index))
+    ! FIXME: unfinished; missing error handling ..
+end subroutine psb_d_rsb_csgetptn
+
+subroutine psb_d_rsb_csput(nz,ia,ja,val,a,imin,imax,jmin,jmax,info,gtl) 
+  use psb_sparse_mod
+  implicit none 
+
+  class(psb_d_rsb_sparse_mat), intent(inout) :: a
+  real(psb_dpk_), intent(in)      :: val(:)
+  integer, intent(in)             :: nz, ia(:), ja(:), imin,imax,jmin,jmax
+  integer, intent(out)            :: info
+  integer, intent(in), optional   :: gtl(:)
+
+
+  Integer            :: err_act
+  character(len=20)  :: name='d_rsb_csput'
+  logical, parameter :: debug=.false.
+  integer            :: nza, i,j,k, nzl, isza, int_err(5)
+    PSBRSB_DEBUG('')
+  info=d_rsb_to_psb_info(rsb_update_elements(a%rsbmptr,val,ia,ja,nz,c_upd_flags))
+end subroutine psb_d_rsb_csput
+
+subroutine psb_d_mv_rsb_to_coo(a,b,info) 
+  use psb_sparse_mod
+  implicit none 
+
+  class(psb_d_rsb_sparse_mat), intent(inout) :: a
+  class(psb_d_coo_sparse_mat), intent(inout)   :: b
+  integer, intent(out)                        :: info
+    PSBRSB_DEBUG('')
+  ! FIXME: use rsb_switch_rsb_matrix_to_coo_sorted !
+  call psb_d_cp_rsb_to_coo(a,b,info)
+  call a%free()
+end subroutine psb_d_mv_rsb_to_coo
+
+    subroutine psb_d_mv_rsb_to_fmt(a,b,info) 
+      class(psb_d_rsb_sparse_mat), intent(inout) :: a
+      class(psb_d_base_sparse_mat), intent(inout)  :: b
+      integer, intent(out)                        :: info
+    PSBRSB_DEBUG('')
+  ! FIXME: could use here rsb_switch_rsb_matrix_to_csr_sorted
+  call psb_d_cp_rsb_to_fmt(a,b,info)
+  call d_rsb_free(a)
+  a%rsbmptr=c_null_ptr
+    end subroutine psb_d_mv_rsb_to_fmt
+  
+subroutine psb_d_mv_rsb_from_fmt(a,b,info) 
+  use psb_sparse_mod
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(inout)  :: a
+  class(psb_d_base_sparse_mat), intent(inout) :: b
+  integer, intent(out)                         :: info
+  ! FIXME: could use here rsb_allocate_rsb_sparse_matrix_from_csr_inplace
+  !if(b%is_sorted()) flags=flags+c_srt_flags 
+  type(psb_d_coo_sparse_mat) :: tmp
+    PSBRSB_DEBUG('')
+  info = psb_success_
+  select type (b)
+  class default
+    call b%mv_to_coo(tmp,info)
+    if (info == psb_success_) call a%mv_from_coo(tmp,info)
+  end select
+end subroutine psb_d_mv_rsb_from_fmt
+
+subroutine psb_d_mv_rsb_from_coo(a,b,info) 
+  use psb_sparse_mod
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(inout) :: a
+  class(psb_d_coo_sparse_mat), intent(inout) :: b
+  integer, intent(out)                        :: info
+    PSBRSB_DEBUG('')
+  ! FIXME: should use rsb_allocate_rsb_sparse_matrix_inplace
+  !if(b%is_sorted()) flags=flags+c_srt_flags 
+  call a%cp_from_coo(b,info)
+  call b%free()
+end subroutine psb_d_mv_rsb_from_coo
+
+subroutine psb_d_rsb_cp_from(a,b)
+  use psb_sparse_mod
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(inout) :: a
+  type(psb_d_rsb_sparse_mat), intent(in)   :: b
+  Integer :: info
+  type(psb_d_coo_sparse_mat) :: tmp
+    PSBRSB_DEBUG('')
+  call b%cp_to_coo(tmp,info)
+  call a%mv_from_coo(tmp,info)
+  call tmp%free()
+end subroutine psb_d_rsb_cp_from
+
+subroutine psb_d_rsb_mv_from(a,b)
+  use psb_sparse_mod
+  implicit none 
+  class(psb_d_rsb_sparse_mat), intent(inout) :: a
+  type(psb_d_rsb_sparse_mat), intent(inout)   :: b
+  Integer :: info
+  type(psb_d_coo_sparse_mat) :: tmp
+    PSBRSB_DEBUG('')
+  call b%mv_to_coo(tmp,info)
+  call a%mv_from_coo(tmp,info)
+end subroutine psb_d_rsb_mv_from
+
+
+#endif
+end module psb_d_rsb_mat_mod
diff --git a/opt/rsb_mod.f03 b/opt/rsb_mod.f03
new file mode 100644
index 00000000..0b0bd155
--- /dev/null
+++ b/opt/rsb_mod.f03
@@ -0,0 +1,812 @@
+module rsb_mod
+  use iso_c_binding
+
+! module constants:
+
+interface
+integer(c_int) function &
+  &rsb_perror&
+  &(errval)&
+  &bind(c,name='rsb_perror')
+use iso_c_binding
+ integer(c_int), value  :: errval
+ end function rsb_perror
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_init&
+  &(io)&
+  &bind(c,name='rsb_init')
+use iso_c_binding
+ type(c_ptr), value  :: io
+ end function rsb_init
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_reinit&
+  &(io)&
+  &bind(c,name='rsb_reinit')
+use iso_c_binding
+ type(c_ptr), value  :: io
+ end function rsb_reinit
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_was_initialized&
+  &()&
+  &bind(c,name='rsb_was_initialized')
+use iso_c_binding
+ end function rsb_was_initialized
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_exit&
+  &()&
+  &bind(c,name='rsb_exit')
+use iso_c_binding
+ end function rsb_exit
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_meminfo&
+  &()&
+  &bind(c,name='rsb_meminfo')
+use iso_c_binding
+ end function rsb_meminfo
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_check_leak&
+  &()&
+  &bind(c,name='rsb_check_leak')
+use iso_c_binding
+ end function rsb_check_leak
+end interface
+
+interface
+type(c_ptr) function &
+  &rsb_allocate_rsb_sparse_matrix_from_csr_const&
+  &(VAc,IAc,JAc,nnz,typecode,m,k,Mb,Kb,flags,errvalp)&
+  &bind(c,name='rsb_allocate_rsb_sparse_matrix_from_csr_const')
+use iso_c_binding
+ real(c_double) :: VAc(*)
+ integer(c_int) :: IAc(*)
+ integer(c_int) :: JAc(*)
+ integer(c_int), value  :: nnz
+ integer(c_int), value  :: typecode
+ integer(c_int), value  :: m
+ integer(c_int), value  :: k
+ integer(c_int), value  :: Mb
+ integer(c_int), value  :: Kb
+ integer(c_int), value  :: flags
+ integer(c_int) :: errvalp
+ end function rsb_allocate_rsb_sparse_matrix_from_csr_const
+end interface
+
+interface
+type(c_ptr) function &
+  &rsb_allocate_rsb_sparse_matrix_from_csr_inplace&
+  &(VA,IA,JA,nnz,typecode,m,k,Mb,Kb,flags,errvalp)&
+  &bind(c,name='rsb_allocate_rsb_sparse_matrix_from_csr_inplace')
+use iso_c_binding
+ real(c_double) :: VA(*)
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ integer(c_int), value  :: nnz
+ integer(c_int), value  :: typecode
+ integer(c_int), value  :: m
+ integer(c_int), value  :: k
+ integer(c_int), value  :: Mb
+ integer(c_int), value  :: Kb
+ integer(c_int), value  :: flags
+ integer(c_int) :: errvalp
+ end function rsb_allocate_rsb_sparse_matrix_from_csr_inplace
+end interface
+
+interface
+type(c_ptr) function &
+  &rsb_allocate_rsb_sparse_matrix_const&
+  &(VAc,IAc,JAc,nnz,typecode,m,k,Mb,Kb,flags,errvalp)&
+  &bind(c,name='rsb_allocate_rsb_sparse_matrix_const')
+use iso_c_binding
+ real(c_double) :: VAc(*)
+ integer(c_int) :: IAc(*)
+ integer(c_int) :: JAc(*)
+ integer(c_int), value  :: nnz
+ integer(c_int), value  :: typecode
+ integer(c_int), value  :: m
+ integer(c_int), value  :: k
+ integer(c_int), value  :: Mb
+ integer(c_int), value  :: Kb
+ integer(c_int), value  :: flags
+ integer(c_int) :: errvalp
+ end function rsb_allocate_rsb_sparse_matrix_const
+end interface
+
+interface
+type(c_ptr) function &
+  &rsb_allocate_rsb_sparse_matrix_inplace&
+  &(VA,IA,JA,nnz,typecode,m,k,Mb,Kb,flags,errvalp)&
+  &bind(c,name='rsb_allocate_rsb_sparse_matrix_inplace')
+use iso_c_binding
+ real(c_double) :: VA(*)
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ integer(c_int), value  :: nnz
+ integer(c_int), value  :: typecode
+ integer(c_int), value  :: m
+ integer(c_int), value  :: k
+ integer(c_int), value  :: Mb
+ integer(c_int), value  :: Kb
+ integer(c_int), value  :: flags
+ integer(c_int) :: errvalp
+ end function rsb_allocate_rsb_sparse_matrix_inplace
+end interface
+
+interface
+type(c_ptr) function &
+  &rsb_free_sparse_matrix&
+  &(matrix)&
+  &bind(c,name='rsb_free_sparse_matrix')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ end function rsb_free_sparse_matrix
+end interface
+
+interface
+type(c_ptr) function &
+  &rsb_clone&
+  &(matrix)&
+  &bind(c,name='rsb_clone')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ end function rsb_clone
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_spmv&
+  &(matrix,x,y,alphap,betap,incx,incy,transa)&
+  &bind(c,name='rsb_spmv')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: x(*)
+ real(c_double) :: y(*)
+ real(c_double) :: alphap
+ real(c_double) :: betap
+ integer(c_int), value  :: incx
+ integer(c_int), value  :: incy
+ integer(c_int), value  :: transa
+ end function rsb_spmv
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_infinity_norm&
+  &(matrix,infinity_norm,transa)&
+  &bind(c,name='rsb_infinity_norm')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: infinity_norm(*)
+ integer(c_int), value  :: transa
+ end function rsb_infinity_norm
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_one_norm&
+  &(matrix,one_norm,transa)&
+  &bind(c,name='rsb_one_norm')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: one_norm(*)
+ integer(c_int), value  :: transa
+ end function rsb_one_norm
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_rows_sums&
+  &(matrix,d)&
+  &bind(c,name='rsb_rows_sums')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: d(*)
+ end function rsb_rows_sums
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_columns_sums&
+  &(matrix,d)&
+  &bind(c,name='rsb_columns_sums')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: d(*)
+ end function rsb_columns_sums
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_absolute_rows_sums&
+  &(matrix,d)&
+  &bind(c,name='rsb_absolute_rows_sums')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: d(*)
+ end function rsb_absolute_rows_sums
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_absolute_columns_sums&
+  &(matrix,d)&
+  &bind(c,name='rsb_absolute_columns_sums')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: d(*)
+ end function rsb_absolute_columns_sums
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_spsv&
+  &(matrix,x,y,alphap,incx,incy,transl)&
+  &bind(c,name='rsb_spsv')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: x(*)
+ real(c_double) :: y(*)
+ real(c_double) :: alphap
+ integer(c_int), value  :: incx
+ integer(c_int), value  :: incy
+ integer(c_int), value  :: transl
+ end function rsb_spsv
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_spsm&
+  &(matrix,b,ldb,nrhs,transt,alphap,betap,order)&
+  &bind(c,name='rsb_spsm')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: b(*)
+ integer(c_int), value  :: ldb
+ integer(c_int), value  :: nrhs
+ integer(c_int), value  :: transt
+ real(c_double) :: alphap
+ real(c_double) :: betap
+ integer(c_int), value  :: order
+ end function rsb_spsm
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_matrix_add_to_dense&
+  &(matrixa,alphap,transa,matrixb,ldb,nr,nc,rowmajor)&
+  &bind(c,name='rsb_matrix_add_to_dense')
+use iso_c_binding
+ type(c_ptr), value  :: matrixa
+ real(c_double) :: alphap
+ integer(c_int), value  :: transa
+ type(c_ptr), value  :: matrixb
+ integer(c_int), value  :: ldb
+ integer(c_int), value  :: nr
+ integer(c_int), value  :: nc
+ integer(c_int), value  :: rowmajor
+ end function rsb_matrix_add_to_dense
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_util_sort_row_major&
+  &(VA,IA,JA,nnz,m,k,typecode,flags)&
+  &bind(c,name='rsb_util_sort_row_major')
+use iso_c_binding
+ real(c_double) :: VA(*)
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ integer(c_int), value  :: nnz
+ integer(c_int), value  :: m
+ integer(c_int), value  :: k
+ integer(c_int), value  :: typecode
+ integer(c_int), value  :: flags
+ end function rsb_util_sort_row_major
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_util_sort_row_major_buffered&
+  &(VA,IA,JA,nnz,m,k,typecode,flags,WA,wb)&
+  &bind(c,name='rsb_util_sort_row_major_buffered')
+use iso_c_binding
+ real(c_double) :: VA(*)
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ integer(c_int), value  :: nnz
+ integer(c_int), value  :: m
+ integer(c_int), value  :: k
+ integer(c_int), value  :: typecode
+ integer(c_int), value  :: flags
+ type(c_ptr), value  :: WA
+ integer(c_int), value  :: wb
+ end function rsb_util_sort_row_major_buffered
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_util_sort_column_major&
+  &(VA,IA,JA,nnz,m,k,typecode,flags)&
+  &bind(c,name='rsb_util_sort_column_major')
+use iso_c_binding
+ real(c_double) :: VA(*)
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ integer(c_int), value  :: nnz
+ integer(c_int), value  :: m
+ integer(c_int), value  :: k
+ integer(c_int), value  :: typecode
+ integer(c_int), value  :: flags
+ end function rsb_util_sort_column_major
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_util_sortcoo&
+  &(VA,IA,JA,nnz,typecode,M_b,K_b,rpntr,cpntr,flags)&
+  &bind(c,name='rsb_util_sortcoo')
+use iso_c_binding
+ real(c_double) :: VA(*)
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ integer(c_int), value  :: nnz
+ integer(c_int), value  :: typecode
+ integer(c_int), value  :: M_b
+ integer(c_int), value  :: K_b
+ type(c_ptr), value  :: rpntr
+ type(c_ptr), value  :: cpntr
+ integer(c_int), value  :: flags
+ end function rsb_util_sortcoo
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_switch_rsb_matrix_to_coo_unsorted&
+  &(matrix,VA,IA,JA,flags)&
+  &bind(c,name='rsb_switch_rsb_matrix_to_coo_unsorted')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: VA(*)
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ integer(c_int), value  :: flags
+ end function rsb_switch_rsb_matrix_to_coo_unsorted
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_switch_rsb_matrix_to_coo_sorted&
+  &(matrix,VA,IA,JA,flags)&
+  &bind(c,name='rsb_switch_rsb_matrix_to_coo_sorted')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: VA(*)
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ integer(c_int), value  :: flags
+ end function rsb_switch_rsb_matrix_to_coo_sorted
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_switch_rsb_matrix_to_csr_sorted&
+  &(matrix,VA,IA,JA,flags)&
+  &bind(c,name='rsb_switch_rsb_matrix_to_csr_sorted')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: VA(*)
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ integer(c_int), value  :: flags
+ end function rsb_switch_rsb_matrix_to_csr_sorted
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_get_coo&
+  &(matrix,VA,IA,JA,flags)&
+  &bind(c,name='rsb_get_coo')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: VA(*)
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ integer(c_int), value  :: flags
+ end function rsb_get_coo
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_get_csr&
+  &(matrix,VA,RP,JA,flags)&
+  &bind(c,name='rsb_get_csr')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: VA(*)
+ type(c_ptr), value  :: RP
+ integer(c_int) :: JA(*)
+ integer(c_int), value  :: flags
+ end function rsb_get_csr
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_getdiag&
+  &(matrix,diagonal)&
+  &bind(c,name='rsb_getdiag')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: diagonal(*)
+ end function rsb_getdiag
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_get_rows_sparse&
+  &(matrix,VA,fr,lr,IA,JA,rnz,flags)&
+  &bind(c,name='rsb_get_rows_sparse')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: VA(*)
+ integer(c_int), value  :: fr
+ integer(c_int), value  :: lr
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ integer(c_int) :: rnz
+ integer(c_int), value  :: flags
+ end function rsb_get_rows_sparse
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_get_block_sparse_pattern&
+  &(matrix,fr,lr,fc,lc,IA,JA,IREN,JREN,rnz,flags)&
+  &bind(c,name='rsb_get_block_sparse_pattern')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ integer(c_int), value  :: fr
+ integer(c_int), value  :: lr
+ integer(c_int), value  :: fc
+ integer(c_int), value  :: lc
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ type(c_ptr), value  :: IREN
+ type(c_ptr), value  :: JREN
+ integer(c_int) :: rnz
+ integer(c_int), value  :: flags
+ end function rsb_get_block_sparse_pattern
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_get_block_sparse&
+  &(matrix,VA,fr,lr,fc,lc,IA,JA,IREN,JREN,rnz,flags)&
+  &bind(c,name='rsb_get_block_sparse')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: VA(*)
+ integer(c_int), value  :: fr
+ integer(c_int), value  :: lr
+ integer(c_int), value  :: fc
+ integer(c_int), value  :: lc
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ type(c_ptr), value  :: IREN
+ type(c_ptr), value  :: JREN
+ integer(c_int) :: rnz
+ integer(c_int), value  :: flags
+ end function rsb_get_block_sparse
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_get_columns_sparse&
+  &(matrix,columns,fc,lc,IA,JA,rnz,flags)&
+  &bind(c,name='rsb_get_columns_sparse')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ type(c_ptr), value  :: columns
+ integer(c_int), value  :: fc
+ integer(c_int), value  :: lc
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ integer(c_int) :: rnz
+ integer(c_int), value  :: flags
+ end function rsb_get_columns_sparse
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_get_matrix_nnz&
+  &(matrix)&
+  &bind(c,name='rsb_get_matrix_nnz')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ end function rsb_get_matrix_nnz
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_get_matrix_n_rows&
+  &(matrix)&
+  &bind(c,name='rsb_get_matrix_n_rows')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ end function rsb_get_matrix_n_rows
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_get_matrix_n_columns&
+  &(matrix)&
+  &bind(c,name='rsb_get_matrix_n_columns')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ end function rsb_get_matrix_n_columns
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_sizeof&
+  &(matrix)&
+  &bind(c,name='rsb_sizeof')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ end function rsb_sizeof
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_get_block_nnz&
+  &(matrix,fr,lr,fc,lc,flags,errvalp)&
+  &bind(c,name='rsb_get_block_nnz')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ integer(c_int), value  :: fr
+ integer(c_int), value  :: lr
+ integer(c_int), value  :: fc
+ integer(c_int), value  :: lc
+ integer(c_int), value  :: flags
+ integer(c_int) :: errvalp
+ end function rsb_get_block_nnz
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_get_rows_nnz&
+  &(matrix,fr,lr,flags,errvalp)&
+  &bind(c,name='rsb_get_rows_nnz')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ integer(c_int), value  :: fr
+ integer(c_int), value  :: lr
+ integer(c_int), value  :: flags
+ integer(c_int) :: errvalp
+ end function rsb_get_rows_nnz
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_assign&
+  &(new_matrix,matrix)&
+  &bind(c,name='rsb_assign')
+use iso_c_binding
+ type(c_ptr), value  :: new_matrix
+ type(c_ptr), value  :: matrix
+ end function rsb_assign
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_transpose&
+  &(matrixp)&
+  &bind(c,name='rsb_transpose')
+use iso_c_binding
+ type(c_ptr), value  :: matrixp
+ end function rsb_transpose
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_htranspose&
+  &(matrixp)&
+  &bind(c,name='rsb_htranspose')
+use iso_c_binding
+ type(c_ptr), value  :: matrixp
+ end function rsb_htranspose
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_elemental_scale&
+  &(matrix,alphap)&
+  &bind(c,name='rsb_elemental_scale')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: alphap
+ end function rsb_elemental_scale
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_elemental_scale_inv&
+  &(matrix,alphap)&
+  &bind(c,name='rsb_elemental_scale_inv')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: alphap
+ end function rsb_elemental_scale_inv
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_elemental_pow&
+  &(matrix,alphap)&
+  &bind(c,name='rsb_elemental_pow')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: alphap
+ end function rsb_elemental_pow
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_update_elements&
+  &(matrix,VA,IA,JA,nnz,flags)&
+  &bind(c,name='rsb_update_elements')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: VA(*)
+ integer(c_int) :: IA(*)
+ integer(c_int) :: JA(*)
+ integer(c_int), value  :: nnz
+ integer(c_int), value  :: flags
+ end function rsb_update_elements
+end interface
+
+interface
+type(c_ptr) function &
+  &rsb_matrix_sum&
+  &(matrixa,alphap,transa,matrixb,betap,transb,errvalp)&
+  &bind(c,name='rsb_matrix_sum')
+use iso_c_binding
+ type(c_ptr), value  :: matrixa
+ real(c_double) :: alphap
+ integer(c_int), value  :: transa
+ type(c_ptr), value  :: matrixb
+ real(c_double) :: betap
+ integer(c_int), value  :: transb
+ integer(c_int) :: errvalp
+ end function rsb_matrix_sum
+end interface
+
+interface
+type(c_ptr) function &
+  &rsb_matrix_mul&
+  &(matrixa,alphap,transa,matrixb,betap,transb,errvalp)&
+  &bind(c,name='rsb_matrix_mul')
+use iso_c_binding
+ type(c_ptr), value  :: matrixa
+ real(c_double) :: alphap
+ integer(c_int), value  :: transa
+ type(c_ptr), value  :: matrixb
+ real(c_double) :: betap
+ integer(c_int), value  :: transb
+ integer(c_int) :: errvalp
+ end function rsb_matrix_mul
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_negation&
+  &(matrix)&
+  &bind(c,name='rsb_negation')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ end function rsb_negation
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_scal&
+  &(matrix,d,transa)&
+  &bind(c,name='rsb_scal')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: d(*)
+ integer(c_int), value  :: transa
+ end function rsb_scal
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_scale_rows&
+  &(matrix,d)&
+  &bind(c,name='rsb_scale_rows')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ real(c_double) :: d(*)
+ end function rsb_scale_rows
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_reinit_matrix&
+  &(matrix)&
+  &bind(c,name='rsb_reinit_matrix')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ end function rsb_reinit_matrix
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_psblas_trans_to_rsb_trans&
+  &(trans)&
+  &bind(c,name='rsb_psblas_trans_to_rsb_trans')
+use iso_c_binding
+ character(c_char), value  :: trans
+ end function rsb_psblas_trans_to_rsb_trans
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_print_matrix_t&
+  &(matrix)&
+  &bind(c,name='rsb_print_matrix_t')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ end function rsb_print_matrix_t
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_print_matrix_unsorted_coo&
+  &(matrix)&
+  &bind(c,name='rsb_print_matrix_unsorted_coo')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ end function rsb_print_matrix_unsorted_coo
+end interface
+
+interface
+integer(c_int) function &
+  &rsb_save_matrix_file_as_matrix_market&
+  &(matrix,filename)&
+  &bind(c,name='rsb_save_matrix_file_as_matrix_market')
+use iso_c_binding
+ type(c_ptr), value  :: matrix
+ type(c_ptr), value  :: filename
+ end function rsb_save_matrix_file_as_matrix_market
+end interface
+
+interface
+type(c_ptr) function &
+  &rsb_load_matrix_file_as_matrix_market&
+  &(filename,flags,typecode,errvalp)&
+  &bind(c,name='rsb_load_matrix_file_as_matrix_market')
+use iso_c_binding
+ type(c_ptr), value  :: filename
+ integer(c_int), value  :: flags
+ integer(c_int), value  :: typecode
+ integer(c_int) :: errvalp
+ end function rsb_load_matrix_file_as_matrix_market
+end interface
+end module rsb_mod
diff --git a/test/newfmt/Makefile b/test/newfmt/Makefile
index c8bc953e..f835e6db 100644
--- a/test/newfmt/Makefile
+++ b/test/newfmt/Makefile
@@ -18,7 +18,7 @@ EXEDIR=./runs
 all: ppde spde
 
 ppde: ppde.o 
-	$(F90LINK)  ppde.o -o ppde $(PSBLAS_LIB) $(LDLIBS)
+	$(F90LINK)  ppde.o -o ppde $(PSBLAS_LIB) $(LDLIBS) -lgomp
 	/bin/mv ppde $(EXEDIR)
 
 
diff --git a/test/newfmt/ppde.f90 b/test/newfmt/ppde.f90
index 808e51a5..6fd6000c 100644
--- a/test/newfmt/ppde.f90
+++ b/test/newfmt/ppde.f90
@@ -66,6 +66,7 @@ program ppde
   use psb_prec_mod
   use psb_krylov_mod
   use psb_d_ell_mat_mod
+  use psb_d_rsb_mat_mod
   implicit none
 
   ! input parameters
@@ -78,9 +79,10 @@ program ppde
   real(psb_dpk_) :: t1, t2, tprec 
 
   ! sparse matrix and preconditioner
-  type(psb_dspmat_type) :: a
+  type(psb_dspmat_type) :: a,bm
   type(psb_dprec_type)  :: prec
   type(psb_d_ell_sparse_mat) :: aell
+  type(psb_d_rsb_sparse_mat) :: arsb
   ! descriptor
   type(psb_desc_type)   :: desc_a
   ! dense matrices
@@ -95,7 +97,7 @@ program ppde
 
   ! other variables
   integer            :: info, i
-  character(len=20)  :: name,ch_err
+  character(len=20)  :: name,ch_err, fname
 
   info=psb_success_
 
@@ -121,7 +123,7 @@ program ppde
   !
   call psb_barrier(ictxt)
   t1 = psb_wtime()
-  call create_matrix(idim,a,b,x,desc_a,ictxt,afmt,info,mold=aell)  
+  call create_matrix(idim,a,b,x,desc_a,ictxt,afmt,info,mold=arsb)  
   call psb_barrier(ictxt)
   t2 = psb_wtime() - t1
   if(info /= psb_success_) then
@@ -130,6 +132,10 @@ program ppde
     call psb_errpush(info,name,a_err=ch_err)
     goto 9999
   end if
+
+  call a%cscnv(bm,info,type='CSR')
+  write(fname,'(a,i2.2,a,i2.2,a)') 'mat',iam,'-',np,'.mtx'
+  call bm%print(fname,head='%Test sparse gen RSB')
   if (iam == psb_root_) write(psb_out_unit,'("Overall matrix creation time : ",es12.5)')t2
   if (iam == psb_root_) write(psb_out_unit,'(" ")')
   !
diff --git a/test/newfmt/runs/ppde.inp b/test/newfmt/runs/ppde.inp
index 43c4a5b2..bf76e9cc 100644
--- a/test/newfmt/runs/ppde.inp
+++ b/test/newfmt/runs/ppde.inp
@@ -1,8 +1,8 @@
 7                      Number of entries below this
 BICGSTAB                 Iterative method BICGSTAB CGS  BICG BICGSTABL RGMRES
-BJAC                   Preconditioner NONE  DIAG  BJAC 
+DIAG                   Preconditioner NONE  DIAG  BJAC 
 CSR                    Storage format for matrix A:  CSR COO JAD 
-060                    Domain size (acutal system is this**3)
+020                    Domain size (acutal system is this**3)
 2                      Stopping criterion
 0100                   MAXIT
 -1                     ITRACE
diff --git a/test/serial/d_matgen.f03 b/test/serial/d_matgen.f03
index e3b83313..439e9cd3 100644
--- a/test/serial/d_matgen.f03
+++ b/test/serial/d_matgen.f03
@@ -427,16 +427,16 @@ contains
    !write (*,*) acxx%val
    !write (*,*) diag
     
-!!$    t1 = psb_wtime()
-!!$    call a_n%cscnv(info,mold=acsr)
-!!$
-!!$    if(info /= psb_success_) then
-!!$      info=psb_err_from_subroutine_
-!!$      ch_err='asb rout.'
-!!$      call psb_errpush(info,name,a_err=ch_err)
-!!$      goto 9999
-!!$    end if
-!!$    tmov = psb_wtime()-t1
+    t1 = psb_wtime()
+    call a_n%cscnv(info,mold=acsr)
+
+    if(info /= psb_success_) then
+      info=psb_err_from_subroutine_
+      ch_err='asb rout.'
+      call psb_errpush(info,name,a_err=ch_err)
+      goto 9999
+    end if
+    tmov = psb_wtime()-t1
 ! !$    call a_n%print(21)
     anorm = a_n%csnmi()
     write(psb_err_unit,*) 'Nrm infinity ',anorm
diff --git a/test/serial/psb_d_rsb_mat_mod.F03 b/test/serial/psb_d_rsb_mat_mod.F03
index 4d1cce61..f30c4ce9 100644
--- a/test/serial/psb_d_rsb_mat_mod.F03
+++ b/test/serial/psb_d_rsb_mat_mod.F03
@@ -86,7 +86,7 @@ module psb_d_rsb_mat_mod
     !PSBRSB_DEBUG('')
     res=-1 ! FIXME
 #ifdef HAVE_LIBRSB
-    res=d_rsb_to_psb_info(rsb_init())
+    res=d_rsb_to_psb_info(rsb_init(c_null_ptr))
 #endif
   end function psv_rsb_mat_init
 
@@ -301,7 +301,7 @@ subroutine psb_d_rsb_csmm(alpha,a,x,beta,y,info,trans)
   end if
   ldx=size(x,1); ldy=size(y,1)
   nc=min(size(x,2),size(y,2) )
-  info=d_rsb_to_psb_info(rsb_spmm(a%rsbmptr,x,y,ldx,ldy,nc,rsb_psblas_trans_to_rsb_trans(trans_),alpha,beta,c_f_order))
+!!$  info=d_rsb_to_psb_info(rsb_spmm(a%rsbmptr,x,y,ldx,ldy,nc,rsb_psblas_trans_to_rsb_trans(trans_),alpha,beta,c_f_order))
 end subroutine psb_d_rsb_csmm
 
 subroutine psb_d_rsb_cssm(alpha,a,x,beta,y,info,trans) 
@@ -322,7 +322,7 @@ subroutine psb_d_rsb_cssm(alpha,a,x,beta,y,info,trans)
   end if
   ldx=size(x,1); ldy=size(y,1)
   nc=min(size(x,2),size(y,2) )
-  info=d_rsb_to_psb_info(rsb_spsm(a%rsbmptr,y,ldy,nc,rsb_psblas_trans_to_rsb_trans(trans_),alpha,beta,c_f_order))
+!!$  info=d_rsb_to_psb_info(rsb_spsm(a%rsbmptr,y,ldy,nc,rsb_psblas_trans_to_rsb_trans(trans_),alpha,beta,c_f_order))
 end subroutine
 
 subroutine psb_d_rsb_rowsum(d,a)