mld2p4:

17 years ago · 1584f7d8a8
parent 3cb6e674c3
commit 1584f7d8a8
34 changed files with 492 additions and 473 deletions
--- a/Make.inc
+++ b/Make.inc
@ -14,6 +14,15 @@ include $(PSBLASDIR)/Make.inc
 #  (compiled with) the compilers specified in the        #
 #  PSBLAS main Make.inc                                  #
 #                                                        #
 #  Examples:                                             #
 #  UMFLIBS=-lumfpack -lamd -L/path/to/UMFPACK            #
 #  UMFFLAGS=-DHave_UMF_ -I/path/to/UMFPACK               #
 #                                                        #
 #  SLULIBS=-lslu   -L/path/to/SuperLU                    #
 #  SLUFLAGS=-DHave_SLU_ -I/path/to/SuperLU               #
 #                                                        #
 #  SLUDISTLIBS=-lslud -L/path/to/SuperLUDist             #
 #  SLUDISTFLAGS=-DHave_SLUDist_ -I/path/to/SuperLUDist   #
 #                                                        #
 ##########################################################
@ -24,8 +33,8 @@ SLUDISTLIBS=-lslud -L/usr/local/SLUDist_2.0/gnu43
 SLUDISTFLAGS=-DHave_SLUDist_ -I/usr/local/SLUDist_2.0/gnu43
-UMFLIBS=-lumfpack -lamd -L/home/sfilippo/LIB/Umfpack_gcc41
+UMFLIBS=-lumfpack -lamd -L/home/sfilippo/LIB/Umfpack_gcc41/
-UMFFLAGS=-DHave_UMF_ -I/home/sfilippo/LIB/Umfpack_gcc41
+UMFFLAGS=-DHave_UMF_ -I/home/sfilippo/LIB/Umfpack_gcc41/
 LDLIBS=$(SLULIBS) $(SLUDISTLIBS) $(UMFLIBS) $(PSBLDLIBS)
--- a/mlprec/mld_basep_bld_mod.f90
+++ b/mlprec/mld_basep_bld_mod.f90
@ -229,7 +229,7 @@ module mld_basep_bld_mod
      type(psb_dspmat_type),intent(in)    :: a
      type(psb_dspmat_type),intent(inout) :: l,u
      type(psb_dspmat_type),intent(in), optional, target :: blck
-      real(kind(1.d0)), intent(inout)     ::  d(:)
+      real(psb_dpk_), intent(inout)     ::  d(:)
    end subroutine mld_dilu0_fact
    subroutine mld_zilu0_fact(ialg,a,l,u,d,info,blck)
      use psb_base_mod
@ -238,7 +238,7 @@ module mld_basep_bld_mod
      type(psb_zspmat_type),intent(in)    :: a
      type(psb_zspmat_type),intent(inout) :: l,u
      type(psb_zspmat_type),intent(in), optional, target :: blck
-      complex(kind(1.d0)), intent(inout)     ::  d(:)
+      complex(psb_dpk_), intent(inout)     ::  d(:)
    end subroutine mld_zilu0_fact
  end interface
@ -250,7 +250,7 @@ module mld_basep_bld_mod
      type(psb_dspmat_type),intent(in)    :: a
      type(psb_dspmat_type),intent(inout) :: l,u
      type(psb_dspmat_type),intent(in), optional, target :: blck
-      real(kind(1.d0)), intent(inout)     ::  d(:)
+      real(psb_dpk_), intent(inout)     ::  d(:)
    end subroutine mld_diluk_fact
    subroutine mld_ziluk_fact(fill_in,ialg,a,l,u,d,info,blck)
      use psb_base_mod
@ -259,7 +259,7 @@ module mld_basep_bld_mod
      type(psb_zspmat_type),intent(in)    :: a
      type(psb_zspmat_type),intent(inout) :: l,u
      type(psb_zspmat_type),intent(in), optional, target :: blck
-      complex(kind(1.d0)), intent(inout)     ::  d(:)
+      complex(psb_dpk_), intent(inout)     ::  d(:)
    end subroutine mld_ziluk_fact
  end interface
@ -267,22 +267,22 @@ module mld_basep_bld_mod
    subroutine mld_dilut_fact(fill_in,thres,a,l,u,d,info,blck)
      use psb_base_mod
      integer, intent(in)                 :: fill_in
-      real(kind(1.d0)), intent(in)        :: thres
+      real(psb_dpk_), intent(in)        :: thres
      integer, intent(out)                :: info
      type(psb_dspmat_type),intent(in)    :: a
      type(psb_dspmat_type),intent(inout) :: l,u
      type(psb_dspmat_type),intent(in), optional, target :: blck
-      real(kind(1.d0)), intent(inout)     ::  d(:)
+      real(psb_dpk_), intent(inout)     ::  d(:)
    end subroutine mld_dilut_fact
    subroutine mld_zilut_fact(fill_in,thres,a,l,u,d,info,blck)
      use psb_base_mod
      integer, intent(in)                 :: fill_in
-      real(kind(1.d0)), intent(in)        :: thres
+      real(psb_dpk_), intent(in)        :: thres
      integer, intent(out)                :: info
      type(psb_zspmat_type),intent(in)    :: a
      type(psb_zspmat_type),intent(inout) :: l,u
      type(psb_zspmat_type),intent(in), optional, target :: blck
-      complex(kind(1.d0)), intent(inout)  ::  d(:)
+      complex(psb_dpk_), intent(inout)  ::  d(:)
    end subroutine mld_zilut_fact
  end interface
 end module mld_basep_bld_mod
--- a/mlprec/mld_daggrmat_smth_asb.F90
+++ b/mlprec/mld_daggrmat_smth_asb.F90
@ -126,7 +126,7 @@ subroutine mld_daggrmat_smth_asb(a,desc_a,ac,desc_ac,p,info)
  logical            :: ml_global_nmb
  integer            :: debug_level, debug_unit
  integer, parameter :: ncmax=16
-  real(kind(1.d0))   :: omega, anorm, tmp, dg
+  real(psb_dpk_)   :: omega, anorm, tmp, dg
  name='mld_aggrmat_smth_asb'
  if(psb_get_errstatus().ne.0) return 
@ -190,7 +190,7 @@ subroutine mld_daggrmat_smth_asb(a,desc_a,ac,desc_ac,p,info)
  if (info /= 0) then 
    info=4025
    call psb_errpush(info,name,i_err=(/nrow,0,0,0,0/),&
-         & a_err='real(kind(1.d0))')
+         & a_err='real(psb_dpk_)')
    goto 9999      
  end if
--- a/mlprec/mld_das_aply.f90
+++ b/mlprec/mld_das_aply.f90
@ -52,16 +52,16 @@
 !
 !
 ! Arguments:
-!   alpha      -  real(kind(0.d0)), input.
+!   alpha      -  real(psb_dpk_), input.
 !                 The scalar alpha.
 !   prec       -  type(mld_dbaseprc_type), input.
 !                 The base preconditioner data structure containing the local part
 !                 of the preconditioner K.
-!   x          -  real(kind(0.d0)), dimension(:), input.
+!   x          -  real(psb_dpk_), dimension(:), input.
 !                 The local part of the vector X.
-!   beta       -  real(kind(0.d0)), input.
+!   beta       -  real(psb_dpk_), input.
 !                 The scalar beta.
-!   y          -  real(kind(0.d0)), dimension(:), input/output.
+!   y          -  real(psb_dpk_), dimension(:), input/output.
 !                 The local part of the vector Y.
 !   desc_data  -  type(psb_desc_type), input.
 !                 The communication descriptor associated to the matrix to be
@ -69,7 +69,7 @@
 !   trans      -  character, optional.
 !                 If trans='N','n' then op(K^(-1)) = K^(-1);
 !                 if trans='T','t' then op(K^(-1)) = K^(-T) (transpose of K^(-1)).
-!   work       -  real(kind(0.d0)), dimension (:), optional, target.
+!   work       -  real(psb_dpk_), dimension (:), optional, target.
 !                 Workspace. Its size must be at least 4*psb_cd_get_local_cols(desc_data).
 !   info       -  integer, output.
 !                 Error code.
@ -84,16 +84,16 @@ subroutine mld_das_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
  ! Arguments
  type(psb_desc_type),intent(in)      :: desc_data
  type(mld_dbaseprc_type), intent(in) :: prec
-  real(kind(0.d0)),intent(in)         :: x(:)
+  real(psb_dpk_),intent(in)         :: x(:)
-  real(kind(0.d0)),intent(inout)      :: y(:)
+  real(psb_dpk_),intent(inout)      :: y(:)
-  real(kind(0.d0)),intent(in)         :: alpha,beta
+  real(psb_dpk_),intent(in)         :: alpha,beta
  character(len=1)                    :: trans
-  real(kind(0.d0)),target             :: work(:)
+  real(psb_dpk_),target             :: work(:)
  integer, intent(out)                :: info
  ! Local variables
  integer :: n_row,n_col, int_err(5), nrow_d
-  real(kind(1.d0)), pointer :: ww(:), aux(:), tx(:),ty(:)
+  real(psb_dpk_), pointer :: ww(:), aux(:), tx(:),ty(:)
  integer           :: ictxt,np,me,isz, err_act
  character(len=20) :: name, ch_err
  character         :: trans_
@ -154,7 +154,7 @@ subroutine mld_das_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
        allocate(ww(isz),tx(isz),ty(isz),stat=info)
        if (info /= 0) then 
          call psb_errpush(4025,name,i_err=(/3*isz,0,0,0,0/),&
-               & a_err='real(kind(1.d0))')
+               & a_err='real(psb_dpk_)')
          goto 9999      
        end if
      else if ((3*isz) <= size(work)) then 
@ -164,7 +164,7 @@ subroutine mld_das_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
        allocate(aux(4*isz),stat=info)
        if (info /= 0) then 
          call psb_errpush(4025,name,i_err=(/4*isz,0,0,0,0/),&
-               & a_err='real(kind(1.d0))')
+               & a_err='real(psb_dpk_)')
          goto 9999      
        end if
      else 
@ -172,7 +172,7 @@ subroutine mld_das_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
             &aux(4*isz),stat=info)
        if (info /= 0) then 
          call psb_errpush(4025,name,i_err=(/4*isz,0,0,0,0/),&
-               & a_err='real(kind(1.d0))')
+               & a_err='real(psb_dpk_)')
          goto 9999      
        end if
--- a/mlprec/mld_dbaseprec_aply.f90
+++ b/mlprec/mld_dbaseprec_aply.f90
@ -56,16 +56,16 @@
 !
 !
 ! Arguments:
-!   alpha      -  real(kind(0.d0)), input.
+!   alpha      -  real(psb_dpk_), input.
 !                 The scalar alpha.
 !   prec       -  type(mld_dbaseprc_type), input.
 !                 The base preconditioner data structure containing the local part
 !                 of the preconditioner K.
-!   x          -  real(kind(0.d0)), dimension(:), input.
+!   x          -  real(psb_dpk_), dimension(:), input.
 !                 The local part of the vector X.
-!   beta       -  real(kind(0.d0)), input.
+!   beta       -  real(psb_dpk_), input.
 !                 The scalar beta.
-!   y          -  real(kind(0.d0)), dimension(:), input/output.
+!   y          -  real(psb_dpk_), dimension(:), input/output.
 !                 The local part of the vector Y.
 !   desc_data  -  type(psb_desc_type), input.
 !                 The communication descriptor associated to the matrix to be
@ -73,7 +73,7 @@
 !   trans      -  character, optional.
 !                 If trans='N','n' then op(K^(-1)) = K^(-1);
 !                 if trans='T','t' then op(K^(-1)) = K^(-T) (transpose of K^(-1)).
-!   work       -  real(kind(0.d0)), dimension (:), optional, target.
+!   work       -  real(psb_dpk_), dimension (:), optional, target.
 !                 Workspace. Its size must be at least 4*psb_cd_get_local_cols(desc_data).
 !   info       -  integer, output.
 !                 Error code.
@ -88,15 +88,15 @@ subroutine mld_dbaseprec_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
  ! Arguments
  type(psb_desc_type),intent(in)      :: desc_data
  type(mld_dbaseprc_type), intent(in) :: prec
-  real(kind(0.d0)),intent(in)         :: x(:)
+  real(psb_dpk_),intent(in)         :: x(:)
-  real(kind(0.d0)),intent(inout)      :: y(:)
+  real(psb_dpk_),intent(inout)      :: y(:)
-  real(kind(0.d0)),intent(in)         :: alpha,beta
+  real(psb_dpk_),intent(in)         :: alpha,beta
  character(len=1)                    :: trans
-  real(kind(0.d0)),target             :: work(:)
+  real(psb_dpk_),target             :: work(:)
  integer, intent(out)                :: info
  ! Local variables
-  real(kind(1.d0)), pointer :: ww(:)
+  real(psb_dpk_), pointer :: ww(:)
  integer           :: ictxt, np, me, err_act
  integer           :: n_row, int_err(5)
  character(len=20) :: name, ch_err
@ -140,7 +140,7 @@ subroutine mld_dbaseprec_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
    else
      allocate(ww(size(x)),stat=info)
      if (info /= 0) then 
-        call psb_errpush(4025,name,i_err=(/size(x),0,0,0,0/),a_err='real(kind(1.d0))')
+        call psb_errpush(4025,name,i_err=(/size(x),0,0,0,0/),a_err='real(psb_dpk_)')
        goto 9999      
      end if
    end if
--- a/mlprec/mld_dilu0_fact.f90
+++ b/mlprec/mld_dilu0_fact.f90
@ -85,7 +85,7 @@
 !               The U factor (except its diagonal) in the incomplete factorization.
 !               Note: its allocation is managed by the calling routine mld_ilu_bld,
 !               hence it cannot be only intent(out).
-!    d       -  real(kind(1.d0)), dimension(:), input/output.
+!    d       -  real(psb_dpk_), dimension(:), input/output.
 !               The inverse of the diagonal entries of the U factor in the incomplete
 !               factorization.
 !               Note: its allocation is managed by the calling routine mld_ilu_bld,
@ -110,7 +110,7 @@ subroutine mld_dilu0_fact(ialg,a,l,u,d,info,blck)
  integer, intent(in)                 :: ialg
  type(psb_dspmat_type),intent(in)    :: a
  type(psb_dspmat_type),intent(inout) :: l,u
-  real(kind(1.d0)), intent(inout)     :: d(:)
+  real(psb_dpk_), intent(inout)     :: d(:)
  integer, intent(out)                :: info
  type(psb_dspmat_type),intent(in), optional, target :: blck
@ -249,10 +249,10 @@ contains
  !               to build an Additive Schwarz base preconditioner with overlap
  !               greater than 0. If the overlap is 0 or the matrix has been
  !               reordered (see mld_fact_bld), then b does not contain any row.
-  !    d       -  real(kind(1.d0)), dimension(:), output.
+  !    d       -  real(psb_dpk_), dimension(:), output.
  !               The inverse of the diagonal entries of the U factor in the
  !               incomplete factorization.
-  !    laspk   -  real(kind(1.d0)), dimension(:), input/output.
+  !    laspk   -  real(psb_dpk_), dimension(:), input/output.
  !               The entries of U are stored according to the CSR format.
  !               The L factor in the incomplete factorization.
  !    lia1    -  integer, dimension(:), input/output.
@ -261,7 +261,7 @@ contains
  !    lia2    -  integer, dimension(:), input/output.
  !               The indices identifying the first nonzero entry of each row
  !               of the L factor in laspk, according to the CSR storage format. 
-  !    uaspk   -  real(kind(1.d0)), dimension(:), input/output.
+  !    uaspk   -  real(psb_dpk_), dimension(:), input/output.
  !               The U factor in the incomplete factorization.
  !               The entries of U are stored according to the CSR format.
  !    uia1    -  integer, dimension(:), input/output.
@ -288,11 +288,11 @@ contains
    integer,intent(inout)            :: m,l1,l2,info
    integer, intent(in)              :: ma,mb
    integer, dimension(:), intent(inout) :: lia1,lia2,uia1,uia2
-    real(kind(1.d0)), dimension(:),intent(inout)  :: laspk,uaspk,d
+    real(psb_dpk_), dimension(:),intent(inout)  :: laspk,uaspk,d
    ! Local variables
    integer :: i,j,k,l,low1,low2,kk,jj,ll, ktrw,err_act
-    real(kind(1.d0)) :: dia,temp
+    real(psb_dpk_) :: dia,temp
    integer, parameter :: nrb=16
    type(psb_dspmat_type) :: trw
    integer             :: int_err(5) 
@ -519,10 +519,10 @@ contains
  !               The column indices of the nonzero entries of the lower triangle
  !               copied in laspk row by row (see mld_dilu0_factint), according
  !               to the CSR storage format.
-  !    laspk   -  real(kind(1.d0)), dimension(:), input/output.
+  !    laspk   -  real(psb_dpk_), dimension(:), input/output.
  !               The array where the entries of the row corresponding to the
  !               lower triangle are copied.
-  !    dia     -  real(kind(1.d0)), output.
+  !    dia     -  real(psb_dpk_), output.
  !               The diagonal entry of the copied row.
  !    l2      -  integer, input/output.
  !               Pointer to the last occupied entry of uaspk.
@ -530,7 +530,7 @@ contains
  !               The column indices of the nonzero entries of the upper triangle
  !               copied in uaspk row by row (see mld_dilu0_factint), according
  !               to the CSR storage format.
-  !    uaspk   -  real(kind(1.d0)), dimension(:), input/output.
+  !    uaspk   -  real(psb_dpk_), dimension(:), input/output.
  !               The array where the entries of the row corresponding to the
  !               upper triangle are copied. 
  !    ktrw    -  integer, input/output.
@ -557,7 +557,7 @@ contains
    integer, intent(in)                  :: i,m,jd,jmin,jmax
    integer, intent(inout)               :: ktrw,l1,l2
    integer, intent(inout)               :: lia1(:), uia1(:)
-    real(kind(1.d0)), intent(inout)      :: laspk(:), uaspk(:), dia
+    real(psb_dpk_), intent(inout)      :: laspk(:), uaspk(:), dia
    ! Local variables
    integer               :: k,j,info,irb
--- a/mlprec/mld_diluk_fact.f90
+++ b/mlprec/mld_diluk_fact.f90
@ -82,7 +82,7 @@
 !               The U factor (except its diagonal) in the incomplete factorization.
 !               Note: its allocation is managed by the calling routine mld_ilu_bld,
 !               hence it cannot be only intent(out).
-!    d       -  real(kind(1.d0)), dimension(:), input/output.
+!    d       -  real(psb_dpk_), dimension(:), input/output.
 !               The inverse of the diagonal entries of the U factor in the incomplete
 !               factorization.
 !               Note: its allocation is managed by the calling routine mld_ilu_bld,
@ -109,7 +109,7 @@ subroutine mld_diluk_fact(fill_in,ialg,a,l,u,d,info,blck)
  type(psb_dspmat_type),intent(in)    :: a
  type(psb_dspmat_type),intent(inout) :: l,u
  type(psb_dspmat_type),intent(in), optional, target :: blck
-  real(kind(1.d0)), intent(inout)     ::  d(:)
+  real(psb_dpk_), intent(inout)     ::  d(:)
  !     Local Variables
  integer   :: l1, l2, m, err_act
@ -237,10 +237,10 @@ contains
  !               to build an Additive Schwarz base preconditioner with overlap
  !               greater than 0. If the overlap is 0 or the matrix has been reordered
  !               (see mld_fact_bld), then b does not contain any row.
-  !    d       -  real(kind(1.d0)), dimension(:), output.
+  !    d       -  real(psb_dpk_), dimension(:), output.
  !               The inverse of the diagonal entries of the U factor in the incomplete
  !               factorization.
-  !    laspk   -  real(kind(1.d0)), dimension(:), input/output.
+  !    laspk   -  real(psb_dpk_), dimension(:), input/output.
  !               The L factor in the incomplete factorization.
  !    lia1    -  integer, dimension(:), input/output.
  !               The column indices of the nonzero entries of the L factor,
@ -248,7 +248,7 @@ contains
  !    lia2    -  integer, dimension(:), input/output.
  !               The indices identifying the first nonzero entry of each row
  !               of the L factor in laspk, according to the CSR storage format. 
-  !    uaspk   -  real(kind(1.d0)), dimension(:), input/output.
+  !    uaspk   -  real(psb_dpk_), dimension(:), input/output.
  !               The U factor in the incomplete factorization.
  !               The entries of U are stored according to the CSR format.
  !    uia1    -  integer, dimension(:), input/output.
@ -276,13 +276,13 @@ contains
    type(psb_dspmat_type), intent(in)            :: a,b
    integer, intent(inout)                       :: m,l1,l2,info
    integer, allocatable, intent(inout)          :: lia1(:),lia2(:),uia1(:),uia2(:)
-    real(kind(1.d0)), allocatable, intent(inout) :: laspk(:),uaspk(:)
+    real(psb_dpk_), allocatable, intent(inout) :: laspk(:),uaspk(:)
-    real(kind(1.d0)), intent(inout)              :: d(:)
+    real(psb_dpk_), intent(inout)              :: d(:)
  ! Local variables
    integer :: ma,mb,i, ktrw,err_act,nidx
    integer, allocatable          :: uplevs(:), rowlevs(:),idxs(:)
-    real(kind(1.d0)), allocatable :: row(:)
+    real(psb_dpk_), allocatable :: row(:)
    type(psb_int_heap) :: heap
    type(psb_dspmat_type) :: trw
    character(len=20), parameter  :: name='mld_diluk_factint'
@ -455,7 +455,7 @@ contains
  !               The maximum valid column index.
  !               The output matrix will contain a clipped copy taken from
  !               a(1:m,jmin:jmax).
-  !    row     -  real(kind(1.d0)), dimension(:), input/output.
+  !    row     -  real(psb_dpk_), dimension(:), input/output.
  !               In input it is the null vector (see mld_iluk_factint and
  !               iluk_copyout). In output it contains the row extracted
  !               from the matrix A. It actually contains a full row, i.e.
@ -493,7 +493,7 @@ contains
    integer, intent(in)                  :: i,m,jmin,jmax
    integer, intent(inout)               :: ktrw,info
    integer, intent(inout)               :: rowlevs(:)
-    real(kind(1.d0)), intent(inout)      :: row(:)
+    real(psb_dpk_), intent(inout)      :: row(:)
    type(psb_int_heap), intent(inout)    :: heap
  ! Local variables
@ -591,7 +591,7 @@ contains
  !    i       -  integer, input.
  !               The local index of the row to which the factorization is
  !               applied.
-  !    row     -  real(kind(1.d0)), dimension(:), input/output.
+  !    row     -  real(psb_dpk_), dimension(:), input/output.
  !               In input it contains the row to which the elimination step
  !               has to be applied. In output it contains the row after the
  !               elimination step. It actually contains a full row, i.e.
@ -608,7 +608,7 @@ contains
  !               in the processed row. In input it contains the indices concerning
  !               the row before the elimination step, while in output it contains
  !               the indices concerning the transformed row.
-  !    d       -  real(kind(1.d0)), input.
+  !    d       -  real(psb_dpk_), input.
  !               The inverse of the diagonal entries of the part of the U factor
  !               above the current row (see iluk_copyout).
  !    uia1    -  integer, dimension(:), input.
@ -621,7 +621,7 @@ contains
  !               the U factor above the current row, stored in uaspk row by row
  !               (see iluk_copyout, called by mld_diluk_factint), according to
  !               the CSR storage format.
-  !    uaspk   -  real(kind(1.d0)), dimension(:), input.
+  !    uaspk   -  real(psb_dpk_), dimension(:), input.
  !               The entries of the U factor above the current row (except the
  !               diagonal ones), stored according to the CSR format.
  !    uplevs  -  integer, dimension(:), input.
@ -651,11 +651,11 @@ contains
    integer, intent(inout)               :: rowlevs(:)
    integer, allocatable, intent(inout)  :: idxs(:)
    integer, intent(inout)               :: uia1(:),uia2(:),uplevs(:)
-    real(kind(1.d0)), intent(inout)      :: row(:), uaspk(:),d(:)
+    real(psb_dpk_), intent(inout)      :: row(:), uaspk(:),d(:)
    ! Local variables
    integer               :: k,j,lrwk,jj,lastk, iret
-    real(kind(1.d0))      :: rwk
+    real(psb_dpk_)      :: rwk
    info = 0
    if (.not.allocated(idxs)) then 
@ -762,7 +762,7 @@ contains
  !               The local index of the row to be copied.
  !    m       -  integer, input.
  !               The number of rows of the local matrix under factorization.
-  !    row     -  real(kind(1.d0)), dimension(:), input/output.
+  !    row     -  real(psb_dpk_), dimension(:), input/output.
  !               It contains, input, the row to be copied, and, in output,
  !               the null vector (the latter is used in the next call to
  !               iluk_copyin in mld_iluk_fact).
@ -792,10 +792,10 @@ contains
  !               The indices identifying the first nonzero entry of each row
  !               of the L factor, copied in laspk row by row (see 
  !               mld_diluk_factint), according to the CSR storage format.
-  !    laspk   -  real(kind(1.d0)), dimension(:), input/output.
+  !    laspk   -  real(psb_dpk_), dimension(:), input/output.
  !               The array where the entries of the row corresponding to the
  !               L factor are copied.
-  !    d       -  real(kind(1.d0)), dimension(:), input/output.
+  !    d       -  real(psb_dpk_), dimension(:), input/output.
  !               The array where the inverse of the diagonal entry of the
  !               row is copied (only d(i) is used by the routine). 
  !    uia1    -  integer, dimension(:), input/output.
@ -806,7 +806,7 @@ contains
  !               The indices identifying the first nonzero entry of each row
  !               of the U factor copied in uaspk row by row (see
  !               mld_dilu_fctint), according to the CSR storage format.
-  !    uaspk   -  real(kind(1.d0)), dimension(:), input/output.
+  !    uaspk   -  real(psb_dpk_), dimension(:), input/output.
  !               The array where the entries of the row corresponding to the
  !               U factor are copied.
  !    uplevs  -  integer, dimension(:), input.
@ -825,8 +825,8 @@ contains
    integer, intent(inout)                       :: l1, l2, info
    integer, intent(inout)                       :: rowlevs(:), idxs(:)
    integer, allocatable, intent(inout)          :: uia1(:), uia2(:), lia1(:), lia2(:),uplevs(:)
-    real(kind(1.d0)), allocatable, intent(inout) :: uaspk(:), laspk(:)
+    real(psb_dpk_), allocatable, intent(inout) :: uaspk(:), laspk(:)
-    real(kind(1.d0)), intent(inout)              :: row(:), d(:)
+    real(psb_dpk_), intent(inout)              :: row(:), d(:)
    ! Local variables
    integer               :: j,isz,err_act,int_err(5),idxp
--- a/mlprec/mld_dilut_fact.f90
+++ b/mlprec/mld_dilut_fact.f90
@ -78,7 +78,7 @@
 !               The U factor (except its diagonal) in the incomplete factorization.
 !               Note: its allocation is managed by the calling routine mld_ilu_bld,
 !               hence it cannot be only intent(out).
-!    d       -  real(kind(1.d0)), dimension(:), input/output.
+!    d       -  real(psb_dpk_), dimension(:), input/output.
 !               The inverse of the diagonal entries of the U factor in the incomplete
 !               factorization.
 !               Note: its allocation is managed by the calling routine mld_ilu_bld,
@ -101,11 +101,11 @@ subroutine mld_dilut_fact(fill_in,thres,a,l,u,d,info,blck)
  ! Arguments
  integer, intent(in)                 :: fill_in
-  real(kind(1.d0)), intent(in)        :: thres
+  real(psb_dpk_), intent(in)        :: thres
  integer, intent(out)                :: info
  type(psb_dspmat_type),intent(in)    :: a
  type(psb_dspmat_type),intent(inout) :: l,u
-  real(kind(1.d0)), intent(inout)     :: d(:)
+  real(psb_dpk_), intent(inout)     :: d(:)
  type(psb_dspmat_type),intent(in), optional, target :: blck
  !     Local Variables
@ -238,10 +238,10 @@ contains
  !               to build an Additive Schwarz base preconditioner with overlap
  !               greater than 0. If the overlap is 0 or the matrix has been reordered
  !               (see mld_fact_bld), then b does not contain any row.
-  !    d       -  real(kind(1.d0)), dimension(:), output.
+  !    d       -  real(psb_dpk_), dimension(:), output.
  !               The inverse of the diagonal entries of the U factor in the incomplete
  !               factorization.
-  !    laspk   -  real(kind(1.d0)), dimension(:), input/output.
+  !    laspk   -  real(psb_dpk_), dimension(:), input/output.
  !               The L factor in the incomplete factorization.
  !    lia1    -  integer, dimension(:), input/output.
  !               The column indices of the nonzero entries of the L factor,
@ -249,7 +249,7 @@ contains
  !    lia2    -  integer, dimension(:), input/output.
  !               The indices identifying the first nonzero entry of each row
  !               of the L factor in laspk, according to the CSR storage format. 
-  !    uaspk   -  real(kind(1.d0)), dimension(:), input/output.
+  !    uaspk   -  real(psb_dpk_), dimension(:), input/output.
  !               The U factor in the incomplete factorization.
  !               The entries of U are stored according to the CSR format.
  !    uia1    -  integer, dimension(:), input/output.
@ -274,18 +274,18 @@ contains
  ! Arguments
    integer, intent(in)                          :: fill_in
-    real(kind(1.d0)), intent(in)                 :: thres
+    real(psb_dpk_), intent(in)                 :: thres
    type(psb_dspmat_type), intent(in)            :: a,b
    integer, intent(inout)                       :: m,l1,l2,info
    integer, allocatable, intent(inout)          :: lia1(:),lia2(:),uia1(:),uia2(:)
-    real(kind(1.d0)), allocatable, intent(inout) :: laspk(:),uaspk(:)
+    real(psb_dpk_), allocatable, intent(inout) :: laspk(:),uaspk(:)
-    real(kind(1.d0)), intent(inout)              :: d(:)
+    real(psb_dpk_), intent(inout)              :: d(:)
    ! Local Variables
    integer :: i, ktrw,err_act,nidx,nlw,nup,jmaxup, ma, mb
-    real(kind(1.d0)) ::  nrmi
+    real(psb_dpk_) ::  nrmi
    integer, allocatable          :: idxs(:)
-    real(kind(1.d0)), allocatable :: row(:)
+    real(psb_dpk_), allocatable :: row(:)
    type(psb_int_heap) :: heap
    type(psb_dspmat_type) :: trw
    character(len=20), parameter  :: name='mld_dilut_factint'
@ -451,9 +451,9 @@ contains
  !    jmaxup  -  integer, output.
  !               The column index of the first entry with maximum absolute
  !               value in the part of the row belonging to the upper triangle
-  !    nrmi    -  real(kind(1.d0)), output.
+  !    nrmi    -  real(psb_dpk_), output.
  !               The 2-norm of the current row.
-  !    row     -  real(kind(1.d0)), dimension(:), input/output.
+  !    row     -  real(psb_dpk_), dimension(:), input/output.
  !               In input it is the null vector (see mld_ilut_factint and
  !               ilut_copyout). In output it contains the row extracted
  !               from the matrix A. It actually contains a full row, i.e.
@ -486,13 +486,13 @@ contains
    type(psb_dspmat_type), intent(inout) :: trw
    integer, intent(in)                  :: i, m,jmin,jmax,jd
    integer, intent(inout)               :: ktrw,nlw,nup,jmaxup,info
-    real(kind(1.d0)), intent(inout)      :: nrmi,row(:)
+    real(psb_dpk_), intent(inout)      :: nrmi,row(:)
    type(psb_int_heap), intent(inout)    :: heap
    integer               :: k,j,irb,kin,nz
    integer, parameter    :: nrb=16
-    real(kind(1.d0))      :: dmaxup
+    real(psb_dpk_)      :: dmaxup
-    real(kind(1.d0)), external    :: dnrm2
+    real(psb_dpk_), external    :: dnrm2
    character(len=20), parameter  :: name='mld_dilut_factint'
    if (psb_get_errstatus() /= 0) return 
@ -627,10 +627,10 @@ contains
  !               The threshold t, i.e. the drop tolerance, in ILU(k,t).
  !    i       -  integer, input.
  !               The local index of the row to which the factorization is applied.
-  !    nrmi    -  real(kind(1.d0)), input.
+  !    nrmi    -  real(psb_dpk_), input.
  !               The 2-norm of the row to which the elimination step has to be
  !               applied.
-  !    row     -  real(kind(1.d0)), dimension(:), input/output.
+  !    row     -  real(psb_dpk_), dimension(:), input/output.
  !               In input it contains the row to which the elimination step
  !               has to be applied. In output it contains the row after the
  !               elimination step. It actually contains a full row, i.e.
@ -641,7 +641,7 @@ contains
  !               the row before the elimination step, while in output it contains
  !               the previous indices plus the ones corresponding to transformed
  !               entries in the 'upper part' that have not been dropped.
-  !    d       -  real(kind(1.d0)), input.
+  !    d       -  real(psb_dpk_), input.
  !               The inverse of the diagonal entries of the part of the U factor
  !               above the current row (see ilut_copyout).
  !    uia1    -  integer, dimension(:), input.
@ -654,7 +654,7 @@ contains
  !               the U factor above the current row, stored in uaspk row by row
  !               (see ilut_copyout, called by mld_dilut_factint), according to
  !               the CSR storage format.
-  !    uaspk   -  real(kind(1.d0)), dimension(:), input.
+  !    uaspk   -  real(psb_dpk_), dimension(:), input.
  !               The entries of the U factor above the current row (except the
  !               diagonal ones), stored according to the CSR format.
  !    nidx    -  integer, output.
@ -678,14 +678,14 @@ contains
    type(psb_int_heap), intent(inout)   :: heap 
    integer, intent(in)                 :: i
    integer, intent(inout)              :: nidx,info
-    real(kind(1.d0)), intent(in)        :: thres,nrmi
+    real(psb_dpk_), intent(in)        :: thres,nrmi
    integer, allocatable, intent(inout) :: idxs(:)
    integer, intent(inout)              :: uia1(:),uia2(:)
-    real(kind(1.d0)), intent(inout)     :: row(:), uaspk(:),d(:)
+    real(psb_dpk_), intent(inout)     :: row(:), uaspk(:),d(:)
    ! Local Variables
    integer               :: k,j,jj,lastk,iret
-    real(kind(1.d0))      :: rwk
+    real(psb_dpk_)      :: rwk
    info  = 0
    call psb_ensure_size(200,idxs,info)
@ -810,9 +810,9 @@ contains
  !    jmaxup  -  integer, input.
  !               The column index of the first entry with maximum absolute
  !               value in the 'upper part' of the row in the initial matrix.
-  !    nrmi    -  real(kind(1.d0)), input.
+  !    nrmi    -  real(psb_dpk_), input.
  !               The 2-norm of the current row in the initial matrix.
-  !    row     -  real(kind(1.d0)), dimension(:), input/output.
+  !    row     -  real(psb_dpk_), dimension(:), input/output.
  !               It contains, input, the row to be copied, and, in output,
  !               the null vector (the latter is used in the next call to
  !               ilut_copyin in mld_ilut_fact).
@ -835,10 +835,10 @@ contains
  !               The indices identifying the first nonzero entry of each row
  !               of the L factor, copied in laspk row by row (see 
  !               mld_dilut_factint), according to the CSR storage format.
-  !    laspk   -  real(kind(1.d0)), dimension(:), input/output.
+  !    laspk   -  real(psb_dpk_), dimension(:), input/output.
  !               The array where the entries of the row corresponding to the
  !               L factor are copied.
-  !    d       -  real(kind(1.d0)), dimension(:), input/output.
+  !    d       -  real(psb_dpk_), dimension(:), input/output.
  !               The array where the inverse of the diagonal entry of the
  !               row is copied (only d(i) is used by the routine). 
  !    uia1    -  integer, dimension(:), input/output.
@ -849,7 +849,7 @@ contains
  !               The indices identifying the first nonzero entry of each row
  !               of the U factor copied in uaspk row by row (see
  !               mld_dilu_fctint), according to the CSR storage format.
-  !    uaspk   -  real(kind(1.d0)), dimension(:), input/output.
+  !    uaspk   -  real(psb_dpk_), dimension(:), input/output.
  !               The array where the entries of the row corresponding to the
  !               U factor are copied.
  !
@ -865,14 +865,14 @@ contains
    integer, intent(in)                         :: idxs(:)
    integer, intent(inout)                      :: l1,l2, info
    integer, allocatable, intent(inout)         :: uia1(:),uia2(:), lia1(:),lia2(:)
-    real(kind(1.d0)), intent(in)                :: thres,nrmi
+    real(psb_dpk_), intent(in)                :: thres,nrmi
-    real(kind(1.d0)),allocatable, intent(inout) :: uaspk(:), laspk(:)
+    real(psb_dpk_),allocatable, intent(inout) :: uaspk(:), laspk(:)
-    real(kind(1.d0)), intent(inout)             :: row(:), d(:)
+    real(psb_dpk_), intent(inout)             :: row(:), d(:)
    ! Local variables
-    real(kind(1.d0)),allocatable :: xw(:)
+    real(psb_dpk_),allocatable :: xw(:)
    integer, allocatable         :: xwid(:), indx(:)
-    real(kind(1.d0))             :: witem
+    real(psb_dpk_)             :: witem
    integer                      :: widx
    integer                      :: k,isz,err_act,int_err(5),idxp, nz
    type(psb_double_idx_heap)    :: heap
@ -899,7 +899,7 @@ contains
    if (info /= 0) then 
      info=4025
      call psb_errpush(info,name,i_err=(/3*nidx,0,0,0,0/),&
-           & a_err='real(kind(1.d0))')
+           & a_err='real(psb_dpk_)')
      goto 9999      
    end if
--- a/mlprec/mld_dmlprec_aply.f90
+++ b/mlprec/mld_dmlprec_aply.f90
@ -76,7 +76,7 @@
 !
 !
 ! Arguments:
-!   alpha      -   real(kind(0.d0)), input.
+!   alpha      -   real(psb_dpk_), input.
 !                  The scalar alpha.
 !   baseprecv  -   type(mld_dbaseprc_type), dimension(:), input.
 !                  The array of base preconditioner data structures containing the
@ -98,7 +98,7 @@
 !                                             It maps vectors (ilev) ---> (ilev-1).
 !         baseprecv(ilev)%av(mld_sm_pr_t_) -  The smoothed prolongator transpose.   
 !                                             It maps vectors (ilev-1) ---> (ilev).
-!      baseprecv(ilev)%d         -  real(kind(1.d0)), dimension(:), allocatable.
+!      baseprecv(ilev)%d         -  real(psb_dpk_), dimension(:), allocatable.
 !                                   The diagonal entries of the U factor in the ILU
 !                                   factorization of A(ilev).
 !      baseprecv(ilev)%desc_data -  type(psb_desc_type).
@ -111,7 +111,7 @@
 !      baseprecv(ilev)%iprcparm  -  integer, dimension(:), allocatable.
 !                                   The integer parameters defining the base
 !                                   preconditioner K(ilev).
-!      baseprecv(ilev)%dprcparm  -  real(kind(1.d0)), dimension(:), allocatable.
+!      baseprecv(ilev)%dprcparm  -  real(psb_dpk_), dimension(:), allocatable.
 !                                   The real parameters defining the base preconditioner
 !                                   K(ilev).
 !      baseprecv(ilev)%perm      -  integer, dimension(:), allocatable.
@ -138,14 +138,14 @@
 !      baseprecv(ilev)%base_desc -  type(psb_desc_type), pointer.
 !                                   Pointer to the communication descriptor associated
 !                                   to the sparse matrix pointed by base_a.  
-!      baseprecv(ilev)%dorig     -  real(kind(1.d0)), dimension(:), allocatable.
+!      baseprecv(ilev)%dorig     -  real(psb_dpk_), dimension(:), allocatable.
 !                                   Diagonal entries of the matrix pointed by base_a.
 !                  
-!   x          -  real(kind(0.d0)), dimension(:), input.
+!   x          -  real(psb_dpk_), dimension(:), input.
 !                 The local part of the vector X.
-!   beta       -  real(kind(0.d0)), input.
+!   beta       -  real(psb_dpk_), input.
 !                 The scalar beta.
-!   y          -  real(kind(0.d0)), dimension(:), input/output.
+!   y          -  real(psb_dpk_), dimension(:), input/output.
 !                 The local part of the vector Y.
 !   desc_data  -  type(psb_desc_type), input.
 !                 The communication descriptor associated to the matrix to be
@ -153,7 +153,7 @@
 !   trans      -  character, optional.
 !                 If trans='N','n' then op(M^(-1)) = M^(-1);
 !                 if trans='T','t' then op(M^(-1)) = M^(-T) (transpose of M^(-1)).
-!   work       -  real(kind(0.d0)), dimension (:), optional, target.
+!   work       -  real(psb_dpk_), dimension (:), optional, target.
 !                 Workspace. Its size must be at least 4*psb_cd_get_local_cols(desc_data).
 !   info       -  integer, output.
 !                 Error code.
@ -174,11 +174,11 @@ subroutine mld_dmlprec_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
  ! Arguments
  type(psb_desc_type),intent(in)      :: desc_data
  type(mld_dbaseprc_type), intent(in) :: baseprecv(:)
-  real(kind(0.d0)),intent(in)         :: alpha,beta
+  real(psb_dpk_),intent(in)         :: alpha,beta
-  real(kind(0.d0)),intent(in)         :: x(:)
+  real(psb_dpk_),intent(in)         :: x(:)
-  real(kind(0.d0)),intent(inout)      :: y(:)
+  real(psb_dpk_),intent(inout)      :: y(:)
  character, intent(in)               :: trans
-  real(kind(0.d0)),target             :: work(:)
+  real(psb_dpk_),target             :: work(:)
  integer, intent(out)                :: info
  ! Local variables
@ -366,11 +366,11 @@ contains
    ! Arguments
    type(psb_desc_type),intent(in)      :: desc_data
    type(mld_dbaseprc_type), intent(in) :: baseprecv(:)
-    real(kind(0.d0)),intent(in)         :: alpha,beta
+    real(psb_dpk_),intent(in)         :: alpha,beta
-    real(kind(0.d0)),intent(in)         :: x(:)
+    real(psb_dpk_),intent(in)         :: x(:)
-    real(kind(0.d0)),intent(inout)      :: y(:)
+    real(psb_dpk_),intent(inout)      :: y(:)
    character, intent(in)               :: trans
-    real(kind(0.d0)),target             :: work(:)
+    real(psb_dpk_),target             :: work(:)
    integer, intent(out)                :: info
    ! Local variables
@ -381,7 +381,7 @@ contains
    character(len=20)  :: name
    type psb_mlprec_wrk_type
-      real(kind(1.d0)), allocatable  :: tx(:), ty(:), x2l(:), y2l(:)
+      real(psb_dpk_), allocatable  :: tx(:), ty(:), x2l(:), y2l(:)
    end type psb_mlprec_wrk_type
    type(psb_mlprec_wrk_type), allocatable  :: mlprec_wrk(:)
@ -414,7 +414,7 @@ contains
    if (info /= 0) then 
      info=4025
      call psb_errpush(info,name,i_err=(/size(x)+size(y),0,0,0,0/),&
-           & a_err='real(kind(1.d0))')
+           & a_err='real(psb_dpk_)')
      goto 9999      
    end if
@ -442,7 +442,7 @@ contains
      if (info /= 0) then 
        info=4025
        call psb_errpush(info,name,i_err=(/2*(nc2l+max(n_row,n_col)),0,0,0,0/),&
-             & a_err='real(kind(1.d0))')
+             & a_err='real(psb_dpk_)')
        goto 9999      
      end if
@ -622,11 +622,11 @@ contains
    ! Arguments
    type(psb_desc_type),intent(in)      :: desc_data
    type(mld_dbaseprc_type), intent(in) :: baseprecv(:)
-    real(kind(0.d0)),intent(in)         :: alpha,beta
+    real(psb_dpk_),intent(in)         :: alpha,beta
-    real(kind(0.d0)),intent(in)         :: x(:)
+    real(psb_dpk_),intent(in)         :: x(:)
-    real(kind(0.d0)),intent(inout)      :: y(:)
+    real(psb_dpk_),intent(inout)      :: y(:)
    character, intent(in)               :: trans
-    real(kind(0.d0)),target             :: work(:)
+    real(psb_dpk_),target             :: work(:)
    integer, intent(out)                :: info
    ! Local variables
@ -637,7 +637,7 @@ contains
    character(len=20)  :: name
    type psb_mlprec_wrk_type
-      real(kind(1.d0)), allocatable  :: tx(:), ty(:), x2l(:), y2l(:)
+      real(psb_dpk_), allocatable  :: tx(:), ty(:), x2l(:), y2l(:)
    end type psb_mlprec_wrk_type
    type(psb_mlprec_wrk_type), allocatable  :: mlprec_wrk(:)
@ -674,7 +674,7 @@ contains
    if (info /= 0) then 
      info=4025
      call psb_errpush(info,name,i_err=(/4*nc2l,0,0,0,0/),&
-           & a_err='real(kind(1.d0))')
+           & a_err='real(psb_dpk_)')
      goto 9999
    end if
@ -726,7 +726,7 @@ contains
      if (info /= 0) then 
        info=4025
        call psb_errpush(info,name,i_err=(/4*nc2l,0,0,0,0/),&
-             & a_err='real(kind(1.d0))')
+             & a_err='real(psb_dpk_)')
        goto 9999      
      end if
@ -900,11 +900,11 @@ contains
    ! Arguments
    type(psb_desc_type),intent(in)      :: desc_data
    type(mld_dbaseprc_type), intent(in) :: baseprecv(:)
-    real(kind(0.d0)),intent(in)         :: alpha,beta
+    real(psb_dpk_),intent(in)         :: alpha,beta
-    real(kind(0.d0)),intent(in)         :: x(:)
+    real(psb_dpk_),intent(in)         :: x(:)
-    real(kind(0.d0)),intent(inout)      :: y(:)
+    real(psb_dpk_),intent(inout)      :: y(:)
    character, intent(in)               :: trans
-    real(kind(0.d0)),target             :: work(:)
+    real(psb_dpk_),target             :: work(:)
    integer, intent(out)                :: info
    ! Local variables
@ -915,7 +915,7 @@ contains
    character(len=20)  :: name
    type psb_mlprec_wrk_type
-      real(kind(1.d0)), allocatable  :: tx(:), ty(:), x2l(:), y2l(:)
+      real(psb_dpk_), allocatable  :: tx(:), ty(:), x2l(:), y2l(:)
    end type psb_mlprec_wrk_type
    type(psb_mlprec_wrk_type), allocatable  :: mlprec_wrk(:)
@ -985,7 +985,7 @@ contains
      if (info /= 0) then 
        info=4025
        call psb_errpush(info,name,i_err=(/4*nc2l,0,0,0,0/),&
-             & a_err='real(kind(1.d0))')
+             & a_err='real(psb_dpk_)')
        goto 9999      
      end if
@ -1210,11 +1210,11 @@ contains
    ! Arguments
    type(psb_desc_type),intent(in)      :: desc_data
    type(mld_dbaseprc_type), intent(in) :: baseprecv(:)
-    real(kind(0.d0)),intent(in)         :: alpha,beta
+    real(psb_dpk_),intent(in)         :: alpha,beta
-    real(kind(0.d0)),intent(in)         :: x(:)
+    real(psb_dpk_),intent(in)         :: x(:)
-    real(kind(0.d0)),intent(inout)      :: y(:)
+    real(psb_dpk_),intent(inout)      :: y(:)
    character, intent(in)               :: trans
-    real(kind(0.d0)),target             :: work(:)
+    real(psb_dpk_),target             :: work(:)
    integer, intent(out)                :: info
    ! Local variables
@ -1225,7 +1225,7 @@ contains
    character(len=20)  :: name
    type psb_mlprec_wrk_type
-      real(kind(1.d0)), allocatable  :: tx(:), ty(:), x2l(:), y2l(:)
+      real(psb_dpk_), allocatable  :: tx(:), ty(:), x2l(:), y2l(:)
    end type psb_mlprec_wrk_type
    type(psb_mlprec_wrk_type), allocatable  :: mlprec_wrk(:)
@ -1262,7 +1262,7 @@ contains
    if (info /= 0) then 
      info=4025
      call psb_errpush(info,name,i_err=(/4*nc2l,0,0,0,0/),&
-           & a_err='real(kind(1.d0))')
+           & a_err='real(psb_dpk_)')
      goto 9999
    end if
@ -1312,7 +1312,7 @@ contains
      if (info /= 0) then 
        info=4025
        call psb_errpush(info,name,i_err=(/4*nc2l,0,0,0,0/),&
-             & a_err='real(kind(1.d0))')
+             & a_err='real(psb_dpk_)')
        goto 9999      
      end if
--- a/mlprec/mld_dprec_aply.f90
+++ b/mlprec/mld_dprec_aply.f90
@ -55,9 +55,9 @@
 !    prec       -  type(mld_dprec_type), input.
 !                  The preconditioner data structure containing the local part
 !                  of the preconditioner to be applied.
-!    x          -  real(kind(0.d0)), dimension(:), input.
+!    x          -  real(psb_dpk_), dimension(:), input.
 !                  The local part of the vector X in Y := op(M^(-1)) * X.
-!    y          -  real(kind(0.d0)), dimension(:), output.
+!    y          -  real(psb_dpk_), dimension(:), output.
 !                  The local part of the vector Y in Y := op(M^(-1)) * X.
 !    desc_data  -  type(psb_desc_type), input.
 !                  The communication descriptor associated to the matrix to be
@ -67,7 +67,7 @@
 !    trans      -  character(len=1), optional.
 !                  If trans='N','n' then op(M^(-1)) = M^(-1);
 !                  if trans='T','t' then op(M^(-1)) = M^(-T) (transpose of M^(-1)).
-!    work       -  real(kind(0.d0)), dimension (:), optional, target.
+!    work       -  real(psb_dpk_), dimension (:), optional, target.
 !                  Workspace. Its size must be at
 !                  least 4*psb_cd_get_local_cols(desc_data).
 !    
@ -82,15 +82,15 @@ subroutine mld_dprec_aply(prec,x,y,desc_data,info,trans,work)
  ! Arguments
  type(psb_desc_type),intent(in)      :: desc_data
  type(mld_dprec_type), intent(in)    :: prec
-  real(kind(0.d0)),intent(in)         :: x(:)
+  real(psb_dpk_),intent(in)         :: x(:)
-  real(kind(0.d0)),intent(inout)      :: y(:)
+  real(psb_dpk_),intent(inout)      :: y(:)
  integer, intent(out)                :: info
  character(len=1), optional          :: trans
-  real(kind(0.d0)), optional, target  :: work(:)
+  real(psb_dpk_), optional, target  :: work(:)
  ! Local variables
  character     :: trans_ 
-  real(kind(1.d0)), pointer :: work_(:)
+  real(psb_dpk_), pointer :: work_(:)
  integer :: ictxt,np,me,err_act,iwsz
  character(len=20)   :: name
@ -114,7 +114,7 @@ subroutine mld_dprec_aply(prec,x,y,desc_data,info,trans,work)
    allocate(work_(iwsz),stat=info)
    if (info /= 0) then 
      call psb_errpush(4025,name,i_err=(/iwsz,0,0,0,0/),&
-           &a_err='real(kind(1.d0))')
+           &a_err='real(psb_dpk_)')
      goto 9999      
    end if
@ -186,7 +186,7 @@ end subroutine mld_dprec_aply
 !    prec       -  type(mld_dprec_type), input.
 !                  The preconditioner data structure containing the local part
 !                  of the preconditioner to be applied.
-!    x          -  real(kind(0.d0)), dimension(:), input/output.
+!    x          -  real(psb_dpk_), dimension(:), input/output.
 !                  The local part of vector X in X := op(M^(-1)) * X.
 !    desc_data  -  type(psb_desc_type), input.
 !                  The communication descriptor associated to the matrix to be
@ -208,13 +208,13 @@ subroutine mld_dprec_aply1(prec,x,desc_data,info,trans)
  ! Arguments
  type(psb_desc_type),intent(in)    :: desc_data
  type(mld_dprec_type), intent(in)  :: prec
-  real(kind(0.d0)),intent(inout)    :: x(:)
+  real(psb_dpk_),intent(inout)    :: x(:)
  integer, intent(out)              :: info
  character(len=1), optional        :: trans
  ! Local variables
  integer       :: ictxt,np,me, err_act
-  real(kind(1.d0)), pointer :: WW(:), w1(:)
+  real(psb_dpk_), pointer :: WW(:), w1(:)
  character(len=20)   :: name
  name='mld_dprec_aply1'
@ -229,7 +229,7 @@ subroutine mld_dprec_aply1(prec,x,desc_data,info,trans)
  if (info /= 0) then 
    info=4025
    call psb_errpush(info,name,i_err=(/2*size(x),0,0,0,0/),&
-         & a_err='real(kind(1.d0))')
+         & a_err='real(psb_dpk_)')
    goto 9999      
  end if
--- a/mlprec/mld_dprecset.f90
+++ b/mlprec/mld_dprecset.f90
@ -519,7 +519,7 @@ end subroutine mld_dprecsetc
 !               The number identifying the parameter to be set.
 !               A mnemonic constant has been associated to each of these
 !               numbers, as reported in MLD2P4 user's guide.
-!    val     -  real(kind(1.d0)), input.
+!    val     -  real(psb_dpk_), input.
 !               The value of the parameter to be set. The list of allowed
 !               values is reported in MLD2P4 user's guide.
 !    info    -  integer, output.
@ -540,7 +540,7 @@ subroutine mld_dprecsetd(p,what,val,info,ilev)
  ! Arguments
  type(mld_dprec_type), intent(inout)    :: p
  integer, intent(in)                    :: what 
-  real(kind(1.d0)), intent(in)           :: val
+  real(psb_dpk_), intent(in)           :: val
  integer, intent(out)                   :: info
  integer, optional, intent(in)          :: ilev
--- a/mlprec/mld_dsp_renum.f90
+++ b/mlprec/mld_dsp_renum.f90
@ -99,7 +99,7 @@ subroutine mld_dsp_renum(a,blck,p,atmp,info)
  integer              :: nztota, nztotb, nztmp, nnr, i,k
  integer, allocatable :: itmp(:), itmp2(:)
  integer              :: ictxt,np,me, err_act
-  real(kind(1.d0))     :: t3,t4
+  real(psb_dpk_)     :: t3,t4
  if (psb_get_errstatus().ne.0) return 
  info=0
--- a/mlprec/mld_dsub_aply.f90
+++ b/mlprec/mld_dsub_aply.f90
@ -108,16 +108,16 @@
 !
 ! Arguments:
 !
-!   alpha      -  real(kind(0.d0)), input.
+!   alpha      -  real(psb_dpk_), input.
 !                 The scalar alpha.
 !   prec       -  type(mld_dbaseprec_type), input.
 !                 The 'base preconditioner' data structure containing the local 
 !                 part of the preconditioner or solver.
-!   x          -  real(kind(0.d0)), dimension(:), input.
+!   x          -  real(psb_dpk_), dimension(:), input.
 !                 The local part of the vector X.
-!   beta       -  real(kind(0.d0)), input.
+!   beta       -  real(psb_dpk_), input.
 !                 The scalar beta.
-!   y          -  real(kind(0.d0)), dimension(:), input/output.
+!   y          -  real(psb_dpk_), dimension(:), input/output.
 !                 The local part of the vector Y.
 !   desc_data  -  type(psb_desc_type), input.
 !                 The communication descriptor associated to the matrix to be
@ -127,7 +127,7 @@
 !                 if trans='T','t' then op(K^(-1)) = K^(-T) (transpose of K^(-1)).
 !                 If prec%iprcparm(smooth_sweeps_) > 1, the value of trans provided
 !                 in input is ignored.
-!   work       -  real(kind(0.d0)), dimension (:), target.
+!   work       -  real(psb_dpk_), dimension (:), target.
 !                 Workspace. Its size must be at least 4*psb_cd_get_local_cols(desc_data).
 !   info       -  integer, output.
 !                 Error code.
@ -142,16 +142,16 @@ subroutine mld_dsub_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
  ! Arguments
  type(psb_desc_type), intent(in)        :: desc_data
  type(mld_dbaseprc_type), intent(in)    :: prec
-  real(kind(0.d0)),intent(in)            :: x(:)
+  real(psb_dpk_),intent(in)            :: x(:)
-  real(kind(0.d0)),intent(inout)         :: y(:)
+  real(psb_dpk_),intent(inout)         :: y(:)
-  real(kind(0.d0)),intent(in)            :: alpha,beta
+  real(psb_dpk_),intent(in)            :: alpha,beta
  character(len=1),intent(in)            :: trans
-  real(kind(0.d0)),target, intent(inout) :: work(:)
+  real(psb_dpk_),target, intent(inout) :: work(:)
  integer, intent(out)                   :: info
  ! Local variables
  integer :: n_row,n_col
-  real(kind(1.d0)), pointer :: ww(:), aux(:), tx(:),ty(:)
+  real(psb_dpk_), pointer :: ww(:), aux(:), tx(:),ty(:)
  integer :: ictxt,np,me,i, err_act
  character(len=20)   :: name
  character           :: trans_
@ -185,7 +185,7 @@ subroutine mld_dsub_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
      if (info /= 0) then 
        info=4025
        call psb_errpush(info,name,i_err=(/4*n_col,0,0,0,0/),&
-             & a_err='real(kind(1.d0))')
+             & a_err='real(psb_dpk_)')
        goto 9999      
      end if
    endif
@ -194,7 +194,7 @@ subroutine mld_dsub_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
    if (info /= 0) then 
      info=4025
      call psb_errpush(info,name,i_err=(/5*n_col,0,0,0,0/),&
-           & a_err='real(kind(1.d0))')
+           & a_err='real(psb_dpk_)')
      goto 9999      
    end if
  endif
@ -225,7 +225,7 @@ subroutine mld_dsub_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
    if (info /= 0) then 
      info=4025
      call psb_errpush(info,name,i_err=(/2*n_col,0,0,0,0/),&
-           & a_err='real(kind(1.d0))')
+           & a_err='real(psb_dpk_)')
      goto 9999      
    end if
--- a/mlprec/mld_dsub_solve.f90
+++ b/mlprec/mld_dsub_solve.f90
@ -86,16 +86,16 @@
 !
 ! Arguments:
 !
-!   alpha      -  real(kind(0.d0)), input.
+!   alpha      -  real(psb_dpk_), input.
 !                 The scalar alpha.
 !   prec       -  type(mld_dbaseprec_type), input.
 !                 The 'base preconditioner' data structure containing the local 
 !                 part of the L and U factors of the matrix A.
-!   x          -  real(kind(0.d0)), dimension(:), input.
+!   x          -  real(psb_dpk_), dimension(:), input.
 !                 The local part of the vector X.
-!   beta       -  real(kind(0.d0)), input.
+!   beta       -  real(psb_dpk_), input.
 !                 The scalar beta.
-!   y          -  real(kind(0.d0)), dimension(:), input/output.
+!   y          -  real(psb_dpk_), dimension(:), input/output.
 !                 The local part of the vector Y.
 !   desc_data  -  type(psb_desc_type), input.
 !                 The communication descriptor associated to the matrix to be
@ -105,7 +105,7 @@
 !                 if trans='T','t' then op(K^(-1)) = K^(-T) (transpose of K^(-1)).
 !                 If prec%iprcparm(smooth_sweeps_) > 1, the value of trans provided
 !                 in input is ignored.
-!   work       -  real(kind(0.d0)), dimension (:), target.
+!   work       -  real(psb_dpk_), dimension (:), target.
 !                 Workspace. Its size must be at least 4*psb_cd_get_local_cols(desc_data).
 !   info       -  integer, output.
 !                 Error code.
@ -120,26 +120,27 @@ subroutine mld_dsub_solve(alpha,prec,x,beta,y,desc_data,trans,work,info)
  ! Arguments
  type(psb_desc_type), intent(in)        :: desc_data
  type(mld_dbaseprc_type), intent(in)    :: prec
-  real(kind(0.d0)),intent(in)            :: x(:)
+  real(psb_dpk_),intent(in)            :: x(:)
-  real(kind(0.d0)),intent(inout)         :: y(:)
+  real(psb_dpk_),intent(inout)         :: y(:)
-  real(kind(0.d0)),intent(in)            :: alpha,beta
+  real(psb_dpk_),intent(in)            :: alpha,beta
  character(len=1),intent(in)            :: trans
-  real(kind(0.d0)),target, intent(inout) :: work(:)
+  real(psb_dpk_),target, intent(inout) :: work(:)
  integer, intent(out)                   :: info
  ! Local variables
  integer :: n_row,n_col
-  real(kind(1.d0)), pointer :: ww(:), aux(:), tx(:),ty(:)
+  real(psb_dpk_), pointer :: ww(:), aux(:), tx(:),ty(:)
  integer :: ictxt,np,me,i, err_act
  character(len=20)   :: name
  character           :: trans_
  interface 
    subroutine mld_dumf_solve(flag,m,x,b,n,ptr,info)
      use psb_base_mod
      integer, intent(in)  :: flag,m,n,ptr
      integer, intent(out) :: info
-      real(kind(1.d0)), intent(in)    :: b(*)
+      real(psb_dpk_), intent(in)    :: b(*)
-      real(kind(1.d0)), intent(inout) :: x(*)
+      real(psb_dpk_), intent(inout) :: x(*)
    end subroutine mld_dumf_solve
  end interface
@ -172,7 +173,7 @@ subroutine mld_dsub_solve(alpha,prec,x,beta,y,desc_data,trans,work,info)
      if (info /= 0) then 
        info=4025
        call psb_errpush(info,name,i_err=(/4*n_col,0,0,0,0/),&
-             & a_err='real(kind(1.d0))')
+             & a_err='real(psb_dpk_)')
        goto 9999      
      end if
    endif
@ -181,7 +182,7 @@ subroutine mld_dsub_solve(alpha,prec,x,beta,y,desc_data,trans,work,info)
    if (info /= 0) then 
      info=4025
      call psb_errpush(info,name,i_err=(/5*n_col,0,0,0,0/),&
-           & a_err='real(kind(1.d0))')
+           & a_err='real(psb_dpk_)')
      goto 9999      
    end if
  endif
--- a/mlprec/mld_inner_mod.f90
+++ b/mlprec/mld_inner_mod.f90
@ -45,11 +45,11 @@ module mld_inner_mod
      use mld_prec_type
      type(psb_desc_type),intent(in)      :: desc_data
      type(mld_dbaseprc_type), intent(in) :: prec
-      real(kind(0.d0)),intent(in)         :: x(:)
+      real(psb_dpk_),intent(in)         :: x(:)
-      real(kind(0.d0)),intent(inout)      :: y(:)
+      real(psb_dpk_),intent(inout)      :: y(:)
-      real(kind(0.d0)),intent(in)         :: alpha,beta
+      real(psb_dpk_),intent(in)         :: alpha,beta
      character(len=1)                    :: trans
-      real(kind(0.d0)),target             :: work(:)
+      real(psb_dpk_),target             :: work(:)
      integer, intent(out)                :: info
    end subroutine mld_dbaseprec_aply
    subroutine mld_zbaseprec_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
@ -57,11 +57,11 @@ module mld_inner_mod
      use mld_prec_type
      type(psb_desc_type),intent(in)      :: desc_data
      type(mld_zbaseprc_type), intent(in) :: prec
-      complex(kind(1.d0)),intent(in)      :: x(:)
+      complex(psb_dpk_),intent(in)      :: x(:)
-      complex(kind(1.d0)),intent(inout)   :: y(:)
+      complex(psb_dpk_),intent(inout)   :: y(:)
-      complex(kind(1.d0)),intent(in)      :: alpha,beta
+      complex(psb_dpk_),intent(in)      :: alpha,beta
      character(len=1)                    :: trans
-      complex(kind(1.d0)),target          :: work(:)
+      complex(psb_dpk_),target          :: work(:)
      integer, intent(out)                :: info
    end subroutine mld_zbaseprec_aply
  end interface
@ -72,11 +72,11 @@ module mld_inner_mod
      use mld_prec_type
      type(psb_desc_type),intent(in)      :: desc_data
      type(mld_dbaseprc_type), intent(in) :: prec
-      real(kind(0.d0)),intent(in)         :: x(:)
+      real(psb_dpk_),intent(in)         :: x(:)
-      real(kind(0.d0)),intent(inout)      :: y(:)
+      real(psb_dpk_),intent(inout)      :: y(:)
-      real(kind(0.d0)),intent(in)         :: alpha,beta
+      real(psb_dpk_),intent(in)         :: alpha,beta
      character(len=1)                    :: trans
-      real(kind(0.d0)),target             :: work(:)
+      real(psb_dpk_),target             :: work(:)
      integer, intent(out)                :: info
    end subroutine mld_das_aply
    subroutine mld_zas_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
@ -84,11 +84,11 @@ module mld_inner_mod
      use mld_prec_type
      type(psb_desc_type),intent(in)      :: desc_data
      type(mld_zbaseprc_type), intent(in) :: prec
-      complex(kind(1.d0)),intent(in)      :: x(:)
+      complex(psb_dpk_),intent(in)      :: x(:)
-      complex(kind(1.d0)),intent(inout)   :: y(:)
+      complex(psb_dpk_),intent(inout)   :: y(:)
-      complex(kind(1.d0)),intent(in)      :: alpha,beta
+      complex(psb_dpk_),intent(in)      :: alpha,beta
      character(len=1)                    :: trans
-      complex(kind(1.d0)),target          :: work(:)
+      complex(psb_dpk_),target          :: work(:)
      integer, intent(out)                :: info
    end subroutine mld_zas_aply
  end interface
@ -99,11 +99,11 @@ module mld_inner_mod
      use mld_prec_type
      type(psb_desc_type),intent(in)      :: desc_data
      type(mld_dbaseprc_type), intent(in) :: baseprecv(:)
-      real(kind(0.d0)),intent(in)         :: alpha,beta
+      real(psb_dpk_),intent(in)         :: alpha,beta
-      real(kind(0.d0)),intent(in)         :: x(:)
+      real(psb_dpk_),intent(in)         :: x(:)
-      real(kind(0.d0)),intent(inout)      :: y(:)
+      real(psb_dpk_),intent(inout)      :: y(:)
      character                           :: trans
-      real(kind(0.d0)),target             :: work(:)
+      real(psb_dpk_),target             :: work(:)
      integer, intent(out)                :: info
    end subroutine mld_dmlprec_aply
    subroutine mld_zmlprec_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
@ -111,11 +111,11 @@ module mld_inner_mod
      use mld_prec_type
      type(psb_desc_type),intent(in)      :: desc_data
      type(mld_zbaseprc_type), intent(in) :: baseprecv(:)
-      complex(kind(0.d0)),intent(in)      :: alpha,beta
+      complex(psb_dpk_),intent(in)      :: alpha,beta
-      complex(kind(0.d0)),intent(in)      :: x(:)
+      complex(psb_dpk_),intent(in)      :: x(:)
-      complex(kind(0.d0)),intent(inout)   :: y(:)
+      complex(psb_dpk_),intent(inout)   :: y(:)
      character                           :: trans
-      complex(kind(0.d0)),target          :: work(:)
+      complex(psb_dpk_),target          :: work(:)
      integer, intent(out)                :: info
    end subroutine mld_zmlprec_aply
  end interface
@ -127,11 +127,11 @@ module mld_inner_mod
      use mld_prec_type
      type(psb_desc_type), intent(in)       :: desc_data
      type(mld_dbaseprc_type), intent(in)   :: prec
-      real(kind(0.d0)),intent(in)           :: x(:)
+      real(psb_dpk_),intent(in)           :: x(:)
-      real(kind(0.d0)),intent(inout)        :: y(:)
+      real(psb_dpk_),intent(inout)        :: y(:)
-      real(kind(0.d0)),intent(in)           :: alpha,beta
+      real(psb_dpk_),intent(in)           :: alpha,beta
      character(len=1),intent(in)           :: trans
-      real(kind(0.d0)),target,intent(inout) :: work(:)
+      real(psb_dpk_),target,intent(inout) :: work(:)
      integer, intent(out)                  :: info
    end subroutine mld_dsub_aply
    subroutine mld_zsub_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
@ -139,11 +139,11 @@ module mld_inner_mod
      use mld_prec_type
      type(psb_desc_type), intent(in)          :: desc_data
      type(mld_zbaseprc_type), intent(in)      :: prec
-      complex(kind(0.d0)),intent(in)           :: x(:)
+      complex(psb_dpk_),intent(in)           :: x(:)
-      complex(kind(0.d0)),intent(inout)        :: y(:)
+      complex(psb_dpk_),intent(inout)        :: y(:)
-      complex(kind(0.d0)),intent(in)           :: alpha,beta
+      complex(psb_dpk_),intent(in)           :: alpha,beta
      character(len=1),intent(in)              :: trans
-      complex(kind(0.d0)),target,intent(inout) :: work(:)
+      complex(psb_dpk_),target,intent(inout) :: work(:)
      integer, intent(out)                     :: info
    end subroutine mld_zsub_aply
  end interface
@ -155,11 +155,11 @@ module mld_inner_mod
      use mld_prec_type
      type(psb_desc_type), intent(in)       :: desc_data
      type(mld_dbaseprc_type), intent(in)   :: prec
-      real(kind(0.d0)),intent(in)           :: x(:)
+      real(psb_dpk_),intent(in)           :: x(:)
-      real(kind(0.d0)),intent(inout)        :: y(:)
+      real(psb_dpk_),intent(inout)        :: y(:)
-      real(kind(0.d0)),intent(in)           :: alpha,beta
+      real(psb_dpk_),intent(in)           :: alpha,beta
      character(len=1),intent(in)           :: trans
-      real(kind(0.d0)),target,intent(inout) :: work(:)
+      real(psb_dpk_),target,intent(inout) :: work(:)
      integer, intent(out)                  :: info
    end subroutine mld_dsub_solve
    subroutine mld_zsub_solve(alpha,prec,x,beta,y,desc_data,trans,work,info)
@ -167,11 +167,11 @@ module mld_inner_mod
      use mld_prec_type
      type(psb_desc_type), intent(in)          :: desc_data
      type(mld_zbaseprc_type), intent(in)      :: prec
-      complex(kind(0.d0)),intent(in)           :: x(:)
+      complex(psb_dpk_),intent(in)           :: x(:)
-      complex(kind(0.d0)),intent(inout)        :: y(:)
+      complex(psb_dpk_),intent(inout)        :: y(:)
-      complex(kind(0.d0)),intent(in)           :: alpha,beta
+      complex(psb_dpk_),intent(in)           :: alpha,beta
      character(len=1),intent(in)              :: trans
-      complex(kind(0.d0)),target,intent(inout) :: work(:)
+      complex(psb_dpk_),target,intent(inout) :: work(:)
      integer, intent(out)                     :: info
    end subroutine mld_zsub_solve
  end interface
--- a/mlprec/mld_prec_mod.f90
+++ b/mlprec/mld_prec_mod.f90
@ -81,7 +81,7 @@ module mld_prec_mod
      use mld_prec_type
      type(mld_dprec_type), intent(inout)    :: p
      integer, intent(in)                    :: what 
-      real(kind(1.d0)), intent(in)           :: val
+      real(psb_dpk_), intent(in)           :: val
      integer, intent(out)                   :: info
      integer, optional, intent(in)          :: ilev
    end subroutine mld_dprecsetd
@ -108,7 +108,7 @@ module mld_prec_mod
      use mld_prec_type
      type(mld_zprec_type), intent(inout)    :: p
      integer, intent(in)                    :: what 
-      real(kind(1.d0)), intent(in)           :: val
+      real(psb_dpk_), intent(in)           :: val
      integer, intent(out)                   :: info
      integer, optional, intent(in)          :: ilev
    end subroutine mld_zprecsetd
@ -144,18 +144,18 @@ module mld_prec_mod
      use mld_prec_type
      type(psb_desc_type),intent(in)    :: desc_data
      type(mld_dprec_type), intent(in)  :: prec
-      real(kind(0.d0)),intent(in)       :: x(:)
+      real(psb_dpk_),intent(in)       :: x(:)
-      real(kind(0.d0)),intent(inout)    :: y(:)
+      real(psb_dpk_),intent(inout)    :: y(:)
      integer, intent(out)              :: info
      character(len=1), optional        :: trans
-      real(kind(0.d0)),intent(inout), optional, target :: work(:)
+      real(psb_dpk_),intent(inout), optional, target :: work(:)
    end subroutine mld_dprec_aply
    subroutine mld_dprec_aply1(prec,x,desc_data,info,trans)
      use psb_base_mod
      use mld_prec_type
      type(psb_desc_type),intent(in)    :: desc_data
      type(mld_dprec_type), intent(in)  :: prec
-      real(kind(0.d0)),intent(inout)    :: x(:)
+      real(psb_dpk_),intent(inout)    :: x(:)
      integer, intent(out)              :: info
      character(len=1), optional        :: trans
    end subroutine mld_dprec_aply1
@ -164,18 +164,18 @@ module mld_prec_mod
      use mld_prec_type
      type(psb_desc_type),intent(in)    :: desc_data
      type(mld_zprec_type), intent(in)  :: prec
-      complex(kind(0.d0)),intent(in)    :: x(:)
+      complex(psb_dpk_),intent(in)    :: x(:)
-      complex(kind(0.d0)),intent(inout) :: y(:)
+      complex(psb_dpk_),intent(inout) :: y(:)
      integer, intent(out)              :: info
      character(len=1), optional        :: trans
-      complex(kind(0.d0)),intent(inout), optional, target :: work(:)
+      complex(psb_dpk_),intent(inout), optional, target :: work(:)
    end subroutine mld_zprec_aply
    subroutine mld_zprec_aply1(prec,x,desc_data,info,trans)
      use psb_base_mod
      use mld_prec_type
      type(psb_desc_type),intent(in)    :: desc_data
      type(mld_zprec_type), intent(in)  :: prec
-      complex(kind(0.d0)),intent(inout) :: x(:)
+      complex(psb_dpk_),intent(inout) :: x(:)
      integer, intent(out)              :: info
      character(len=1), optional        :: trans
    end subroutine mld_zprec_aply1
--- a/mlprec/mld_prec_type.f90
+++ b/mlprec/mld_prec_type.f90
@ -63,7 +63,7 @@ module mld_prec_type
  ! blows up on some systems.
  !
  use psb_base_mod, only : psb_dspmat_type, psb_zspmat_type, psb_desc_type,&
-       & psb_inter_desc_type, psb_sizeof
+       & psb_inter_desc_type, psb_sizeof, psb_dpk_
  !
  ! Type: mld_dprec_type, mld_zprec_type
@ -113,7 +113,7 @@ module mld_prec_type
  !      Shouldn't we keep just one of the last two items and handle the transpose
  !      in the Sparse BLAS? Maybe.
  !
-  !   d            -  real(kind(1.d0)), dimension(:), allocatable.
+  !   d            -  real(psb_dpk_), dimension(:), allocatable.
  !                   The diagonal entries of the U factor in the ILU factorization
  !                   of A(ilev).
  !   desc_data    -  type(psb_desc_type).
@ -126,7 +126,7 @@ module mld_prec_type
  !   iprcparm     -  integer, dimension(:), allocatable.
  !                   The integer parameters defining the base preconditioner K(ilev)
  !                   (the iprcparm entries and values are specified below).
-  !   dprcparm     -  real(kind(1.d0)), dimension(:), allocatable.
+  !   dprcparm     -  real(psb_dpk_), dimension(:), allocatable.
  !                   The real parameters defining the base preconditioner K(ilev)
  !                   (the dprcparm entries and values are specified below).
  !   perm         -  integer, dimension(:), allocatable.
@ -149,7 +149,7 @@ module mld_prec_type
  !   base_desc    -  type(psb_desc_type), pointer.
  !                   Pointer to the communication descriptor associated to the sparse
  !                   matrix pointed by base_a. 
-  !   dorig        -  real(kind(1.d0)), dimension(:), allocatable.
+  !   dorig        -  real(psb_dpk_), dimension(:), allocatable.
  !                   Diagonal entries of the matrix pointed by base_a.
  !
  !   Note that when the LU factorization of the matrix A(ilev) is computed instead of
@ -161,15 +161,15 @@ module mld_prec_type
  type mld_dbaseprc_type
    type(psb_dspmat_type), allocatable :: av(:) 
-    real(kind(1.d0)), allocatable      :: d(:)  
+    real(psb_dpk_), allocatable      :: d(:)  
    type(psb_desc_type)                :: desc_data , desc_ac
    integer, allocatable               :: iprcparm(:) 
-    real(kind(1.d0)), allocatable      :: dprcparm(:) 
+    real(psb_dpk_), allocatable      :: dprcparm(:) 
    integer, allocatable               :: perm(:),  invperm(:) 
    integer, allocatable               :: mlia(:), nlaggr(:) 
    type(psb_dspmat_type), pointer     :: base_a    => null() 
    type(psb_desc_type), pointer       :: base_desc => null() 
-    real(kind(1.d0)), allocatable      :: dorig(:) 
+    real(psb_dpk_), allocatable      :: dorig(:) 
    type(psb_inter_desc_type)          :: map_desc
  end type mld_dbaseprc_type
@ -180,15 +180,15 @@ module mld_prec_type
  type mld_zbaseprc_type
    type(psb_zspmat_type), allocatable :: av(:) 
-    complex(kind(1.d0)), allocatable   :: d(:)  
+    complex(psb_dpk_), allocatable   :: d(:)  
    type(psb_desc_type)                :: desc_data , desc_ac
    integer, allocatable               :: iprcparm(:) 
-    real(kind(1.d0)), allocatable      :: dprcparm(:) 
+    real(psb_dpk_), allocatable      :: dprcparm(:) 
    integer, allocatable               :: perm(:),  invperm(:) 
    integer, allocatable               :: mlia(:), nlaggr(:) 
    type(psb_zspmat_type), pointer     :: base_a    => null() 
    type(psb_desc_type), pointer       :: base_desc => null() 
-    complex(kind(1.d0)), allocatable   :: dorig(:)
+    complex(psb_dpk_), allocatable   :: dorig(:)
    type(psb_inter_desc_type)          :: map_desc
  end type mld_zbaseprc_type
@ -409,10 +409,10 @@ contains
      end if
    end if
-    if (allocated(prec%dprcparm)) val = val + 8 * size(prec%dprcparm)
+    if (allocated(prec%dprcparm)) val = val + psb_sizeof_dp * size(prec%dprcparm)
-    if (allocated(prec%d))        val = val + 8 * size(prec%d)
+    if (allocated(prec%d))        val = val + psb_sizeof_dp * size(prec%d)
-    if (allocated(prec%perm))     val = val + 4 * size(prec%perm)
+    if (allocated(prec%perm))     val = val + psb_sizeof_int * size(prec%perm)
-    if (allocated(prec%invperm))  val = val + 4 * size(prec%invperm)
+    if (allocated(prec%invperm))  val = val + psb_sizeof_int * size(prec%invperm)
                                  val = val + psb_sizeof(prec%desc_data)
    if (allocated(prec%av))  then 
      do i=1,size(prec%av)
@ -446,10 +446,10 @@ contains
      end if
    end if
-    if (allocated(prec%dprcparm)) val = val + 8 * size(prec%dprcparm)
+    if (allocated(prec%dprcparm)) val = val + psb_sizeof_dp * size(prec%dprcparm)
-    if (allocated(prec%d))        val = val + 16 * size(prec%d)
+    if (allocated(prec%d))        val = val + 2 * psb_sizeof_dp * size(prec%d)
-    if (allocated(prec%perm))     val = val + 4 * size(prec%perm)
+    if (allocated(prec%perm))     val = val + psb_sizeof_int * size(prec%perm)
-    if (allocated(prec%invperm))  val = val + 4 * size(prec%invperm)
+    if (allocated(prec%invperm))  val = val + psb_sizeof_int * size(prec%invperm)
                                  val = val + psb_sizeof(prec%desc_data)
    if (allocated(prec%av))  then 
      do i=1,size(prec%av)
@ -857,14 +857,14 @@ contains
  end function is_legal_ml_lev
  function is_legal_omega(ip)
    use psb_base_mod
-    real(kind(1.d0)), intent(in) :: ip
+    real(psb_dpk_), intent(in) :: ip
    logical             :: is_legal_omega
    is_legal_omega = ((ip>=0.0d0).and.(ip<=2.0d0))
    return
  end function is_legal_omega
  function is_legal_fact_thrs(ip)
    use psb_base_mod
-    real(kind(1.d0)), intent(in) :: ip
+    real(psb_dpk_), intent(in) :: ip
    logical             :: is_legal_fact_thrs
    is_legal_fact_thrs = (ip>=0.0d0)
@ -894,12 +894,13 @@ contains
  subroutine mld_dcheck_def(ip,name,id,is_legal)
    use psb_base_mod
-    real(kind(1.d0)), intent(inout) :: ip
+    real(psb_dpk_), intent(inout) :: ip
-    real(kind(1.d0)), intent(in)    :: id
+    real(psb_dpk_), intent(in)    :: id
    character(len=*), intent(in) :: name
    interface 
      function is_legal(i)
-        real(kind(1.d0)), intent(in) :: i
+        use psb_base_mod
        real(psb_dpk_), intent(in) :: i
        logical             :: is_legal
      end function is_legal
    end interface
--- a/mlprec/mld_zaggrmat_smth_asb.F90
+++ b/mlprec/mld_zaggrmat_smth_asb.F90
@ -126,7 +126,7 @@ subroutine mld_zaggrmat_smth_asb(a,desc_a,ac,desc_ac,p,info)
  logical            :: ml_global_nmb
  integer            :: debug_level, debug_unit
  integer, parameter :: ncmax=16
-  real(kind(1.d0))   :: omega, anorm, tmp, dg
+  real(psb_dpk_)   :: omega, anorm, tmp, dg
  name='mld_aggrmat_smth_asb'
  if(psb_get_errstatus().ne.0) return 
@ -190,7 +190,7 @@ subroutine mld_zaggrmat_smth_asb(a,desc_a,ac,desc_ac,p,info)
  if (info /= 0) then 
    info=4025
    call psb_errpush(info,name,i_err=(/nrow,0,0,0,0/),&
-         & a_err='real(kind(1.d0))')
+         & a_err='real(psb_dpk_)')
    goto 9999      
  end if
--- a/mlprec/mld_zas_aply.f90
+++ b/mlprec/mld_zas_aply.f90
@ -52,16 +52,16 @@
 !
 !
 ! Arguments:
-!   alpha      -  real(kind(0.d0)), input.
+!   alpha      -  real(psb_dpk_), input.
 !                 The scalar alpha.
 !   prec       -  type(mld_dbaseprc_type), input.
 !                 The base preconditioner data structure containing the local part
 !                 of the preconditioner K.
-!   x          -  real(kind(0.d0)), dimension(:), input.
+!   x          -  real(psb_dpk_), dimension(:), input.
 !                 The local part of the vector X.
-!   beta       -  real(kind(0.d0)), input.
+!   beta       -  real(psb_dpk_), input.
 !                 The scalar beta.
-!   y          -  real(kind(0.d0)), dimension(:), input/output.
+!   y          -  real(psb_dpk_), dimension(:), input/output.
 !                 The local part of the vector Y.
 !   desc_data  -  type(psb_desc_type), input.
 !                 The communication descriptor associated to the matrix to be
@ -69,7 +69,7 @@
 !   trans      -  character, optional.
 !                 If trans='N','n' then op(K^(-1)) = K^(-1);
 !                 if trans='T','t' then op(K^(-1)) = K^(-T) (transpose of K^(-1)).
-!   work       -  real(kind(0.d0)), dimension (:), optional, target.
+!   work       -  real(psb_dpk_), dimension (:), optional, target.
 !                 Workspace. Its size must be at least 4*psb_cd_get_local_cols(desc_data).
 !   info       -  integer, output.
 !                 Error code.
@ -84,16 +84,16 @@ subroutine mld_zas_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
  ! Arguments
  type(psb_desc_type),intent(in)      :: desc_data
  type(mld_zbaseprc_type), intent(in) :: prec
-  complex(kind(0.d0)),intent(in)         :: x(:)
+  complex(psb_dpk_),intent(in)         :: x(:)
-  complex(kind(0.d0)),intent(inout)      :: y(:)
+  complex(psb_dpk_),intent(inout)      :: y(:)
-  complex(kind(0.d0)),intent(in)         :: alpha,beta
+  complex(psb_dpk_),intent(in)         :: alpha,beta
  character(len=1)                    :: trans
-  complex(kind(0.d0)),target             :: work(:)
+  complex(psb_dpk_),target             :: work(:)
  integer, intent(out)                :: info
  ! Local variables
  integer :: n_row,n_col, int_err(5), nrow_d
-  complex(kind(1.d0)), pointer :: ww(:), aux(:), tx(:),ty(:)
+  complex(psb_dpk_), pointer :: ww(:), aux(:), tx(:),ty(:)
  integer           :: ictxt,np,me,isz, err_act
  character(len=20) :: name, ch_err
  character         :: trans_
@ -154,7 +154,7 @@ subroutine mld_zas_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
        allocate(ww(isz),tx(isz),ty(isz),stat=info)
        if (info /= 0) then 
          call psb_errpush(4025,name,i_err=(/3*isz,0,0,0,0/),&
-               & a_err='complex(kind(1.d0))')
+               & a_err='complex(psb_dpk_)')
          goto 9999      
        end if
      else if ((3*isz) <= size(work)) then 
@ -164,7 +164,7 @@ subroutine mld_zas_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
        allocate(aux(4*isz),stat=info)
        if (info /= 0) then 
          call psb_errpush(4025,name,i_err=(/4*isz,0,0,0,0/),&
-               & a_err='complex(kind(1.d0))')
+               & a_err='complex(psb_dpk_)')
          goto 9999      
        end if
      else 
@ -172,7 +172,7 @@ subroutine mld_zas_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
             &aux(4*isz),stat=info)
        if (info /= 0) then 
          call psb_errpush(4025,name,i_err=(/4*isz,0,0,0,0/),&
-               & a_err='complex(kind(1.d0))')
+               & a_err='complex(psb_dpk_)')
          goto 9999      
        end if
--- a/mlprec/mld_zbaseprec_aply.f90
+++ b/mlprec/mld_zbaseprec_aply.f90
@ -56,16 +56,16 @@
 !
 !
 ! Arguments:
-!   alpha      -  complex(kind(0.d0)), input.
+!   alpha      -  complex(psb_dpk_), input.
 !                 The scalar alpha.
 !   prec       -  type(mld_zbaseprc_type), input.
 !                 The base preconditioner data structure containing the local part
 !                 of the preconditioner K.
-!   x          -  complex(kind(0.d0)), dimension(:), input.
+!   x          -  complex(psb_dpk_), dimension(:), input.
 !                 The local part of the vector X.
-!   beta       -  complex(kind(0.d0)), input.
+!   beta       -  complex(psb_dpk_), input.
 !                 The scalar beta.
-!   y          -  complex(kind(0.d0)), dimension(:), input/output.
+!   y          -  complex(psb_dpk_), dimension(:), input/output.
 !                 The local part of the vector Y.
 !   desc_data  -  type(psb_desc_type), input.
 !                 The communication descriptor associated to the matrix to be
@ -73,7 +73,7 @@
 !   trans      -  character, optional.
 !                 If trans='N','n' then op(K^(-1)) = K^(-1);
 !                 if trans='T','t' then op(K^(-1)) = K^(-T) (transpose of K^(-1)).
-!   work       -  real(kind(0.d0)), dimension (:), optional, target.
+!   work       -  real(psb_dpk_), dimension (:), optional, target.
 !                 Workspace. Its size must be at least 4*psb_cd_get_local_cols(desc_data).
 !   info       -  integer, output.
 !                 Error code.
@ -88,15 +88,15 @@ subroutine mld_zbaseprec_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
  ! Arguments
  type(psb_desc_type),intent(in)      :: desc_data
  type(mld_zbaseprc_type), intent(in) :: prec
-  complex(kind(0.d0)),intent(in)      :: x(:)
+  complex(psb_dpk_),intent(in)      :: x(:)
-  complex(kind(0.d0)),intent(inout)   :: y(:)
+  complex(psb_dpk_),intent(inout)   :: y(:)
-  complex(kind(0.d0)),intent(in)      :: alpha,beta
+  complex(psb_dpk_),intent(in)      :: alpha,beta
  character(len=1)                    :: trans
-  complex(kind(0.d0)),target          :: work(:)
+  complex(psb_dpk_),target          :: work(:)
  integer, intent(out)                :: info
  ! Local variables
-  complex(kind(1.d0)), pointer :: ww(:)
+  complex(psb_dpk_), pointer :: ww(:)
  integer           :: ictxt, np, me, err_act
  integer           :: n_row, int_err(5)
  character(len=20) :: name, ch_err
@ -140,7 +140,7 @@ subroutine mld_zbaseprec_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
    else
      allocate(ww(size(x)),stat=info)
      if (info /= 0) then 
-        call psb_errpush(4025,name,i_err=(/size(x),0,0,0,0/),a_err='complex(kind(1.d0))')
+        call psb_errpush(4025,name,i_err=(/size(x),0,0,0,0/),a_err='complex(psb_dpk_)')
        goto 9999      
      end if
    end if
--- a/mlprec/mld_zilu0_fact.f90
+++ b/mlprec/mld_zilu0_fact.f90
@ -85,7 +85,7 @@
 !               The U factor (except its diagonal) in the incomplete factorization.
 !               Note: its allocation is managed by the calling routine mld_ilu_bld,
 !               hence it cannot be only intent(out).
-!    d       -  complex(kind(1.d0)), dimension(:), input/output.
+!    d       -  complex(psb_dpk_), dimension(:), input/output.
 !               The inverse of the diagonal entries of the U factor in the incomplete
 !               factorization.
 !               Note: its allocation is managed by the calling routine mld_ilu_bld,
@ -110,7 +110,7 @@ subroutine mld_zilu0_fact(ialg,a,l,u,d,info,blck)
  integer, intent(in)                 :: ialg
  type(psb_zspmat_type),intent(in)    :: a
  type(psb_zspmat_type),intent(inout) :: l,u
-  complex(kind(1.d0)), intent(inout)     :: d(:)
+  complex(psb_dpk_), intent(inout)     :: d(:)
  integer, intent(out)                :: info
  type(psb_zspmat_type),intent(in), optional, target :: blck
@ -249,10 +249,10 @@ contains
  !               to build an Additive Schwarz base preconditioner with overlap
  !               greater than 0. If the overlap is 0 or the matrix has been
  !               reordered (see mld_fact_bld), then b does not contain any row.
-  !    d       -  complex(kind(1.d0)), dimension(:), output.
+  !    d       -  complex(psb_dpk_), dimension(:), output.
  !               The inverse of the diagonal entries of the U factor in the
  !               incomplete factorization.
-  !    laspk   -  complex(kind(1.d0)), dimension(:), input/output.
+  !    laspk   -  complex(psb_dpk_), dimension(:), input/output.
  !               The entries of U are stored according to the CSR format.
  !               The L factor in the incomplete factorization.
  !    lia1    -  integer, dimension(:), input/output.
@ -261,7 +261,7 @@ contains
  !    lia2    -  integer, dimension(:), input/output.
  !               The indices identifying the first nonzero entry of each row
  !               of the L factor in laspk, according to the CSR storage format. 
-  !    uaspk   -  complex(kind(1.d0)), dimension(:), input/output.
+  !    uaspk   -  complex(psb_dpk_), dimension(:), input/output.
  !               The U factor in the incomplete factorization.
  !               The entries of U are stored according to the CSR format.
  !    uia1    -  integer, dimension(:), input/output.
@ -288,11 +288,11 @@ contains
    integer,intent(inout)                            :: m,l1,l2,info
    integer, intent(in)                              :: ma,mb
    integer, dimension(:), intent(inout)             :: lia1,lia2,uia1,uia2
-    complex(kind(1.d0)), dimension(:), intent(inout) :: laspk,uaspk,d
+    complex(psb_dpk_), dimension(:), intent(inout) :: laspk,uaspk,d
    ! Local variables
    integer :: i,j,k,l,low1,low2,kk,jj,ll, ktrw,err_act
-    complex(kind(1.d0)) :: dia,temp
+    complex(psb_dpk_) :: dia,temp
    integer, parameter :: nrb=16
    type(psb_zspmat_type) :: trw
    integer             :: int_err(5) 
@ -519,10 +519,10 @@ contains
  !               The column indices of the nonzero entries of the lower triangle
  !               copied in laspk row by row (see mld_zilu0_factint), according
  !               to the CSR storage format.
-  !    laspk   -  complex(kind(1.d0)), dimension(:), input/output.
+  !    laspk   -  complex(psb_dpk_), dimension(:), input/output.
  !               The array where the entries of the row corresponding to the
  !               lower triangle are copied.
-  !    dia     -  complex(kind(1.d0)), output.
+  !    dia     -  complex(psb_dpk_), output.
  !               The diagonal entry of the copied row.
  !    l2      -  integer, input/output.
  !               Pointer to the last occupied entry of uaspk.
@ -530,7 +530,7 @@ contains
  !               The column indices of the nonzero entries of the upper triangle
  !               copied in uaspk row by row (see mld_zilu0_factint), according
  !               to the CSR storage format.
-  !    uaspk   -  complex(kind(1.d0)), dimension(:), input/output.
+  !    uaspk   -  complex(psb_dpk_), dimension(:), input/output.
  !               The array where the entries of the row corresponding to the
  !               upper triangle are copied. 
  !    ktrw    -  integer, input/output.
@ -557,7 +557,7 @@ contains
    integer, intent(in)                  :: i,m,jd,jmin,jmax
    integer, intent(inout)               :: ktrw,l1,l2
    integer, intent(inout)               :: lia1(:), uia1(:)
-    complex(kind(1.d0)), intent(inout)   :: laspk(:), uaspk(:), dia
+    complex(psb_dpk_), intent(inout)   :: laspk(:), uaspk(:), dia
    ! Local variables
    integer               :: k,j,info,irb
--- a/mlprec/mld_ziluk_fact.f90
+++ b/mlprec/mld_ziluk_fact.f90
@ -82,7 +82,7 @@
 !               The U factor (except its diagonal) in the incomplete factorization.
 !               Note: its allocation is managed by the calling routine mld_ilu_bld,
 !               hence it cannot be only intent(out).
-!    d       -  complex(kind(1.d0)), dimension(:), input/output.
+!    d       -  complex(psb_dpk_), dimension(:), input/output.
 !               The inverse of the diagonal entries of the U factor in the incomplete
 !               factorization.
 !               Note: its allocation is managed by the calling routine mld_ilu_bld,
@ -109,7 +109,7 @@ subroutine mld_ziluk_fact(fill_in,ialg,a,l,u,d,info,blck)
  type(psb_zspmat_type),intent(in)    :: a
  type(psb_zspmat_type),intent(inout) :: l,u
  type(psb_zspmat_type),intent(in), optional, target :: blck
-  complex(kind(1.d0)), intent(inout)     ::  d(:)
+  complex(psb_dpk_), intent(inout)     ::  d(:)
  !     Local Variables
  integer   :: l1, l2, m, err_act
@ -237,10 +237,10 @@ contains
  !               to build an Additive Schwarz base preconditioner with overlap
  !               greater than 0. If the overlap is 0 or the matrix has been reordered
  !               (see mld_fact_bld), then b does not contain   any row.
-  !    d       -  complex(kind(1.d0)), dimension(:), output.
+  !    d       -  complex(psb_dpk_), dimension(:), output.
  !               The inverse of the diagonal entries of the U factor in the incomplete
  !               factorization.
-  !    laspk   -  complex(kind(1.d0)), dimension(:), input/output.
+  !    laspk   -  complex(psb_dpk_), dimension(:), input/output.
  !               The L factor in the incomplete factorization.
  !    lia1    -  integer, dimension(:), input/output.
  !               The column indices of the nonzero entries of the L factor,
@ -248,7 +248,7 @@ contains
  !    lia2    -  integer, dimension(:), input/output.
  !               The indices identifying the first nonzero entry of each row
  !               of the L factor in laspk, according to the CSR storage format. 
-  !    uaspk   -  complex(kind(1.d0)), dimension(:), input/output.
+  !    uaspk   -  complex(psb_dpk_), dimension(:), input/output.
  !               The U factor in the incomplete factorization.
  !               The entries of U are stored according to the CSR format.
  !    uia1    -  integer, dimension(:), input/output.
@ -276,13 +276,13 @@ contains
    type(psb_zspmat_type), intent(in)               :: a,b
    integer, intent(inout)                          :: m,l1,l2,info
    integer, allocatable, intent(inout)             :: lia1(:),lia2(:),uia1(:),uia2(:)
-    complex(kind(1.d0)), allocatable, intent(inout) :: laspk(:),uaspk(:)
+    complex(psb_dpk_), allocatable, intent(inout) :: laspk(:),uaspk(:)
-    complex(kind(1.d0)), intent(inout)              :: d(:)
+    complex(psb_dpk_), intent(inout)              :: d(:)
  ! Local variables
    integer :: ma,mb,i, ktrw,err_act,nidx
    integer, allocatable            :: uplevs(:), rowlevs(:),idxs(:)
-    complex(kind(1.d0)), allocatable :: row(:)
+    complex(psb_dpk_), allocatable :: row(:)
    type(psb_int_heap) :: heap
    logical,parameter  :: debug=.false.
    type(psb_zspmat_type) :: trw
@ -456,7 +456,7 @@ contains
  !               The maximum valid column index.
  !               The output matrix will contain a clipped copy taken from
  !               a(1:m,jmin:jmax).
-  !    row     -  complex(kind(1.d0)), dimension(:), input/output.
+  !    row     -  complex(psb_dpk_), dimension(:), input/output.
  !               In input it is the null vector (see mld_iluk_factint and
  !               iluk_copyout). In output it contains the row extracted
  !               from the matrix A. It actually contains a full row, i.e.
@ -494,7 +494,7 @@ contains
    integer, intent(in)                  :: i,m,jmin,jmax
    integer, intent(inout)               :: ktrw,info
    integer, intent(inout)               :: rowlevs(:)
-    complex(kind(1.d0)), intent(inout)   :: row(:)
+    complex(psb_dpk_), intent(inout)   :: row(:)
    type(psb_int_heap), intent(inout)    :: heap
  ! Local variables
@ -592,7 +592,7 @@ contains
  !    i       -  integer, input.
  !               The local index of the row to which the factorization is
  !               applied.
-  !    row     -  complex(kind(1.d0)), dimension(:), input/output.
+  !    row     -  complex(psb_dpk_), dimension(:), input/output.
  !               In input it contains the row to which the elimination step
  !               has to be applied. In output it contains the row after the
  !               elimination step. It actually contains a full row, i.e.
@ -609,7 +609,7 @@ contains
  !               in the processed row. In input it contains the indices concerning
  !               the row before the elimination step, while in output it contains
  !               the indices concerning the transformed row.
-  !    d       -  complex(kind(1.d0)), input.
+  !    d       -  complex(psb_dpk_), input.
  !               The inverse of the diagonal entries of the part of the U factor
  !               above the current row (see iluk_copyout).
  !    uia1    -  integer, dimension(:), input.
@ -622,7 +622,7 @@ contains
  !               the U factor above the current row, stored in uaspk row by row
  !               (see iluk_copyout, called by mld_ziluk_factint), according to
  !               the CSR storage format.
-  !    uaspk   -  complex(kind(1.d0)), dimension(:), input.
+  !    uaspk   -  complex(psb_dpk_), dimension(:), input.
  !               The entries of the U factor above the current row (except the
  !               diagonal ones), stored according to the CSR format.
  !    uplevs  -  integer, dimension(:), input.
@ -652,11 +652,11 @@ contains
    integer, intent(inout)               :: rowlevs(:)
    integer, allocatable, intent(inout)  :: idxs(:)
    integer, intent(inout)               :: uia1(:),uia2(:),uplevs(:)
-    complex(kind(1.d0)), intent(inout)   :: row(:), uaspk(:),d(:)
+    complex(psb_dpk_), intent(inout)   :: row(:), uaspk(:),d(:)
    ! Local variables
    integer               :: k,j,lrwk,jj,lastk, iret
-    complex(kind(1.d0))   :: rwk
+    complex(psb_dpk_)   :: rwk
    info = 0
    if (.not.allocated(idxs)) then 
@ -762,7 +762,7 @@ contains
  !               The local index of the row to be copied.
  !    m       -  integer, input.
  !               The number of rows of the local matrix under factorization.
-  !    row     -  complex(kind(1.d0)), dimension(:), input/output.
+  !    row     -  complex(psb_dpk_), dimension(:), input/output.
  !               It contains, input, the row to be copied, and, in output,
  !               the null vector (the latter is used in the next call to
  !               iluk_copyin in mld_iluk_fact).
@ -792,10 +792,10 @@ contains
  !               The indices identifying the first nonzero entry of each row
  !               of the L factor, copied in laspk row by row (see 
  !               mld_ziluk_factint), according to the CSR storage format.
-  !    laspk   -  complex(kind(1.d0)), dimension(:), input/output.
+  !    laspk   -  complex(psb_dpk_), dimension(:), input/output.
  !               The array where the entries of the row corresponding to the
  !               L factor are copied.
-  !    d       -  complex(kind(1.d0)), dimension(:), input/output.
+  !    d       -  complex(psb_dpk_), dimension(:), input/output.
  !               The array where the inverse of the diagonal entry of the
  !               row is copied (only d(i) is used by the routine). 
  !    uia1    -  integer, dimension(:), input/output.
@ -806,7 +806,7 @@ contains
  !               The indices identifying the first nonzero entry of each row
  !               of the U factor copied in uaspk row by row (see
  !               mld_zilu_fctint), according to the CSR storage format.
-  !    uaspk   -  complex(kind(1.d0)), dimension(:), input/output.
+  !    uaspk   -  complex(psb_dpk_), dimension(:), input/output.
  !               The array where the entries of the row corresponding to the
  !               U factor are copied.
  !    uplevs  -  integer, dimension(:), input.
@ -825,8 +825,8 @@ contains
    integer, intent(inout)                          :: l1, l2, info
    integer, intent(inout)                          :: rowlevs(:), idxs(:)
    integer, allocatable, intent(inout)             :: uia1(:), uia2(:), lia1(:), lia2(:),uplevs(:)
-    complex(kind(1.d0)), allocatable, intent(inout) :: uaspk(:), laspk(:)
+    complex(psb_dpk_), allocatable, intent(inout) :: uaspk(:), laspk(:)
-    complex(kind(1.d0)), intent(inout)              :: row(:), d(:)
+    complex(psb_dpk_), intent(inout)              :: row(:), d(:)
    ! Local variables
    integer               :: j,isz,err_act,int_err(5),idxp
--- a/mlprec/mld_zilut_fact.f90
+++ b/mlprec/mld_zilut_fact.f90
@ -78,7 +78,7 @@
 !               The U factor (except its diagonal) in the incomplete factorization.
 !               Note: its allocation is managed by the calling routine mld_ilu_bld,
 !               hence it cannot be only intent(out).
-!    d       -  complex(kind(1.d0)), dimension(:), input/output.
+!    d       -  complex(psb_dpk_), dimension(:), input/output.
 !               The inverse of the diagonal entries of the U factor in the incomplete
 !               factorization.
 !               Note: its allocation is managed by the calling routine mld_ilu_bld,
@ -101,11 +101,11 @@ subroutine mld_zilut_fact(fill_in,thres,a,l,u,d,info,blck)
  ! Arguments
  integer, intent(in)                 :: fill_in
-  real(kind(1.d0)), intent(in)        :: thres
+  real(psb_dpk_), intent(in)        :: thres
  integer, intent(out)                :: info
  type(psb_zspmat_type),intent(in)    :: a
  type(psb_zspmat_type),intent(inout) :: l,u
-  complex(kind(1.d0)), intent(inout)     :: d(:)
+  complex(psb_dpk_), intent(inout)     :: d(:)
  type(psb_zspmat_type),intent(in), optional, target :: blck
  !     Local Variables
@ -238,10 +238,10 @@ contains
  !               to build an Additive Schwarz base preconditioner with overlap
  !               greater than 0. If the overlap is 0 or the matrix has been reordered
  !               (see mld_fact_bld), then b does not contain   any row.
-  !    d       -  complex(kind(1.d0)), dimension(:), output.
+  !    d       -  complex(psb_dpk_), dimension(:), output.
  !               The inverse of the diagonal entries of the U factor in the incomplete
  !               factorization.
-  !    laspk   -  complex(kind(1.d0)), dimension(:), input/output.
+  !    laspk   -  complex(psb_dpk_), dimension(:), input/output.
  !               The L factor in the incomplete factorization.
  !    lia1    -  integer, dimension(:), input/output.
  !               The column indices of the nonzero entries of the L factor,
@ -249,7 +249,7 @@ contains
  !    lia2    -  integer, dimension(:), input/output.
  !               The indices identifying the first nonzero entry of each row
  !               of the L factor in laspk, according to the CSR storage format. 
-  !    uaspk   -  complex(kind(1.d0)), dimension(:), input/output.
+  !    uaspk   -  complex(psb_dpk_), dimension(:), input/output.
  !               The U factor in the incomplete factorization.
  !               The entries of U are stored according to the CSR format.
  !    uia1    -  integer, dimension(:), input/output.
@ -274,18 +274,18 @@ contains
  ! Arguments
    integer, intent(in)                             :: fill_in
-    real(kind(1.d0)), intent(in)                    :: thres
+    real(psb_dpk_), intent(in)                    :: thres
    type(psb_zspmat_type), intent(in)               :: a,b
    integer, intent(inout)                          :: m,l1,l2,info
    integer, allocatable, intent(inout)             :: lia1(:),lia2(:),uia1(:),uia2(:)
-    complex(kind(1.d0)), allocatable, intent(inout) :: laspk(:),uaspk(:)
+    complex(psb_dpk_), allocatable, intent(inout) :: laspk(:),uaspk(:)
-    complex(kind(1.d0)), intent(inout)              :: d(:)
+    complex(psb_dpk_), intent(inout)              :: d(:)
    ! Local Variables
    integer :: i, ktrw,err_act,nidx,nlw,nup,jmaxup, ma, mb
-    real(kind(1.d0)) ::  nrmi
+    real(psb_dpk_) ::  nrmi
    integer, allocatable          :: idxs(:)
-    complex(kind(1.d0)), allocatable :: row(:)
+    complex(psb_dpk_), allocatable :: row(:)
    type(psb_int_heap) :: heap
    type(psb_zspmat_type) :: trw
    character(len=20), parameter  :: name='mld_zilut_factint'
@ -451,9 +451,9 @@ contains
  !    jmaxup  -  integer, output.
  !               The column index of the first entry with maximum absolute
  !               value in the part of the row belonging to the upper triangle
-  !    nrmi    -  real(kind(1.d0)), output.
+  !    nrmi    -  real(psb_dpk_), output.
  !               The 2-norm of the current row.
-  !    row     -  complex(kind(1.d0)), dimension(:), input/output.
+  !    row     -  complex(psb_dpk_), dimension(:), input/output.
  !               In input it is the null vector (see mld_ilut_factint and
  !               ilut_copyout). In output it contains the row extracted
  !               from the matrix A. It actually contains a full row, i.e.
@ -486,14 +486,14 @@ contains
    type(psb_zspmat_type), intent(inout) :: trw
    integer, intent(in)                  :: i, m,jmin,jmax,jd
    integer, intent(inout)               :: ktrw,nlw,nup,jmaxup,info
-    real(kind(1.d0)), intent(inout)      :: nrmi
+    real(psb_dpk_), intent(inout)      :: nrmi
-    complex(kind(1.d0)), intent(inout)      :: row(:)
+    complex(psb_dpk_), intent(inout)      :: row(:)
    type(psb_int_heap), intent(inout)    :: heap
    integer               :: k,j,irb,kin,nz
    integer, parameter    :: nrb=16
-    real(kind(1.d0))      :: dmaxup
+    real(psb_dpk_)      :: dmaxup
-    real(kind(1.d0)), external    :: dznrm2
+    real(psb_dpk_), external    :: dznrm2
    character(len=20), parameter  :: name='mld_zilut_factint'
    if (psb_get_errstatus() /= 0) return 
@ -628,10 +628,10 @@ contains
  !               The threshold t, i.e. the drop tolerance, in ILU(k,t).
  !    i       -  integer, input.
  !               The local index of the row to which the factorization is applied.
-  !    nrmi    -  real(kind(1.d0)), input.
+  !    nrmi    -  real(psb_dpk_), input.
  !               The 2-norm of the row to which the elimination step has to be
  !               applied.
-  !    row     -  complex(kind(1.d0)), dimension(:), input/output.
+  !    row     -  complex(psb_dpk_), dimension(:), input/output.
  !               In input it contains the row to which the elimination step
  !               has to be applied. In output it contains the row after the
  !               elimination step. It actually contains a full row, i.e.
@ -642,7 +642,7 @@ contains
  !               the row before the elimination step, while in output it contains
  !               the previous indices plus the ones corresponding to transformed
  !               entries in the 'upper part' that have not been dropped.
-  !    d       -  complex(kind(1.d0)), input.
+  !    d       -  complex(psb_dpk_), input.
  !               The inverse of the diagonal entries of the part of the U factor
  !               above the current row (see ilut_copyout).
  !    uia1    -  integer, dimension(:), input.
@ -655,7 +655,7 @@ contains
  !               the U factor above the current row, stored in uaspk row by row
  !               (see ilut_copyout, called by mld_zilut_factint), according to
  !               the CSR storage format.
-  !    uaspk   -  complex(kind(1.d0)), dimension(:), input.
+  !    uaspk   -  complex(psb_dpk_), dimension(:), input.
  !               The entries of the U factor above the current row (except the
  !               diagonal ones), stored according to the CSR format.
  !    nidx    -  integer, output.
@ -679,14 +679,14 @@ contains
    type(psb_int_heap), intent(inout)   :: heap 
    integer, intent(in)                 :: i
    integer, intent(inout)              :: nidx,info
-    real(kind(1.d0)), intent(in)        :: thres,nrmi
+    real(psb_dpk_), intent(in)        :: thres,nrmi
    integer, allocatable, intent(inout) :: idxs(:)
    integer, intent(inout)              :: uia1(:),uia2(:)
-    complex(kind(1.d0)), intent(inout)     :: row(:), uaspk(:),d(:)
+    complex(psb_dpk_), intent(inout)     :: row(:), uaspk(:),d(:)
    ! Local Variables
    integer               :: k,j,jj,lastk, iret
-    complex(kind(1.d0))   :: rwk
+    complex(psb_dpk_)   :: rwk
    info  = 0
    call psb_ensure_size(200,idxs,info)
@ -811,9 +811,9 @@ contains
  !    jmaxup  -  integer, input.
  !               The column index of the first entry with maximum absolute
  !               value in the 'upper part' of the row in the initial matrix.
-  !    nrmi    -  real(kind(1.d0)), input.
+  !    nrmi    -  real(psb_dpk_), input.
  !               The 2-norm of the current row in the initial matrix.
-  !    row     -  complex(kind(1.d0)), dimension(:), input/output.
+  !    row     -  complex(psb_dpk_), dimension(:), input/output.
  !               It contains, input, the row to be copied, and, in output,
  !               the null vector (the latter is used in the next call to
  !               ilut_copyin in mld_ilut_fact).
@ -836,10 +836,10 @@ contains
  !               The indices identifying the first nonzero entry of each row
  !               of the L factor, copied in laspk row by row (see 
  !               mld_zilut_factint), according to the CSR storage format.
-  !    laspk   -  complex(kind(1.d0)), dimension(:), input/output.
+  !    laspk   -  complex(psb_dpk_), dimension(:), input/output.
  !               The array where the entries of the row corresponding to the
  !               L factor are copied.
-  !    d       -  complex(kind(1.d0)), dimension(:), input/output.
+  !    d       -  complex(psb_dpk_), dimension(:), input/output.
  !               The array where the inverse of the diagonal entry of the
  !               row is copied (only d(i) is used by the routine). 
  !    uia1    -  integer, dimension(:), input/output.
@ -850,7 +850,7 @@ contains
  !               The indices identifying the first nonzero entry of each row
  !               of the U factor copied in uaspk row by row (see
  !               mld_zilu_fctint), according to the CSR storage format.
-  !    uaspk   -  complex(kind(1.d0)), dimension(:), input/output.
+  !    uaspk   -  complex(psb_dpk_), dimension(:), input/output.
  !               The array where the entries of the row corresponding to the
  !               U factor are copied.
  !
@ -866,14 +866,14 @@ contains
    integer, intent(in)                            :: idxs(:)
    integer, intent(inout)                         :: l1,l2, info
    integer, allocatable, intent(inout)            :: uia1(:),uia2(:), lia1(:),lia2(:)
-    real(kind(1.d0)), intent(in)                   :: thres,nrmi
+    real(psb_dpk_), intent(in)                   :: thres,nrmi
-    complex(kind(1.d0)),allocatable, intent(inout) :: uaspk(:), laspk(:)
+    complex(psb_dpk_),allocatable, intent(inout) :: uaspk(:), laspk(:)
-    complex(kind(1.d0)), intent(inout)             :: row(:), d(:)
+    complex(psb_dpk_), intent(inout)             :: row(:), d(:)
    ! Local variables
-    complex(kind(1.d0)),allocatable :: xw(:)
+    complex(psb_dpk_),allocatable :: xw(:)
    integer, allocatable         :: xwid(:), indx(:)
-    complex(kind(1.d0))             :: witem
+    complex(psb_dpk_)             :: witem
    integer                      :: widx
    integer                      :: k,isz,err_act,int_err(5),idxp, nz
    type(psb_dcomplex_idx_heap)  :: heap
@ -900,7 +900,7 @@ contains
    if (info /= 0) then 
      info=4025
      call psb_errpush(info,name,i_err=(/3*nidx,0,0,0,0/),&
-           & a_err='complex(kind(1.d0))')
+           & a_err='complex(psb_dpk_)')
      goto 9999      
    end if
--- a/mlprec/mld_zmlprec_aply.f90
+++ b/mlprec/mld_zmlprec_aply.f90
@ -76,7 +76,7 @@
 !
 !
 ! Arguments:
-!  alpha       -   complex(kind(0.d0)), input.
+!  alpha       -   complex(psb_dpk_), input.
 !                  The scalar alpha.
 !   baseprecv  -   type(mld_zbaseprc_type), dimension(:), input.
 !                  The array of base preconditioner data structures containing the
@ -98,7 +98,7 @@
 !                                             It maps vectors (ilev) ---> (ilev-1).
 !         baseprecv(ilev)%av(mld_sm_pr_t_) -  The smoothed prolongator transpose.   
 !                                             It maps vectors (ilev-1) ---> (ilev).
-!      baseprecv(ilev)%d         -  complex(kind(1.d0)), dimension(:), allocatable.
+!      baseprecv(ilev)%d         -  complex(psb_dpk_), dimension(:), allocatable.
 !                                   The diagonal entries of the U factor in the ILU
 !                                   factorization of A(ilev).
 !      baseprecv(ilev)%desc_data -  type(psb_desc_type).
@ -111,7 +111,7 @@
 !      baseprecv(ilev)%iprcparm  -  integer, dimension(:), allocatable.
 !                                   The integer parameters defining the base
 !                                   preconditioner K(ilev).
-!      baseprecv(ilev)%dprcparm  -  complex(kind(1.d0)), dimension(:), allocatable.
+!      baseprecv(ilev)%dprcparm  -  complex(psb_dpk_), dimension(:), allocatable.
 !                                   The real parameters defining the base preconditioner
 !                                   K(ilev).
 !      baseprecv(ilev)%perm      -  integer, dimension(:), allocatable.
@ -138,14 +138,14 @@
 !      baseprecv(ilev)%base_desc -  type(psb_desc_type), pointer.
 !                                   Pointer to the communication descriptor associated
 !                                   to the sparse matrix pointed by base_a.  
-!      baseprecv(ilev)%dorig     -  complex(kind(1.d0)), dimension(:), allocatable.
+!      baseprecv(ilev)%dorig     -  complex(psb_dpk_), dimension(:), allocatable.
 !                                   Diagonal entries of the matrix pointed by base_a.
 !                  
-!   x          -  complex(kind(0.d0)), dimension(:), input.
+!   x          -  complex(psb_dpk_), dimension(:), input.
 !                 The local part of the vector X.
-!   beta       -  complex(kind(0.d0)), input.
+!   beta       -  complex(psb_dpk_), input.
 !                 The scalar beta.
-!   y          -  complex(kind(0.d0)), dimension(:), input/output.
+!   y          -  complex(psb_dpk_), dimension(:), input/output.
 !                 The local part of the vector Y.
 !   desc_data  -  type(psb_desc_type), input.
 !                 The communication descriptor associated to the matrix to be
@ -153,7 +153,7 @@
 !   trans      -  character, optional.
 !                 If trans='N','n' then op(M^(-1)) = M^(-1);
 !                 if trans='T','t' then op(M^(-1)) = M^(-T) (transpose of M^(-1)).
-!   work       -  complex(kind(0.d0)), dimension (:), optional, target.
+!   work       -  complex(psb_dpk_), dimension (:), optional, target.
 !                 Workspace. Its size must be at least 4*psb_cd_get_local_cols(desc_data).
 !   info       -  integer, output.
 !                 Error code.
@ -174,11 +174,11 @@ subroutine mld_zmlprec_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
  ! Arguments
  type(psb_desc_type),intent(in)      :: desc_data
  type(mld_zbaseprc_type), intent(in) :: baseprecv(:)
-  complex(kind(1.d0)),intent(in)      :: alpha,beta
+  complex(psb_dpk_),intent(in)      :: alpha,beta
-  complex(kind(1.d0)),intent(in)      :: x(:)
+  complex(psb_dpk_),intent(in)      :: x(:)
-  complex(kind(1.d0)),intent(inout)   :: y(:)
+  complex(psb_dpk_),intent(inout)   :: y(:)
  character, intent(in)               :: trans
-  complex(kind(1.d0)),target          :: work(:)
+  complex(psb_dpk_),target          :: work(:)
  integer, intent(out)                :: info
  ! Local variables
@ -367,11 +367,11 @@ contains
    ! Arguments
    type(psb_desc_type),intent(in)      :: desc_data
    type(mld_zbaseprc_type), intent(in) :: baseprecv(:)
-    complex(kind(1.d0)),intent(in)      :: alpha,beta
+    complex(psb_dpk_),intent(in)      :: alpha,beta
-    complex(kind(1.d0)),intent(in)      :: x(:)
+    complex(psb_dpk_),intent(in)      :: x(:)
-    complex(kind(1.d0)),intent(inout)   :: y(:)
+    complex(psb_dpk_),intent(inout)   :: y(:)
    character, intent(in)               :: trans
-    complex(kind(1.d0)),target          :: work(:)
+    complex(psb_dpk_),target          :: work(:)
    integer, intent(out)                :: info
    ! Local variables
@ -382,7 +382,7 @@ contains
    character(len=20)  :: name
    type psb_mlprec_wrk_type
-      complex(kind(1.d0)), allocatable  :: tx(:),ty(:),x2l(:),y2l(:)
+      complex(psb_dpk_), allocatable  :: tx(:),ty(:),x2l(:),y2l(:)
    end type psb_mlprec_wrk_type
    type(psb_mlprec_wrk_type), allocatable :: mlprec_wrk(:)
@ -415,7 +415,7 @@ contains
    if (info /= 0) then 
      info=4025
      call psb_errpush(info,name,i_err=(/size(x)+size(y),0,0,0,0/),&
-           & a_err='real(kind(1.d0))')
+           & a_err='real(psb_dpk_)')
      goto 9999      
    end if
@ -444,7 +444,7 @@ contains
      if (info /= 0) then 
        info=4025
        call psb_errpush(info,name,i_err=(/2*(nc2l+max(n_row,n_col)),0,0,0,0/),&
-             & a_err='real(kind(1.d0))')
+             & a_err='real(psb_dpk_)')
        goto 9999      
      end if
@ -625,11 +625,11 @@ contains
    ! Arguments
    type(psb_desc_type),intent(in)      :: desc_data
    type(mld_zbaseprc_type), intent(in) :: baseprecv(:)
-    complex(kind(1.d0)),intent(in)      :: alpha,beta
+    complex(psb_dpk_),intent(in)      :: alpha,beta
-    complex(kind(1.d0)),intent(in)      :: x(:)
+    complex(psb_dpk_),intent(in)      :: x(:)
-    complex(kind(1.d0)),intent(inout)   :: y(:)
+    complex(psb_dpk_),intent(inout)   :: y(:)
    character, intent(in)               :: trans
-    complex(kind(1.d0)),target          :: work(:)
+    complex(psb_dpk_),target          :: work(:)
    integer, intent(out)                :: info
    ! Local variables
@ -640,7 +640,7 @@ contains
    character(len=20)  :: name
    type psb_mlprec_wrk_type
-      complex(kind(1.d0)), allocatable  :: tx(:),ty(:),x2l(:),y2l(:)
+      complex(psb_dpk_), allocatable  :: tx(:),ty(:),x2l(:),y2l(:)
    end type psb_mlprec_wrk_type
    type(psb_mlprec_wrk_type), allocatable :: mlprec_wrk(:)
@ -677,7 +677,7 @@ contains
    if (info /= 0) then 
      info=4025
      call psb_errpush(info,name,i_err=(/4*nc2l,0,0,0,0/),&
-           & a_err='real(kind(1.d0))')
+           & a_err='real(psb_dpk_)')
      goto 9999
    end if
@ -729,7 +729,7 @@ contains
      if (info /= 0) then 
        info=4025
        call psb_errpush(info,name,i_err=(/4*nc2l,0,0,0,0/),&
-             & a_err='real(kind(1.d0))')
+             & a_err='real(psb_dpk_)')
        goto 9999      
      end if
@ -904,11 +904,11 @@ contains
    ! Arguments
    type(psb_desc_type),intent(in)      :: desc_data
    type(mld_zbaseprc_type), intent(in) :: baseprecv(:)
-    complex(kind(1.d0)),intent(in)      :: alpha,beta
+    complex(psb_dpk_),intent(in)      :: alpha,beta
-    complex(kind(1.d0)),intent(in)      :: x(:)
+    complex(psb_dpk_),intent(in)      :: x(:)
-    complex(kind(1.d0)),intent(inout)   :: y(:)
+    complex(psb_dpk_),intent(inout)   :: y(:)
    character, intent(in)               :: trans
-    complex(kind(1.d0)),target          :: work(:)
+    complex(psb_dpk_),target          :: work(:)
    integer, intent(out)                :: info
    ! Local variables
@ -919,7 +919,7 @@ contains
    character(len=20)  :: name
    type psb_mlprec_wrk_type
-      complex(kind(1.d0)), allocatable  :: tx(:),ty(:),x2l(:),y2l(:)
+      complex(psb_dpk_), allocatable  :: tx(:),ty(:),x2l(:),y2l(:)
    end type psb_mlprec_wrk_type
    type(psb_mlprec_wrk_type), allocatable :: mlprec_wrk(:)
@ -989,7 +989,7 @@ contains
      if (info /= 0) then 
        info=4025
        call psb_errpush(info,name,i_err=(/4*nc2l,0,0,0,0/),&
-             & a_err='real(kind(1.d0))')
+             & a_err='real(psb_dpk_)')
        goto 9999      
      end if
@ -1214,11 +1214,11 @@ contains
    ! Arguments
    type(psb_desc_type),intent(in)      :: desc_data
    type(mld_zbaseprc_type), intent(in) :: baseprecv(:)
-    complex(kind(1.d0)),intent(in)      :: alpha,beta
+    complex(psb_dpk_),intent(in)      :: alpha,beta
-    complex(kind(1.d0)),intent(in)      :: x(:)
+    complex(psb_dpk_),intent(in)      :: x(:)
-    complex(kind(1.d0)),intent(inout)   :: y(:)
+    complex(psb_dpk_),intent(inout)   :: y(:)
    character, intent(in)               :: trans
-    complex(kind(1.d0)),target          :: work(:)
+    complex(psb_dpk_),target          :: work(:)
    integer, intent(out)                :: info
    ! Local variables
@ -1229,7 +1229,7 @@ contains
    character(len=20)  :: name
    type psb_mlprec_wrk_type
-      complex(kind(1.d0)), allocatable  :: tx(:),ty(:),x2l(:),y2l(:)
+      complex(psb_dpk_), allocatable  :: tx(:),ty(:),x2l(:),y2l(:)
    end type psb_mlprec_wrk_type
    type(psb_mlprec_wrk_type), allocatable :: mlprec_wrk(:)
@ -1266,7 +1266,7 @@ contains
    if (info /= 0) then 
      info=4025
      call psb_errpush(info,name,i_err=(/4*nc2l,0,0,0,0/),&
-           & a_err='real(kind(1.d0))')
+           & a_err='real(psb_dpk_)')
      goto 9999
    end if
@ -1316,7 +1316,7 @@ contains
      if (info /= 0) then 
        info=4025
        call psb_errpush(info,name,i_err=(/4*nc2l,0,0,0,0/),&
-             & a_err='real(kind(1.d0))')
+             & a_err='real(psb_dpk_)')
        goto 9999      
      end if
--- a/mlprec/mld_zprec_aply.f90
+++ b/mlprec/mld_zprec_aply.f90
@ -55,9 +55,9 @@
 !    prec       -  type(mld_zprec_type), input.
 !                  The preconditioner data structure containing the local part
 !                  of the preconditioner to be applied.
-!    x          -  complex(kind(0.d0)), dimension(:), input.
+!    x          -  complex(psb_dpk_), dimension(:), input.
 !                  The local part of the vector X in Y := op(M^(-1)) * X.
-!    y          -  complex(kind(0.d0)), dimension(:), output.
+!    y          -  complex(psb_dpk_), dimension(:), output.
 !                  The local part of the vector Y in Y := op(M^(-1)) * X.
 !    desc_data  -  type(psb_desc_type), input.
 !                  The communication descriptor associated to the matrix to be
@ -67,7 +67,7 @@
 !    trans      -  character(len=1), optional.
 !                  If trans='N','n' then op(M^(-1)) = M^(-1);
 !                  if trans='T','t' then op(M^(-1)) = M^(-T) (transpose of M^(-1)).
-!    work       -  complex(kind(0.d0)), dimension (:), optional, target.
+!    work       -  complex(psb_dpk_), dimension (:), optional, target.
 !                  Workspace. Its size must be at
 !                  least 4*psb_cd_get_local_cols(desc_data).
 !    
@ -82,15 +82,15 @@ subroutine mld_zprec_aply(prec,x,y,desc_data,info,trans,work)
  ! Arguments
  type(psb_desc_type),intent(in)      :: desc_data
  type(mld_zprec_type), intent(in)    :: prec
-  complex(kind(0.d0)),intent(in)      :: x(:)
+  complex(psb_dpk_),intent(in)      :: x(:)
-  complex(kind(0.d0)),intent(inout)   :: y(:)
+  complex(psb_dpk_),intent(inout)   :: y(:)
  integer, intent(out)                :: info
  character(len=1), optional          :: trans
-  complex(kind(0.d0)), optional, target  :: work(:)
+  complex(psb_dpk_), optional, target  :: work(:)
  ! Local variables
  character     :: trans_ 
-  complex(kind(1.d0)), pointer :: work_(:)
+  complex(psb_dpk_), pointer :: work_(:)
  integer :: ictxt,np,me,err_act,iwsz
  character(len=20)   :: name
@ -114,7 +114,7 @@ subroutine mld_zprec_aply(prec,x,y,desc_data,info,trans,work)
    allocate(work_(iwsz),stat=info)
    if (info /= 0) then 
      call psb_errpush(4025,name,i_err=(/iwsz,0,0,0,0/),&
-           & a_err='complex(kind(1.d0))')
+           & a_err='complex(psb_dpk_)')
      goto 9999      
    end if
@ -187,7 +187,7 @@ end subroutine mld_zprec_aply
 !    prec       -  type(mld_zprec_type), input.
 !                  The preconditioner data structure containing the local part
 !                  of the preconditioner to be applied.
-!    x          -  complex(kind(0.d0)), dimension(:), input/output.
+!    x          -  complex(psb_dpk_), dimension(:), input/output.
 !                  The local part of vector X in X := op(M^(-1)) * X.
 !    desc_data  -  type(psb_desc_type), input.
 !                  The communication descriptor associated to the matrix to be
@ -209,13 +209,13 @@ subroutine mld_zprec_aply1(prec,x,desc_data,info,trans)
  ! Arguments
  type(psb_desc_type),intent(in)    :: desc_data
  type(mld_zprec_type), intent(in)  :: prec
-  complex(kind(0.d0)),intent(inout) :: x(:)
+  complex(psb_dpk_),intent(inout) :: x(:)
  integer, intent(out)              :: info
  character(len=1), optional        :: trans
  ! Local variables
  integer :: ictxt,np,me, err_act
-  complex(kind(1.d0)), pointer :: WW(:), w1(:)
+  complex(psb_dpk_), pointer :: WW(:), w1(:)
  character(len=20)   :: name
  name='mld_zprec_aply1'
@ -230,7 +230,7 @@ subroutine mld_zprec_aply1(prec,x,desc_data,info,trans)
  if (info /= 0) then 
    info=4025
    call psb_errpush(info,name,i_err=(/2*size(x),0,0,0,0/),&
-         & a_err='complex(kind(1.d0))')
+         & a_err='complex(psb_dpk_)')
    goto 9999      
  end if
--- a/mlprec/mld_zprecset.f90
+++ b/mlprec/mld_zprecset.f90
@ -519,7 +519,7 @@ end subroutine mld_zprecsetc
 !               The number identifying the parameter to be set.
 !               A mnemonic constant has been associated to each of these
 !               numbers, as reported in MLD2P4 user's guide.
-!    val     -  real(kind(1.d0)), input.
+!    val     -  real(psb_dpk_), input.
 !               The value of the parameter to be set. The list of allowed
 !               values is reported in MLD2P4 user's guide.
 !    info    -  integer, output.
@ -540,7 +540,7 @@ subroutine mld_zprecsetd(p,what,val,info,ilev)
  ! Arguments
  type(mld_zprec_type), intent(inout)    :: p
  integer, intent(in)                    :: what 
-  real(kind(1.d0)), intent(in)           :: val
+  real(psb_dpk_), intent(in)           :: val
  integer, intent(out)                   :: info
  integer, optional, intent(in)          :: ilev
--- a/mlprec/mld_zsp_renum.f90
+++ b/mlprec/mld_zsp_renum.f90
@ -99,7 +99,7 @@ subroutine mld_zsp_renum(a,blck,p,atmp,info)
  integer              :: nztota, nztotb, nztmp, nnr, i,k
  integer, allocatable :: itmp(:), itmp2(:)
  integer              :: ictxt,np,me, err_act
-  real(kind(1.d0))     :: t3,t4
+  real(psb_dpk_)     :: t3,t4
  if (psb_get_errstatus().ne.0) return 
  info=0
--- a/mlprec/mld_zsub_aply.f90
+++ b/mlprec/mld_zsub_aply.f90
@ -108,16 +108,16 @@
 !
 ! Arguments:
 !
-!   alpha      -  complex(kind(0.d0)), input.
+!   alpha      -  complex(psb_dpk_), input.
 !                 The scalar alpha.
 !   prec       -  type(mld_zbaseprec_type), input.
 !                 The 'base preconditioner' data structure containing the local 
 !                 part of the preconditioner or solver.
-!   x          -  complex(kind(0.d0)), dimension(:), input.
+!   x          -  complex(psb_dpk_), dimension(:), input.
 !                 The local part of the vector X.
-!   beta       -  complex(kind(0.d0)), input.
+!   beta       -  complex(psb_dpk_), input.
 !                 The scalar beta.
-!   y          -  complex(kind(0.d0)), dimension(:), input/output.
+!   y          -  complex(psb_dpk_), dimension(:), input/output.
 !                 The local part of the vector Y.
 !   desc_data  -  type(psb_desc_type), input.
 !                 The communication descriptor associated to the matrix to be
@ -128,7 +128,7 @@
 !                 if trans='C','c' then op(K^(-1)) = K^(-C) (transpose conjugate of K^(-1)).
 !                 If prec%iprcparm(smooth_sweeps_) > 1, the value of trans provided
 !                 in input is ignored.
-!   work       -  complex(kind(0.d0)), dimension (:), target.
+!   work       -  complex(psb_dpk_), dimension (:), target.
 !                 Workspace. Its size must be at least 4*psb_cd_get_local_cols(desc_data).
 !   info       -  integer, output.
 !                 Error code.
@ -143,16 +143,16 @@ subroutine mld_zsub_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
  ! Arguments
  type(psb_desc_type), intent(in)           :: desc_data
  type(mld_zbaseprc_type), intent(in)       :: prec
-  complex(kind(0.d0)),intent(in)            :: x(:)
+  complex(psb_dpk_),intent(in)            :: x(:)
-  complex(kind(0.d0)),intent(inout)         :: y(:)
+  complex(psb_dpk_),intent(inout)         :: y(:)
-  complex(kind(0.d0)),intent(in)            :: alpha,beta
+  complex(psb_dpk_),intent(in)            :: alpha,beta
  character(len=1), intent(in)              :: trans
-  complex(kind(0.d0)),target, intent(inout) :: work(:)
+  complex(psb_dpk_),target, intent(inout) :: work(:)
  integer, intent(out)                      :: info
  ! Local variables
  integer :: n_row,n_col
-  complex(kind(1.d0)), pointer :: ww(:), aux(:), tx(:),ty(:)
+  complex(psb_dpk_), pointer :: ww(:), aux(:), tx(:),ty(:)
  integer :: ictxt,np,me,i, err_act
  character(len=20)   :: name
  character           :: trans_
@ -186,7 +186,7 @@ subroutine mld_zsub_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
      if (info /= 0) then 
        info=4025
        call psb_errpush(info,name,i_err=(/4*n_col,0,0,0,0/),&
-             & a_err='complex(kind(1.d0))')
+             & a_err='complex(psb_dpk_)')
        goto 9999      
      end if
    endif
@ -195,7 +195,7 @@ subroutine mld_zsub_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
    if (info /= 0) then 
      info=4025
      call psb_errpush(info,name,i_err=(/5*n_col,0,0,0,0/),&
-           & a_err='complex(kind(1.d0))')
+           & a_err='complex(psb_dpk_)')
      goto 9999      
    end if
  endif
@ -225,7 +225,7 @@ subroutine mld_zsub_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
    if (info /= 0) then 
      info=4025
      call psb_errpush(info,name,i_err=(/2*n_col,0,0,0,0/),&
-           & a_err='complex(kind(1.d0))')
+           & a_err='complex(psb_dpk_)')
      goto 9999      
    end if
--- a/mlprec/mld_zsub_solve.f90
+++ b/mlprec/mld_zsub_solve.f90
@ -86,16 +86,16 @@
 !
 ! Arguments:
 !
-!   alpha      -  complex(kind(0.d0)), input.
+!   alpha      -  complex(psb_dpk_), input.
 !                 The scalar alpha.
 !   prec       -  type(mld_zbaseprec_type), input.
 !                 The 'base preconditioner' data structure containing the local 
 !                 part of the L and U factors of the matrix A.
-!   x          -  complex(kind(0.d0)), dimension(:), input.
+!   x          -  complex(psb_dpk_), dimension(:), input.
 !                 The local part of the vector X.
-!   beta       -  complex(kind(0.d0)), input.
+!   beta       -  complex(psb_dpk_), input.
 !                 The scalar beta.
-!   y          -  complex(kind(0.d0)), dimension(:), input/output.
+!   y          -  complex(psb_dpk_), dimension(:), input/output.
 !                 The local part of the vector Y.
 !   desc_data  -  type(psb_desc_type), input.
 !                 The communication descriptor associated to the matrix to be
@ -106,7 +106,7 @@
 !                 if trans='C','c' then op(K^(-1)) = K^(-C) (transpose conjugate of K^(-1)).
 !                 If prec%iprcparm(smooth_sweeps_) > 1, the value of trans provided
 !                 in input is ignored.
-!   work       -  complex(kind(0.d0)), dimension (:), target.
+!   work       -  complex(psb_dpk_), dimension (:), target.
 !                 Workspace. Its size must be at least 4*psb_cd_get_local_cols(desc_data).
 !   info       -  integer, output.
 !                 Error code.
@ -121,26 +121,27 @@ subroutine mld_zsub_solve(alpha,prec,x,beta,y,desc_data,trans,work,info)
  ! Arguments
  type(psb_desc_type), intent(in)           :: desc_data
  type(mld_zbaseprc_type), intent(in)       :: prec
-  complex(kind(0.d0)),intent(in)            :: x(:)
+  complex(psb_dpk_),intent(in)            :: x(:)
-  complex(kind(0.d0)),intent(inout)         :: y(:)
+  complex(psb_dpk_),intent(inout)         :: y(:)
-  complex(kind(0.d0)),intent(in)            :: alpha,beta
+  complex(psb_dpk_),intent(in)            :: alpha,beta
  character(len=1), intent(in)              :: trans
-  complex(kind(0.d0)),target, intent(inout) :: work(:)
+  complex(psb_dpk_),target, intent(inout) :: work(:)
  integer, intent(out)                      :: info
  ! Local variables
  integer :: n_row,n_col
-  complex(kind(1.d0)), pointer :: ww(:), aux(:), tx(:),ty(:)
+  complex(psb_dpk_), pointer :: ww(:), aux(:), tx(:),ty(:)
  integer :: ictxt,np,me,i, err_act
  character(len=20)   :: name
  character           :: trans_
  interface 
    subroutine mld_zumf_solve(flag,m,x,b,n,ptr,info)
      use psb_base_mod
      integer, intent(in)  :: flag,m,n,ptr
      integer, intent(out) :: info
-      complex(kind(1.d0)), intent(in)    :: b(*)
+      complex(psb_dpk_), intent(in)    :: b(*)
-      complex(kind(1.d0)), intent(inout) :: x(*)
+      complex(psb_dpk_), intent(inout) :: x(*)
    end subroutine mld_zumf_solve
  end interface
@ -173,7 +174,7 @@ subroutine mld_zsub_solve(alpha,prec,x,beta,y,desc_data,trans,work,info)
      if (info /= 0) then 
        info=4025
        call psb_errpush(info,name,i_err=(/4*n_col,0,0,0,0/),&
-             & a_err='complex(kind(1.d0))')
+             & a_err='complex(psb_dpk_)')
        goto 9999      
      end if
    endif
@ -182,7 +183,7 @@ subroutine mld_zsub_solve(alpha,prec,x,beta,y,desc_data,trans,work,info)
    if (info /= 0) then 
      info=4025
      call psb_errpush(info,name,i_err=(/5*n_col,0,0,0,0/),&
-           & a_err='complex(kind(1.d0))')
+           & a_err='complex(psb_dpk_)')
      goto 9999      
    end if
  endif
--- a/test/fileread/df_bench.f90
+++ b/test/fileread/df_bench.f90
@ -22,10 +22,10 @@ program df_bench
  ! dense matrices
-  real(kind(1.d0)), allocatable, target        ::  aux_b(:,:), d(:)
+  real(psb_dpk_), allocatable, target        ::  aux_b(:,:), d(:)
-  real(kind(1.d0)), allocatable , save :: b_col(:), x_col(:), r_col(:), &
+  real(psb_dpk_), allocatable , save :: b_col(:), x_col(:), r_col(:), &
       & x_col_glob(:), r_col_glob(:)
-  real(kind(1.d0)), pointer  :: b_col_glob(:)
+  real(psb_dpk_), pointer  :: b_col_glob(:)
  ! communications data structure
  type(psb_desc_type):: desc_a
@ -36,7 +36,7 @@ program df_bench
  ! solver paramters
  integer            :: iter, itmax, ierr, itrace, ircode, ipart,nlev,&
       & methd, istopc, iprec, ml, irnum, irst, ntry, nmat, ilev,ipsize,asize,cdsize
-  real(kind(1.d0))   :: err, eps
+  real(psb_dpk_)   :: err, eps
  character(len=5)   :: afmt
  character(len=20)  :: name
@ -45,7 +45,7 @@ program df_bench
  ! other variables
  integer            :: i,info,j,m_problem, nm, nt
  integer            :: internal, m,ii,nnzero, nprecs, pp
-  real(kind(1.d0)) :: t1, t2, tprec, r_amax, b_amax,&
+  real(psb_dpk_) :: t1, t2, tprec, r_amax, b_amax,&
       &scale,resmx,resmxp, mttot, mtslv, mtprec
  integer :: nrhs, nrow, n_row, dim, nv, ne
  integer, allocatable  :: ipv(:), neigh(:), ivg(:)
--- a/test/fileread/getp.f90
+++ b/test/fileread/getp.f90
@ -21,7 +21,7 @@ contains
    type(precdata),allocatable  :: precs(:)
    integer             :: iret, istopc,itmax,itrace,ipart,nmat,nprecs,irst,irnum,ntry
    character(len=1024) :: charbuf
-    real(kind(1.d0))    :: eps, omega,thr1,thr2
+    real(psb_dpk_)    :: eps, omega,thr1,thr2
    character           :: afmt*5, lv1*10, lv2*10, pdescr*40
    integer             :: iam, nm, np, i, idx
    integer, parameter  :: npparms=14
--- a/test/fileread/precdata.f90
+++ b/test/fileread/precdata.f90
@ -1,5 +1,5 @@
 module precd
-
+  use psb_base_mod, only : psb_dpk_
  type precdata
     character(len=10)  :: lv1, lv2    ! First and second level prec type
     integer            :: nlev        !
@ -8,7 +8,7 @@ module precd
     integer            :: prol        ! prolongation over application of as
     integer            :: ftype1      ! Factorization type: ILU, SuperLU, UMFPACK. 
     integer            :: fill1       ! Fill-in for factorization 1
-     real(kind(1.d0))   :: thr1        ! Threshold for fact. 1 ILU(T)
+     real(psb_dpk_)   :: thr1        ! Threshold for fact. 1 ILU(T)
     integer            :: mltype      ! additive or multiplicative 2nd level prec
     integer            :: aggr        ! local or global aggregation
     integer            :: smthkind    ! smoothing type
@ -17,9 +17,9 @@ module precd
     integer            :: glbsmth     ! global smoothing
     integer            :: ftype2      ! Factorization type: ILU, SuperLU, UMFPACK. 
     integer            :: fill2       ! Fill-in for factorization 1
-     real(kind(1.d0))   :: thr2        ! Threshold for fact. 1 ILU(T)
+     real(psb_dpk_)   :: thr2        ! Threshold for fact. 1 ILU(T)
     integer            :: jswp        ! Jacobi sweeps
-     real(kind(1.d0))   :: omega       ! smoother omega
+     real(psb_dpk_)   :: omega       ! smoother omega
     character(len=40)  :: descr       ! verbose description of the prec
  end type precdata
--- a/test/fileread/zf_bench.f90
+++ b/test/fileread/zf_bench.f90
@ -22,10 +22,10 @@ program zf_bench
  integer               :: igsmth, matop, novr
  ! dense matrices
-  complex(kind(1.d0)), allocatable, target        ::  aux_b(:,:), d(:)
+  complex(psb_dpk_), allocatable, target        ::  aux_b(:,:), d(:)
-  complex(kind(1.d0)), allocatable , save :: b_col(:), x_col(:), r_col(:), &
+  complex(psb_dpk_), allocatable , save :: b_col(:), x_col(:), r_col(:), &
       & x_col_glob(:), r_col_glob(:)
-  complex(kind(1.d0)), pointer  :: b_col_glob(:)
+  complex(psb_dpk_), pointer  :: b_col_glob(:)
  ! communications data structure
  type(psb_desc_type):: desc_a
@ -37,7 +37,7 @@ program zf_bench
  ! solver paramters
  integer            :: iter, itmax, ierr, itrace, ircode, ipart,nlev,&
       & methd, istopc, iprec, ml, irnum, irst, ntry, nmat, ilev,ipsize,asize,cdsize
-  real(kind(1.d0))   :: err, eps
+  real(psb_dpk_)   :: err, eps
  character(len=5)   :: afmt
  character(len=20)  :: name
@ -46,7 +46,7 @@ program zf_bench
  ! other variables
  integer            :: i,info,j,m_problem, nm, nt
  integer            :: internal, m,ii,nnzero, nprecs, pp
-  real(kind(1.d0)) :: t1, t2, tprec, r_amax, b_amax,&
+  real(psb_dpk_) :: t1, t2, tprec, r_amax, b_amax,&
       &scale,resmx,resmxp, mttot, mtslv, mtprec
  integer :: nrhs, nrow, n_row, dim, nv, ne
  integer, allocatable  :: ipv(:), neigh(:), ivg(:)
--- a/test/pargen/ppde.f90
+++ b/test/pargen/ppde.f90
@ -90,8 +90,8 @@ program ppde
  integer   :: idim
  ! miscellaneous 
-  real(kind(1.d0)), parameter :: one = 1.d0
+  real(psb_dpk_), parameter :: one = 1.d0
-  real(kind(1.d0)) :: t1, t2, tprec 
+  real(psb_dpk_) :: t1, t2, tprec 
  ! sparse matrix and preconditioner
  type(psb_dspmat_type) :: a
@ -99,13 +99,13 @@ program ppde
  ! descriptor
  type(psb_desc_type)   :: desc_a
  ! dense matrices
-  real(kind(1.d0)), allocatable :: b(:), x(:)
+  real(psb_dpk_), allocatable :: b(:), x(:)
  ! blacs parameters
  integer            :: ictxt, iam, np
  ! solver parameters
  integer            :: iter, itmax,itrace, istopc, irst
-  real(kind(1.d0))   :: err, eps
+  real(psb_dpk_)   :: err, eps
  type precdata
     character(len=10)  :: lv1, lvn    ! First level(s) and last  level prec type
@ -115,7 +115,7 @@ program ppde
     integer            :: prol        ! prolongation over application of as
     integer            :: ftype1      ! Factorization type: ILU, SuperLU, UMFPACK. 
     integer            :: fill1       ! Fill-in for factorization 1
-     real(kind(1.d0))   :: thr1        ! Threshold for fact. 1 ILU(T)
+     real(psb_dpk_)   :: thr1        ! Threshold for fact. 1 ILU(T)
     integer            :: mltype      ! additive or multiplicative 2nd level prec
     integer            :: aggr        ! local or global aggregation
     integer            :: smthkind    ! smoothing type
@ -124,9 +124,9 @@ program ppde
     integer            :: glbsmth     ! global smoothing
     integer            :: ftype2      ! Factorization type: ILU, SuperLU, UMFPACK. 
     integer            :: fill2       ! Fill-in for factorization 1
-     real(kind(1.d0))   :: thr2        ! Threshold for fact. 1 ILU(T)
+     real(psb_dpk_)   :: thr2        ! Threshold for fact. 1 ILU(T)
     integer            :: jswp        ! Jacobi sweeps
-     real(kind(1.d0))   :: omega       ! smoother omega
+     real(psb_dpk_)   :: omega       ! smoother omega
     character(len=20)  :: descr       ! verbose description of the prec
  end type precdata
  type(precdata)     :: prectype
@ -471,7 +471,7 @@ contains
    implicit none
    integer                        :: idim
    integer, parameter             :: nbmax=10
-    real(kind(1.d0)), allocatable  :: b(:),xv(:)
+    real(psb_dpk_), allocatable  :: b(:),xv(:)
    type(psb_desc_type)            :: desc_a
    integer                        :: ictxt, info
    character                      :: afmt*5
@ -485,21 +485,21 @@ contains
      end subroutine parts
    end interface   ! local variables
    type(psb_dspmat_type)    :: a
-    real(kind(1.d0))         :: zt(nbmax),glob_x,glob_y,glob_z
+    real(psb_dpk_)         :: zt(nbmax),glob_x,glob_y,glob_z
    integer                  :: m,n,nnz,glob_row
    integer                  :: x,y,z,ia,indx_owner
    integer                  :: np, iam
    integer                  :: element
    integer                  :: nv, inv
    integer, allocatable     :: irow(:),icol(:)
-    real(kind(1.d0)), allocatable :: val(:)
+    real(psb_dpk_), allocatable :: val(:)
    integer, allocatable     :: prv(:)
    ! deltah dimension of each grid cell
    ! deltat discretization time
-    real(kind(1.d0))         :: deltah
+    real(psb_dpk_)         :: deltah
-    real(kind(1.d0)),parameter   :: rhs=0.d0,one=1.d0,zero=0.d0
+    real(psb_dpk_),parameter   :: rhs=0.d0,one=1.d0,zero=0.d0
-    real(kind(1.d0))   :: t1, t2, t3, tins, tasb
+    real(psb_dpk_)   :: t1, t2, t3, tins, tasb
-    real(kind(1.d0))   :: a1, a2, a3, a4, b1, b2, b3 
+    real(psb_dpk_)   :: a1, a2, a3, a4, b1, b2, b3 
    external           :: a1, a2, a3, a4, b1, b2, b3
    integer            :: err_act
    ! common area
@ -750,38 +750,45 @@ end program ppde
 ! functions parametrizing the differential equation 
 !  
 function a1(x,y,z)
-  real(kind(1.d0)) :: a1
+  use psb_base_mod, only : psb_dpk_
-  real(kind(1.d0)) :: x,y,z
+  real(psb_dpk_) :: a1
  real(psb_dpk_) :: x,y,z
  a1=1.d0
 end function a1
 function a2(x,y,z)
-  real(kind(1.d0)) ::  a2
+  use psb_base_mod, only : psb_dpk_
-  real(kind(1.d0)) :: x,y,z
+  real(psb_dpk_) ::  a2
  real(psb_dpk_) :: x,y,z
  a2=2.d1*y
 end function a2
 function a3(x,y,z)
-  real(kind(1.d0)) ::  a3
+  use psb_base_mod, only : psb_dpk_
-  real(kind(1.d0)) :: x,y,z      
+  real(psb_dpk_) ::  a3
  real(psb_dpk_) :: x,y,z      
  a3=1.d0
 end function a3
 function a4(x,y,z)
-  real(kind(1.d0)) ::  a4
+  use psb_base_mod, only : psb_dpk_
-  real(kind(1.d0)) :: x,y,z      
+  real(psb_dpk_) ::  a4
  real(psb_dpk_) :: x,y,z      
  a4=1.d0
 end function a4
 function b1(x,y,z)
-  real(kind(1.d0)) ::  b1   
+  use psb_base_mod, only : psb_dpk_
-  real(kind(1.d0)) :: x,y,z
+  real(psb_dpk_) ::  b1   
  real(psb_dpk_) :: x,y,z
  b1=1.d0
 end function b1
 function b2(x,y,z)
-  real(kind(1.d0)) ::  b2
+  use psb_base_mod, only : psb_dpk_
-  real(kind(1.d0)) :: x,y,z
+  real(psb_dpk_) ::  b2
  real(psb_dpk_) :: x,y,z
  b2=1.d0
 end function b2
 function b3(x,y,z)
-  real(kind(1.d0)) ::  b3
+  use psb_base_mod, only : psb_dpk_
-  real(kind(1.d0)) :: x,y,z
+  real(psb_dpk_) ::  b3
  real(psb_dpk_) :: x,y,z
  b3=1.d0
 end function b3