diff --git a/mlprec/mld_daggrmat_raw_asb.F90 b/mlprec/mld_daggrmat_raw_asb.F90
index a28e70c6..791eba0d 100644
--- a/mlprec/mld_daggrmat_raw_asb.F90
+++ b/mlprec/mld_daggrmat_raw_asb.F90
@@ -99,6 +99,7 @@ subroutine mld_daggrmat_raw_asb(a,desc_a,ac,desc_ac,p,info)
   character(len=20) :: name
   type(psb_dspmat_type)  :: b
   integer, pointer :: nzbr(:), idisp(:)
+  type(psb_dspmat_type), pointer  :: am1,am2
   integer :: nrow, nglob, ncol, ntaggr, nzac, ip, ndx,&
        & naggr, nzt, naggrm1, i
 
@@ -116,6 +117,13 @@ subroutine mld_daggrmat_raw_asb(a,desc_a,ac,desc_ac,p,info)
   nrow  = psb_cd_get_local_rows(desc_a)
   ncol  = psb_cd_get_local_cols(desc_a)
 
+
+  am2 => p%av(mld_sm_pr_t_)
+  am1 => p%av(mld_sm_pr_)
+  call psb_nullify_sp(am1)
+  call psb_nullify_sp(am2)
+
+
   naggr  = p%nlaggr(me+1)
   ntaggr = sum(p%nlaggr)
   allocate(nzbr(np), idisp(np),stat=info)
@@ -140,6 +148,28 @@ subroutine mld_daggrmat_raw_asb(a,desc_a,ac,desc_ac,p,info)
     goto 9999
   end if
 
+  if (p%iprcparm(mld_coarse_mat_) == mld_repl_mat_) then
+    call psb_sp_all(ncol,ntaggr,am1,ncol,info)
+  else
+    call psb_sp_all(ncol,naggr,am1,ncol,info)
+  end if
+
+  if (info /= 0) then
+    call psb_errpush(4010,name,a_err='spall')
+    goto 9999
+  end if
+
+  do i=1,nrow
+    am1%aspk(i) = done
+    am1%ia1(i)  = i
+    am1%ia2(i)  = p%mlia(i)  
+  end do
+  am1%infoa(psb_nnz_) = nrow
+
+  call psb_spcnv(am1,info,afmt='csr',dupl=psb_dupl_add_)
+  call psb_transp(am1,am2)
+
+
   call psb_sp_clip(a,b,info,jmax=nrow)
   if(info /= 0) then
     call psb_errpush(4010,name,a_err='spclip')
diff --git a/mlprec/mld_dmlprec_aply.f90 b/mlprec/mld_dmlprec_aply.f90
index 58629841..08ad9686 100644
--- a/mlprec/mld_dmlprec_aply.f90
+++ b/mlprec/mld_dmlprec_aply.f90
@@ -436,8 +436,7 @@ contains
       nc2l  = psb_cd_get_local_cols(baseprecv(ilev)%base_desc)
       nr2l  = psb_cd_get_local_rows(baseprecv(ilev)%base_desc)
       allocate(mlprec_wrk(ilev)%x2l(nc2l),mlprec_wrk(ilev)%y2l(nc2l),&
-           & mlprec_wrk(ilev)%tx(max(n_row,n_col)),&
-           & mlprec_wrk(ilev)%ty(max(n_row,n_col)), stat=info)
+           & stat=info)
       if (info /= 0) then 
         info=4025
         call psb_errpush(info,name,i_err=(/2*(nc2l+max(n_row,n_col)),0,0,0,0/),&
@@ -445,33 +444,15 @@ contains
         goto 9999      
       end if
 
-      mlprec_wrk(ilev)%x2l(:) = dzero
-      mlprec_wrk(ilev)%y2l(:) = dzero
-      mlprec_wrk(ilev)%tx(1:n_row) = mlprec_wrk(ilev-1)%x2l(1:n_row) 
-      mlprec_wrk(ilev)%tx(n_row+1:max(n_row,n_col)) = dzero
-      mlprec_wrk(ilev)%ty(:) = dzero
-
       ismth = baseprecv(ilev)%iprcparm(mld_aggr_kind_)
       icm   = baseprecv(ilev)%iprcparm(mld_coarse_mat_)
-      if (ismth  /= mld_no_smooth_) then 
-        !
-        ! Apply the smoothed prolongator transpose
-        !
-        call psb_halo(mlprec_wrk(ilev-1)%x2l,baseprecv(ilev-1)%base_desc,&
-             &  info,work=work)         
-        if (info == 0) call psb_csmm(done,baseprecv(ilev)%av(mld_sm_pr_t_),&
-             & mlprec_wrk(ilev-1)%x2l,dzero,mlprec_wrk(ilev)%x2l,info)
-      else
-        !
-        ! Apply the raw aggregation map transpose (take a shortcut)
-        !
-        do i=1,n_row
-          mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) = &
-               &  mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) + &
-               &  mlprec_wrk(ilev-1)%x2l(i)
-        end do
+      
+      ! Apply prolongator transpose, i.e. restriction
+      call psb_forward_map(done,mlprec_wrk(ilev-1)%x2l,&
+           & dzero,mlprec_wrk(ilev)%x2l,&
+           & baseprecv(ilev)%map_desc,info,work=work)
 
-      end if
+      
       if (info /=0) then
         call psb_errpush(4001,name,a_err='Error during restriction')
         goto 9999
@@ -510,22 +491,13 @@ contains
       ismth = baseprecv(ilev)%iprcparm(mld_aggr_kind_)
       icm   = baseprecv(ilev)%iprcparm(mld_coarse_mat_)
 
-      if (ismth  /= mld_no_smooth_) then 
-        !
-        ! Apply the smoothed prolongator
-        !
-        call psb_csmm(done,baseprecv(ilev)%av(mld_sm_pr_),mlprec_wrk(ilev)%y2l,&
-             & done,mlprec_wrk(ilev-1)%y2l,info)
-      else
-        !
-        ! Apply the raw aggregation map (take a shortcut)
-        !
-        do i=1, n_row
-          mlprec_wrk(ilev-1)%y2l(i) = mlprec_wrk(ilev-1)%y2l(i) + &
-               &   mlprec_wrk(ilev)%y2l(baseprecv(ilev)%mlia(i))
-        enddo
+      !
+      ! Apply prolongator
+      !  
+      call psb_backward_map(done,mlprec_wrk(ilev)%y2l,&
+           & done,mlprec_wrk(ilev-1)%y2l,&
+           & baseprecv(ilev)%map_desc,info,work=work)
 
-      end if
       if (info /=0) then
         call psb_errpush(4001,name,a_err='Error during prolongation')
         goto 9999
@@ -704,9 +676,7 @@ contains
       goto 9999
     end if
 
-    mlprec_wrk(1)%y2l(:) = dzero
     mlprec_wrk(1)%x2l(:) = x
-
     !
     ! STEP 2
     !
@@ -758,29 +728,11 @@ contains
         goto 9999      
       end if
 
-      mlprec_wrk(ilev)%x2l(:) = dzero
-      mlprec_wrk(ilev)%y2l(:) = dzero
-      mlprec_wrk(ilev)%tx(:)  = dzero
-
-      if (ismth  /= mld_no_smooth_) then 
-        !
-        ! Apply the smoothed prolongator transpose
-        !
-        call psb_halo(mlprec_wrk(ilev-1)%tx,baseprecv(ilev-1)%base_desc,&
-             & info,work=work) 
-        if (info == 0) call psb_csmm(done,baseprecv(ilev)%av(mld_sm_pr_t_),&
-             & mlprec_wrk(ilev-1)%tx,dzero,mlprec_wrk(ilev)%x2l,info)
-      else
-        !
-        ! Apply the raw aggregation map transpose (take a shortcut)
-        !
-        mlprec_wrk(ilev)%x2l = dzero
-        do i=1,n_row
-          mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) = &
-               & mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) + &
-               &  mlprec_wrk(ilev-1)%tx(i)
-        end do
-      end if
+      ! Apply prolongator transpose, i.e. restriction      
+      call psb_forward_map(done,mlprec_wrk(ilev-1)%tx,&
+           & dzero,mlprec_wrk(ilev)%x2l,&
+           & baseprecv(ilev)%map_desc,info,work=work)
+      
       if (info /=0) then
         call psb_errpush(4001,name,a_err='Error during restriction')
         goto 9999
@@ -826,24 +778,13 @@ contains
       ismth = baseprecv(ilev+1)%iprcparm(mld_aggr_kind_)
       n_row = psb_cd_get_local_rows(baseprecv(ilev)%base_desc)
 
-      if (ismth  /= mld_no_smooth_) then 
-        !
-        ! Apply the smoothed prolongator
-        !
-        if (ismth == mld_smooth_prol_) &
-             & call psb_halo(mlprec_wrk(ilev+1)%y2l,&
-             & baseprecv(ilev+1)%base_desc,info,work=work) 
-        if (info == 0) call psb_csmm(done,baseprecv(ilev+1)%av(mld_sm_pr_),&
-             & mlprec_wrk(ilev+1)%y2l,done,mlprec_wrk(ilev)%y2l,info)
-      else
-        !
-        ! Apply the raw aggregation map (take a shortcut)
-        !
-        do i=1, n_row
-          mlprec_wrk(ilev)%y2l(i) = mlprec_wrk(ilev)%y2l(i) + &
-               & mlprec_wrk(ilev+1)%y2l(baseprecv(ilev+1)%mlia(i))
-        enddo
-      end if
+      !
+      ! Apply prolongator
+      !  
+      call psb_backward_map(done,mlprec_wrk(ilev+1)%y2l,&
+           & done,mlprec_wrk(ilev)%y2l,&
+           & baseprecv(ilev+1)%map_desc,info,work=work)
+
       if (info /=0) then
         call psb_errpush(4001,name,a_err='Error during prolongation')
         goto 9999
@@ -1010,9 +951,6 @@ contains
 
     allocate(mlprec_wrk(1)%x2l(nc2l),mlprec_wrk(1)%y2l(nc2l), &
          & mlprec_wrk(1)%tx(nc2l), stat=info)
-    mlprec_wrk(1)%x2l(:) = dzero
-    mlprec_wrk(1)%y2l(:) = dzero
-    mlprec_wrk(1)%tx(:)  = dzero
 
     call psb_geaxpby(done,x,dzero,mlprec_wrk(1)%tx,&
          & baseprecv(1)%base_desc,info)
@@ -1049,32 +987,11 @@ contains
         goto 9999      
       end if
 
-      mlprec_wrk(ilev)%x2l(:) = dzero
-      mlprec_wrk(ilev)%y2l(:) = dzero
-      mlprec_wrk(ilev)%tx(:)  = dzero
-
-      if (ismth  /= mld_no_smooth_) then 
-        !
-        ! Apply the smoothed prolongator transpose
-        !
-        if (debug_level >= psb_debug_inner_) &
-             & write(debug_unit,*) me,' ',trim(name), ' up sweep ', ilev
-
-        call psb_halo(mlprec_wrk(ilev-1)%x2l,&
-             &  baseprecv(ilev-1)%base_desc,info,work=work) 
-        if (info == 0) call psb_csmm(done,baseprecv(ilev)%av(mld_sm_pr_t_),&
-             & mlprec_wrk(ilev-1)%x2l,dzero,mlprec_wrk(ilev)%x2l,info)
-      else
-        !
-        ! Apply the raw aggregation map transpose (take a shortcut)
-        !
-        do i=1,n_row
-          mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) = &
-               & mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) + &
-               & mlprec_wrk(ilev-1)%x2l(i)
-        end do
-
-      end if
+      ! Apply prolongator transpose, i.e. restriction
+      call psb_forward_map(done,mlprec_wrk(ilev-1)%x2l,&
+           & dzero,mlprec_wrk(ilev)%x2l,&
+           & baseprecv(ilev)%map_desc,info,work=work)
+      
       if (info /=0) then
         call psb_errpush(4001,name,a_err='Error during restriction')
         goto 9999
@@ -1135,25 +1052,13 @@ contains
       ismth = baseprecv(ilev+1)%iprcparm(mld_aggr_kind_)
       n_row = psb_cd_get_local_rows(baseprecv(ilev)%base_desc)
 
-      if (ismth  /= mld_no_smooth_) then
-        !
-        ! Apply the smoothed prolongator
-        !  
-        if (ismth == mld_smooth_prol_) &
-             & call psb_halo(mlprec_wrk(ilev+1)%y2l,baseprecv(ilev+1)%base_desc,&
-             &  info,work=work) 
-        if (info == 0) call psb_csmm(done,baseprecv(ilev+1)%av(mld_sm_pr_),&
-             & mlprec_wrk(ilev+1)%y2l, dzero,mlprec_wrk(ilev)%y2l,info)
-      else
-        !
-        ! Apply the raw aggregation map (take a shortcut)
-        !
-        mlprec_wrk(ilev)%y2l(:) = dzero
-        do i=1, n_row
-          mlprec_wrk(ilev)%y2l(i) = mlprec_wrk(ilev)%y2l(i) + &
-               & mlprec_wrk(ilev+1)%y2l(baseprecv(ilev+1)%mlia(i))
-        enddo
-      end if
+      !
+      ! Apply prolongator
+      !  
+      call psb_backward_map(done,mlprec_wrk(ilev+1)%y2l,&
+           & dzero,mlprec_wrk(ilev)%y2l,&
+           & baseprecv(ilev+1)%map_desc,info,work=work)
+
       if (info /=0) then
         call psb_errpush(4001,name,a_err='Error during prolongation')
         goto 9999
@@ -1358,10 +1263,6 @@ contains
            & a_err='real(kind(1.d0))')
       goto 9999
     end if
-    mlprec_wrk(1)%x2l(:) = dzero
-    mlprec_wrk(1)%y2l(:) = dzero
-    mlprec_wrk(1)%tx(:)  = dzero
-    mlprec_wrk(1)%ty(:)  = dzero
 
     call psb_geaxpby(done,x,dzero,mlprec_wrk(1)%x2l,&
          & baseprecv(1)%base_desc,info)
@@ -1403,8 +1304,8 @@ contains
       nr2l  = psb_cd_get_local_rows(baseprecv(ilev)%base_desc)
       ismth = baseprecv(ilev)%iprcparm(mld_aggr_kind_)
       icm   = baseprecv(ilev)%iprcparm(mld_coarse_mat_)
-      allocate(mlprec_wrk(ilev)%ty(nc2l),mlprec_wrk(ilev)%y2l(nc2l),&
-           &   mlprec_wrk(ilev)%x2l(nc2l), stat=info)
+      allocate(mlprec_wrk(ilev)%tx(nc2l),mlprec_wrk(ilev)%ty(nc2l),&
+           &  mlprec_wrk(ilev)%y2l(nc2l),mlprec_wrk(ilev)%x2l(nc2l), stat=info)
 
       if (info /= 0) then 
         info=4025
@@ -1413,30 +1314,11 @@ contains
         goto 9999      
       end if
 
-      mlprec_wrk(ilev)%x2l(:) = dzero
-      mlprec_wrk(ilev)%y2l(:) = dzero
-      mlprec_wrk(ilev)%tx(:)  = dzero
-      mlprec_wrk(ilev)%ty(:)  = dzero
-
-      if (ismth  /= mld_no_smooth_) then 
-        !
-        ! Apply the smoothed prolongator transpose
-        !
-        call psb_halo(mlprec_wrk(ilev-1)%ty,baseprecv(ilev-1)%base_desc,&
-             & info,work=work) 
-        if (info == 0) call psb_csmm(done,baseprecv(ilev)%av(mld_sm_pr_t_),&
-             & mlprec_wrk(ilev-1)%ty,dzero,mlprec_wrk(ilev)%x2l,info)
-      else
-        !
-        ! Apply the raw aggregation map transpose (take a shortcut)
-        !
-        mlprec_wrk(ilev)%x2l = dzero
-        do i=1,n_row
-          mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) = &
-               & mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) + &
-               &  mlprec_wrk(ilev-1)%ty(i)
-        end do
-      end if
+      ! Apply prolongator transpose, i.e. restriction
+      call psb_forward_map(done,mlprec_wrk(ilev-1)%ty,&
+           & dzero,mlprec_wrk(ilev)%x2l,&
+           & baseprecv(ilev)%map_desc,info,work=work)
+      
       if (info /=0) then
         call psb_errpush(4001,name,a_err='Error during restriction')
         goto 9999
@@ -1485,24 +1367,13 @@ contains
       ismth = baseprecv(ilev+1)%iprcparm(mld_aggr_kind_)
       n_row = psb_cd_get_local_rows(baseprecv(ilev)%base_desc)
 
-      if (ismth  /= mld_no_smooth_) then
-        !
-        ! Apply the smoothed prolongator
-        ! 
-        if (ismth == mld_smooth_prol_) &
-             & call psb_halo(mlprec_wrk(ilev+1)%y2l,baseprecv(ilev+1)%base_desc,&
-             &  info,work=work) 
-        if (info == 0) call psb_csmm(done,baseprecv(ilev+1)%av(mld_sm_pr_),&
-             & mlprec_wrk(ilev+1)%y2l, done,mlprec_wrk(ilev)%y2l,info)
-      else
-        !
-        ! Apply the raw aggregation map (take a shortcut)
-        !
-        do i=1, n_row
-          mlprec_wrk(ilev)%y2l(i) = mlprec_wrk(ilev)%y2l(i) + &
-               & mlprec_wrk(ilev+1)%y2l(baseprecv(ilev+1)%mlia(i))
-        enddo
-      end if
+      !
+      ! Apply prolongator
+      !  
+      call psb_backward_map(done,mlprec_wrk(ilev+1)%y2l,&
+           & done,mlprec_wrk(ilev)%y2l,&
+           & baseprecv(ilev+1)%map_desc,info,work=work)
+
       if (info /=0 ) then
         call psb_errpush(4001,name,a_err='Error during restriction')
         goto 9999
diff --git a/mlprec/mld_dmlprec_bld.f90 b/mlprec/mld_dmlprec_bld.f90
index 1012f25d..2cc2efd0 100644
--- a/mlprec/mld_dmlprec_bld.f90
+++ b/mlprec/mld_dmlprec_bld.f90
@@ -153,6 +153,11 @@ subroutine mld_dmlprec_bld(a,desc_a,p,info)
   p%base_a => p%av(mld_ac_)
   if (info==0) call psb_cdtransfer(desc_ac,p%desc_ac,info)
 
+  p%map_desc = psb_inter_desc(psb_map_aggr_,desc_a,&
+       & p%desc_ac,p%av(mld_sm_pr_t_),p%av(mld_sm_pr_))
+  ! The two matrices from p%av() have been copied, may free them.
+  if (info == 0) call psb_sp_free(p%av(mld_sm_pr_t_),info)
+  if (info == 0) call psb_sp_free(p%av(mld_sm_pr_),info)
   if (info /= 0) then 
     call psb_errpush(4010,name,a_err='psb_cdtransfer')
     goto 9999
diff --git a/mlprec/mld_prec_type.f90 b/mlprec/mld_prec_type.f90
index 45a5664f..83c05950 100644
--- a/mlprec/mld_prec_type.f90
+++ b/mlprec/mld_prec_type.f90
@@ -61,7 +61,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_sizeof
+       & psb_inter_desc_type, psb_sizeof
 
   !
   ! Type: mld_dprec_type, mld_zprec_type
@@ -168,7 +168,7 @@ module mld_prec_type
     type(psb_dspmat_type), pointer     :: base_a    => null() 
     type(psb_desc_type), pointer       :: base_desc => null() 
     real(kind(1.d0)), allocatable      :: dorig(:) 
-
+    type(psb_inter_desc_type)          :: map_desc
   end type mld_dbaseprc_type
 
   type mld_dprec_type
@@ -187,7 +187,7 @@ module mld_prec_type
     type(psb_zspmat_type), pointer     :: base_a    => null() 
     type(psb_desc_type), pointer       :: base_desc => null() 
     complex(kind(1.d0)), allocatable   :: dorig(:)
-
+    type(psb_inter_desc_type)          :: map_desc
   end type mld_zbaseprc_type
 
   type mld_zprec_type
@@ -417,6 +417,7 @@ contains
         val = val + psb_sizeof(prec%av(i))
       end do
     end if
+    val = val + psb_sizeof(prec%map_desc) 
 
     mld_dbaseprc_sizeof = val 
     
@@ -453,6 +454,7 @@ contains
         val = val + psb_sizeof(prec%av(i))
       end do
     end if
+    val = val + psb_sizeof(prec%map_desc) 
     
     mld_zbaseprc_sizeof = val 
     
diff --git a/mlprec/mld_zaggrmat_raw_asb.F90 b/mlprec/mld_zaggrmat_raw_asb.F90
index 9ee8880f..7c4b65c1 100644
--- a/mlprec/mld_zaggrmat_raw_asb.F90
+++ b/mlprec/mld_zaggrmat_raw_asb.F90
@@ -99,6 +99,7 @@ subroutine mld_zaggrmat_raw_asb(a,desc_a,ac,desc_ac,p,info)
   character(len=20) :: name
   type(psb_zspmat_type)  :: b
   integer, pointer :: nzbr(:), idisp(:)
+  type(psb_zspmat_type), pointer  :: am1,am2
   integer :: nrow, nglob, ncol, ntaggr, nzac, ip, ndx,&
        & naggr, nzt,naggrm1, i
 
@@ -116,6 +117,13 @@ subroutine mld_zaggrmat_raw_asb(a,desc_a,ac,desc_ac,p,info)
   nrow  = psb_cd_get_local_rows(desc_a)
   ncol  = psb_cd_get_local_cols(desc_a)
 
+
+  am2 => p%av(mld_sm_pr_t_)
+  am1 => p%av(mld_sm_pr_)
+  call psb_nullify_sp(am1)
+  call psb_nullify_sp(am2)
+
+
   naggr  = p%nlaggr(me+1)
   ntaggr = sum(p%nlaggr)
   allocate(nzbr(np), idisp(np),stat=info)
@@ -140,6 +148,28 @@ subroutine mld_zaggrmat_raw_asb(a,desc_a,ac,desc_ac,p,info)
     goto 9999
   end if
 
+  if (p%iprcparm(mld_coarse_mat_) == mld_repl_mat_) then
+    call psb_sp_all(ncol,ntaggr,am1,ncol,info)
+  else
+    call psb_sp_all(ncol,naggr,am1,ncol,info)
+  end if
+
+  if (info /= 0) then
+    call psb_errpush(4010,name,a_err='spall')
+    goto 9999
+  end if
+
+  do i=1,nrow
+    am1%aspk(i) = zone
+    am1%ia1(i)  = i
+    am1%ia2(i)  = p%mlia(i)  
+  end do
+  am1%infoa(psb_nnz_) = nrow
+
+  call psb_spcnv(am1,info,afmt='csr',dupl=psb_dupl_add_)
+  call psb_transp(am1,am2)
+
+
   call psb_sp_clip(a,b,info,jmax=nrow)
   if(info /= 0) then
     call psb_errpush(4010,name,a_err='spclip')
diff --git a/mlprec/mld_zmlprec_aply.f90 b/mlprec/mld_zmlprec_aply.f90
index f4c5ab7f..d9f310b2 100644
--- a/mlprec/mld_zmlprec_aply.f90
+++ b/mlprec/mld_zmlprec_aply.f90
@@ -228,7 +228,7 @@ subroutine mld_zmlprec_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
     select case(baseprecv(2)%iprcparm(mld_smooth_pos_))
 
     case(mld_post_smooth_)
-      
+
       select case (trans_) 
       case('N')
         call mlt_post_ml_aply(alpha,baseprecv,x,beta,y,desc_data,trans_,work,info)
@@ -438,8 +438,7 @@ contains
       nc2l  = psb_cd_get_local_cols(baseprecv(ilev)%base_desc)
       nr2l  = psb_cd_get_local_rows(baseprecv(ilev)%base_desc)
       allocate(mlprec_wrk(ilev)%x2l(nc2l),mlprec_wrk(ilev)%y2l(nc2l),&
-           & mlprec_wrk(ilev)%tx(max(n_row,n_col)),&
-           & mlprec_wrk(ilev)%ty(max(n_row,n_col)), stat=info)
+           & stat=info)
       if (info /= 0) then 
         info=4025
         call psb_errpush(info,name,i_err=(/2*(nc2l+max(n_row,n_col)),0,0,0,0/),&
@@ -447,33 +446,15 @@ contains
         goto 9999      
       end if
 
-      mlprec_wrk(ilev)%x2l(:) = zzero
-      mlprec_wrk(ilev)%y2l(:) = zzero
-      mlprec_wrk(ilev)%tx(1:n_row) = mlprec_wrk(ilev-1)%x2l(1:n_row) 
-      mlprec_wrk(ilev)%tx(n_row+1:max(n_row,n_col)) = zzero
-      mlprec_wrk(ilev)%ty(:) = zzero
-
       ismth = baseprecv(ilev)%iprcparm(mld_aggr_kind_)
       icm   = baseprecv(ilev)%iprcparm(mld_coarse_mat_)
-      if (ismth  /= mld_no_smooth_) then 
-        !
-        ! Apply the smoothed prolongator transpose
-        !
-        call psb_halo(mlprec_wrk(ilev-1)%x2l,baseprecv(ilev-1)%base_desc,&
-             &  info,work=work) 
-        if (info == 0)  call psb_csmm(zone,baseprecv(ilev)%av(mld_sm_pr_t_),&
-             & mlprec_wrk(ilev-1)%x2l,zzero,mlprec_wrk(ilev)%x2l,info)
-      else
-        !
-        ! Apply the raw aggregation map transpose (take a shortcut)
-        !
-        do i=1,n_row
-          mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) = &
-               &  mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) + &
-               &  mlprec_wrk(ilev-1)%x2l(i)
-        end do
+      
+      ! Apply prolongator transpose, i.e. restriction
+      call psb_forward_map(zone,mlprec_wrk(ilev-1)%x2l,&
+           & zzero,mlprec_wrk(ilev)%x2l,&
+           & baseprecv(ilev)%map_desc,info,work=work)
 
-      end if
+      
       if (info /=0) then
         call psb_errpush(4001,name,a_err='Error during restriction')
         goto 9999
@@ -512,22 +493,13 @@ contains
       ismth = baseprecv(ilev)%iprcparm(mld_aggr_kind_)
       icm   = baseprecv(ilev)%iprcparm(mld_coarse_mat_)
 
-      if (ismth  /= mld_no_smooth_) then 
-        !
-        ! Apply the smoothed prolongator
-        !
-        call psb_csmm(zone,baseprecv(ilev)%av(mld_sm_pr_),mlprec_wrk(ilev)%y2l,&
-             & zone,mlprec_wrk(ilev-1)%y2l,info)
-      else
-        !
-        ! Apply the raw aggregation map (take a shortcut)
-        !
-        do i=1, n_row
-          mlprec_wrk(ilev-1)%y2l(i) = mlprec_wrk(ilev-1)%y2l(i) + &
-               &   mlprec_wrk(ilev)%y2l(baseprecv(ilev)%mlia(i))
-        enddo
+      !
+      ! Apply prolongator
+      !  
+      call psb_backward_map(zone,mlprec_wrk(ilev)%y2l,&
+           & zone,mlprec_wrk(ilev-1)%y2l,&
+           & baseprecv(ilev)%map_desc,info,work=work)
 
-      end if
       if (info /=0) then
         call psb_errpush(4001,name,a_err='Error during prolongation')
         goto 9999
@@ -707,9 +679,7 @@ contains
       goto 9999
     end if
 
-    mlprec_wrk(1)%y2l(:) = zzero
     mlprec_wrk(1)%x2l(:) = x
-
     !
     ! STEP 2
     !
@@ -761,29 +731,11 @@ contains
         goto 9999      
       end if
 
-      mlprec_wrk(ilev)%x2l(:) = zzero
-      mlprec_wrk(ilev)%y2l(:) = zzero
-      mlprec_wrk(ilev)%tx(:)  = zzero
-
-      if (ismth  /= mld_no_smooth_) then 
-        !
-        ! Apply the smoothed prolongator transpose
-        !
-        call psb_halo(mlprec_wrk(ilev-1)%tx,baseprecv(ilev-1)%base_desc,&
-             & info,work=work) 
-        if (info == 0) call psb_csmm(zone,baseprecv(ilev)%av(mld_sm_pr_t_),&
-             & mlprec_wrk(ilev-1)%tx,zzero,mlprec_wrk(ilev)%x2l,info)
-      else
-        !
-        ! Apply the raw aggregation map transpose (take a shortcut)
-        !
-        mlprec_wrk(ilev)%x2l = zzero
-        do i=1,n_row
-          mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) = &
-               & mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) + &
-               &  mlprec_wrk(ilev-1)%tx(i)
-        end do
-      end if
+      ! Apply prolongator transpose, i.e. restriction      
+      call psb_forward_map(zone,mlprec_wrk(ilev-1)%tx,&
+           & zzero,mlprec_wrk(ilev)%x2l,&
+           & baseprecv(ilev)%map_desc,info,work=work)
+      
       if (info /=0) then
         call psb_errpush(4001,name,a_err='Error during restriction')
         goto 9999
@@ -829,24 +781,13 @@ contains
       ismth = baseprecv(ilev+1)%iprcparm(mld_aggr_kind_)
       n_row = psb_cd_get_local_rows(baseprecv(ilev)%base_desc)
 
-      if (ismth  /= mld_no_smooth_) then 
-        !
-        ! Apply the smoothed prolongator
-        !
-        if (ismth == mld_smooth_prol_) &
-             & call psb_halo(mlprec_wrk(ilev+1)%y2l,&
-             & baseprecv(ilev+1)%base_desc,info,work=work) 
-        if (info == 0) call psb_csmm(zone,baseprecv(ilev+1)%av(mld_sm_pr_),&
-             & mlprec_wrk(ilev+1)%y2l,zone,mlprec_wrk(ilev)%y2l,info)
-      else
-        !
-        ! Apply the raw aggregation map (take a shortcut)
-        !
-        do i=1, n_row
-          mlprec_wrk(ilev)%y2l(i) = mlprec_wrk(ilev)%y2l(i) + &
-               & mlprec_wrk(ilev+1)%y2l(baseprecv(ilev+1)%mlia(i))
-        enddo
-      end if
+      !
+      ! Apply prolongator
+      !  
+      call psb_backward_map(zone,mlprec_wrk(ilev+1)%y2l,&
+           & zone,mlprec_wrk(ilev)%y2l,&
+           & baseprecv(ilev+1)%map_desc,info,work=work)
+
       if (info /=0) then
         call psb_errpush(4001,name,a_err='Error during prolongation')
         goto 9999
@@ -1014,9 +955,6 @@ contains
 
     allocate(mlprec_wrk(1)%x2l(nc2l),mlprec_wrk(1)%y2l(nc2l), &
          & mlprec_wrk(1)%tx(nc2l), stat=info)
-    mlprec_wrk(1)%x2l(:) = zzero
-    mlprec_wrk(1)%y2l(:) = zzero
-    mlprec_wrk(1)%tx(:)  = zzero
 
     call psb_geaxpby(zone,x,zzero,mlprec_wrk(1)%tx,&
          & baseprecv(1)%base_desc,info)
@@ -1053,31 +991,11 @@ contains
         goto 9999      
       end if
 
-      mlprec_wrk(ilev)%x2l(:) = zzero
-      mlprec_wrk(ilev)%y2l(:) = zzero
-      mlprec_wrk(ilev)%tx(:)  = zzero
-      if (ismth  /= mld_no_smooth_) then 
-        !
-        ! Apply the smoothed prolongator transpose
-        !
-        if (debug_level >= psb_debug_inner_) &
-             & write(debug_unit,*) me,' ',trim(name), ' up sweep ', ilev
-
-        call psb_halo(mlprec_wrk(ilev-1)%x2l,&
-             &  baseprecv(ilev-1)%base_desc,info,work=work) 
-        if (info == 0) call psb_csmm(zone,baseprecv(ilev)%av(mld_sm_pr_t_),&
-             & mlprec_wrk(ilev-1)%x2l,zzero,mlprec_wrk(ilev)%x2l,info)
-      else
-        !
-        ! Apply the raw aggregation map transpose (take a shortcut)
-        !
-        do i=1,n_row
-          mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) = &
-               & mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) + &
-               & mlprec_wrk(ilev-1)%x2l(i)
-        end do
-
-      end if
+      ! Apply prolongator transpose, i.e. restriction
+      call psb_forward_map(zone,mlprec_wrk(ilev-1)%x2l,&
+           & zzero,mlprec_wrk(ilev)%x2l,&
+           & baseprecv(ilev)%map_desc,info,work=work)
+      
       if (info /=0) then
         call psb_errpush(4001,name,a_err='Error during restriction')
         goto 9999
@@ -1137,25 +1055,13 @@ contains
       ismth = baseprecv(ilev+1)%iprcparm(mld_aggr_kind_)
       n_row = psb_cd_get_local_rows(baseprecv(ilev)%base_desc)
 
-      if (ismth  /= mld_no_smooth_) then
-        !
-        ! Apply the smoothed prolongator
-        !  
-        if (ismth == mld_smooth_prol_) &
-             & call psb_halo(mlprec_wrk(ilev+1)%y2l,baseprecv(ilev+1)%base_desc,&
-             &  info,work=work) 
-        if (info == 0) call psb_csmm(zone,baseprecv(ilev+1)%av(mld_sm_pr_),&
-             & mlprec_wrk(ilev+1)%y2l, zzero,mlprec_wrk(ilev)%y2l,info)
-      else
-        !
-        ! Apply the raw aggregation map (take a shortcut)
-        !
-        mlprec_wrk(ilev)%y2l(:) = zzero
-        do i=1, n_row
-          mlprec_wrk(ilev)%y2l(i) = mlprec_wrk(ilev)%y2l(i) + &
-               & mlprec_wrk(ilev+1)%y2l(baseprecv(ilev+1)%mlia(i))
-        enddo
-      end if
+      !
+      ! Apply prolongator
+      !  
+      call psb_backward_map(zone,mlprec_wrk(ilev+1)%y2l,&
+           & zzero,mlprec_wrk(ilev)%y2l,&
+           & baseprecv(ilev+1)%map_desc,info,work=work)
+
       if (info /=0) then
         call psb_errpush(4001,name,a_err='Error during prolongation')
         goto 9999
@@ -1361,11 +1267,6 @@ contains
            & a_err='real(kind(1.d0))')
       goto 9999
     end if
-    mlprec_wrk(1)%x2l(:) = zzero
-    mlprec_wrk(1)%y2l(:) = zzero
-    mlprec_wrk(1)%tx(:)  = zzero
-    mlprec_wrk(1)%ty(:)  = zzero
-
 
     call psb_geaxpby(zone,x,zzero,mlprec_wrk(1)%x2l,&
          & baseprecv(1)%base_desc,info)
@@ -1407,8 +1308,8 @@ contains
       nr2l  = psb_cd_get_local_rows(baseprecv(ilev)%base_desc)
       ismth = baseprecv(ilev)%iprcparm(mld_aggr_kind_)
       icm   = baseprecv(ilev)%iprcparm(mld_coarse_mat_)
-      allocate(mlprec_wrk(ilev)%ty(nc2l),mlprec_wrk(ilev)%y2l(nc2l),&
-           &   mlprec_wrk(ilev)%x2l(nc2l), stat=info)
+      allocate(mlprec_wrk(ilev)%tx(nc2l),mlprec_wrk(ilev)%ty(nc2l),&
+           &  mlprec_wrk(ilev)%y2l(nc2l),mlprec_wrk(ilev)%x2l(nc2l), stat=info)
 
       if (info /= 0) then 
         info=4025
@@ -1417,30 +1318,11 @@ contains
         goto 9999      
       end if
 
-      mlprec_wrk(ilev)%x2l(:) = zzero
-      mlprec_wrk(ilev)%y2l(:) = zzero
-      mlprec_wrk(ilev)%tx(:)  = zzero
-      mlprec_wrk(ilev)%ty(:)  = zzero
-
-      if (ismth  /= mld_no_smooth_) then 
-        !
-        ! Apply the smoothed prolongator transpose
-        !
-        call psb_halo(mlprec_wrk(ilev-1)%ty,baseprecv(ilev-1)%base_desc,&
-             & info,work=work) 
-        if (info == 0) call psb_csmm(zone,baseprecv(ilev)%av(mld_sm_pr_t_),&
-             & mlprec_wrk(ilev-1)%ty,zzero,mlprec_wrk(ilev)%x2l,info)
-      else
-        !
-        ! Apply the raw aggregation map transpose (take a shortcut)
-        !
-        mlprec_wrk(ilev)%x2l = zzero
-        do i=1,n_row
-          mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) = &
-               & mlprec_wrk(ilev)%x2l(baseprecv(ilev)%mlia(i)) + &
-               &  mlprec_wrk(ilev-1)%ty(i)
-        end do
-      end if
+      ! Apply prolongator transpose, i.e. restriction
+      call psb_forward_map(zone,mlprec_wrk(ilev-1)%ty,&
+           & zzero,mlprec_wrk(ilev)%x2l,&
+           & baseprecv(ilev)%map_desc,info,work=work)
+      
       if (info /=0) then
         call psb_errpush(4001,name,a_err='Error during restriction')
         goto 9999
@@ -1489,24 +1371,13 @@ contains
       ismth = baseprecv(ilev+1)%iprcparm(mld_aggr_kind_)
       n_row = psb_cd_get_local_rows(baseprecv(ilev)%base_desc)
 
-      if (ismth  /= mld_no_smooth_) then
-        !
-        ! Apply the smoothed prolongator
-        ! 
-        if (ismth == mld_smooth_prol_) &
-             & call psb_halo(mlprec_wrk(ilev+1)%y2l,baseprecv(ilev+1)%base_desc,&
-             &  info,work=work) 
-        if (info == 0) call psb_csmm(zone,baseprecv(ilev+1)%av(mld_sm_pr_),&
-             & mlprec_wrk(ilev+1)%y2l,zone,mlprec_wrk(ilev)%y2l,info)
-      else
-        !
-        ! Apply the raw aggregation map (take a shortcut)
-        !
-        do i=1, n_row
-          mlprec_wrk(ilev)%y2l(i) = mlprec_wrk(ilev)%y2l(i) + &
-               & mlprec_wrk(ilev+1)%y2l(baseprecv(ilev+1)%mlia(i))
-        enddo
-      end if
+      !
+      ! Apply prolongator
+      !  
+      call psb_backward_map(zone,mlprec_wrk(ilev+1)%y2l,&
+           & zone,mlprec_wrk(ilev)%y2l,&
+           & baseprecv(ilev+1)%map_desc,info,work=work)
+
       if (info /=0 ) then
         call psb_errpush(4001,name,a_err='Error during restriction')
         goto 9999
diff --git a/mlprec/mld_zmlprec_bld.f90 b/mlprec/mld_zmlprec_bld.f90
index ffaeda90..105d8bbd 100644
--- a/mlprec/mld_zmlprec_bld.f90
+++ b/mlprec/mld_zmlprec_bld.f90
@@ -153,6 +153,11 @@ subroutine mld_zmlprec_bld(a,desc_a,p,info)
   p%base_a => p%av(mld_ac_)
   if (info==0) call psb_cdtransfer(desc_ac,p%desc_ac,info)
 
+  p%map_desc = psb_inter_desc(psb_map_aggr_,desc_a,&
+       & p%desc_ac,p%av(mld_sm_pr_t_),p%av(mld_sm_pr_))
+  ! The two matrices from p%av() have been copied, may free them.
+  if (info == 0) call psb_sp_free(p%av(mld_sm_pr_t_),info)
+  if (info == 0) call psb_sp_free(p%av(mld_sm_pr_),info)
   if (info /= 0) then 
     call psb_errpush(4010,name,a_err='psb_cdtransfer')
     goto 9999