diff --git a/cbind/Makefile b/cbind/Makefile
new file mode 100644
index 00000000..95202fe9
--- /dev/null
+++ b/cbind/Makefile
@@ -0,0 +1,27 @@
+include ../Make.inc
+
+HERE=.
+LIBDIR=../lib
+INCDIR=../include
+MODDIR=../modules/
+LIBNAME=$(CBINDLIBNAME)
+LIBNAME=libmld_cbind.a	
+
+lib: mlprecd
+#	/bin/cp -p $(CPUPDFLAG) $(HERE)/$(LIBNAME) $(LIBDIR)#
+#	/bin/cp -p $(CPUPDFLAG) *.h $(INCDIR)
+#	/bin/cp -p $(CPUPDFLAG) *$(.mod) $(MODDIR)	
+
+
+mlprecd:
+	cd mlprec && $(MAKE) lib LIBNAME=$(LIBNAME) 
+
+
+clean: 
+	cd mlprec && $(MAKE) clean
+	cd krylov && $(MAKE) clean
+
+
+veryclean: clean
+	cd test/pargen && $(MAKE) clean
+	/bin/rm -f $(HERE)/$(LIBNAME) $(LIBMOD) *$(.mod) *.h
diff --git a/cbind/mlprec/Makefile b/cbind/mlprec/Makefile
new file mode 100644
index 00000000..d51a6316
--- /dev/null
+++ b/cbind/mlprec/Makefile
@@ -0,0 +1,34 @@
+include ../../Make.inc
+LIBDIR=../../lib
+INCDIR=../../include
+MODDIR=../../modules
+HERE=.
+
+CINCLUDES=-I. -I$(LIBDIR) -I$(PSBLAS_INCDIR)
+FINCLUDES=$(FMFLAG)$(HERE) $(FMFLAG)$(INCDIR) $(FMFLAG)$(MODDIR) $(PSBLAS_INCLUDES)
+
+
+OBJS=mld_prec_cbind_mod.o mld_dprec_cbind_mod.o  mld_c_dprec.o
+CMOD=mld_cbind.h mld_c_dprec.h mld_const.h
+
+
+LIBMOD=mld_prec_cbind_mod$(.mod) 
+LOCAL_MODS=$(LIBMOD)
+#LIBNAME=$(CPRECLIBNAME)	
+
+
+lib:  $(OBJS) $(CMOD)
+	$(AR) $(HERE)/$(LIBNAME) $(OBJS)
+	$(RANLIB) $(HERE)/$(LIBNAME)
+	/bin/cp -p $(HERE)/$(LIBNAME) $(LIBDIR)
+	/bin/cp -p $(LIBMOD) $(CMOD) $(INCDIR)
+
+mld_prec_cbind_mod.o: mld_dprec_cbind_mod.o 
+#mld_prec_cbind_mod.o: psb_prec_cbind_mod.o
+veryclean: clean
+	/bin/rm -f $(HERE)/$(LIBNAME)
+
+clean:
+	/bin/rm -f $(OBJS) $(LOCAL_MODS)
+
+veryclean: clean
diff --git a/cbind/mlprec/mld_c_dprec.c b/cbind/mlprec/mld_c_dprec.c
new file mode 100644
index 00000000..bf056619
--- /dev/null
+++ b/cbind/mlprec/mld_c_dprec.c
@@ -0,0 +1,21 @@
+#include <stdlib.h>
+#include "mld_c_dprec.h"
+
+mld_c_dprec* mld_c_new_dprec()
+{
+  mld_c_dprec* temp;
+  
+  temp=(mld_c_dprec *) malloc(sizeof(mld_c_dprec));
+  temp->dprec=NULL;
+  return(temp);
+}
+
+
+int mld_c_delete_dprec(mld_c_dprec* p)
+{
+  int iret;
+  iret=mld_c_dprecfree(p);
+  if (iret ==0) free(p);
+  return(iret);
+}
+
diff --git a/cbind/mlprec/mld_c_dprec.h b/cbind/mlprec/mld_c_dprec.h
new file mode 100644
index 00000000..9113f333
--- /dev/null
+++ b/cbind/mlprec/mld_c_dprec.h
@@ -0,0 +1,43 @@
+#ifndef MLD_C_DPREC_
+#define MLD_C_DPREC_
+
+#include "mld_const.h"
+#include "psb_base_cbind.h"
+#include "psb_prec_cbind.h"
+#include "psb_krylov_cbind.h"
+
+/* Object handle related routines */
+/* Note:  mld_get_XXX_handle returns:  <= 0  unsuccessful */
+/*                                     >0    valid handle */
+#ifdef __cplusplus
+extern "C" {
+#endif
+  typedef struct MLD_C_DPREC {
+    void *dprec;
+  } mld_c_dprec; 
+  
+  mld_c_dprec* mld_c_new_dprec();
+  psb_i_t mld_c_delete_dprec(mld_c_dprec* p);
+ 
+  psb_i_t mld_c_dprecinit(mld_c_dprec *ph, const char *ptype);
+  psb_i_t mld_c_dprecseti(mld_c_dprec *ph, const char *what, psb_i_t val);
+  psb_i_t mld_c_dprecsetc(mld_c_dprec *ph, const char *what, const char *val);
+  psb_i_t mld_c_dprecsetr(mld_c_dprec *ph, const char *what, double val);
+  psb_i_t mld_c_dprecbld(psb_c_dspmat *ah, psb_c_descriptor *cdh, mld_c_dprec *ph);
+  psb_i_t mld_c_dhierarchy_build(psb_c_dspmat *ah, psb_c_descriptor *cdh, mld_c_dprec *ph);
+  psb_i_t mld_c_dsmoothers_build(psb_c_dspmat *ah, psb_c_descriptor *cdh, mld_c_dprec *ph);
+  psb_i_t mld_c_dprecfree(mld_c_dprec *ph);
+  psb_i_t mld_c_dprecbld_opt(psb_c_dspmat *ah, psb_c_descriptor *cdh, 
+			  mld_c_dprec *ph, const char *afmt);
+
+
+  psb_i_t mld_c_dkrylov(const char *method, psb_c_dspmat *ah, mld_c_dprec *ph, 
+		  psb_c_dvector *bh, psb_c_dvector *xh,
+		  psb_c_descriptor *cdh, psb_c_SolverOptions *opt);
+
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif
diff --git a/cbind/mlprec/mld_cbind.h b/cbind/mlprec/mld_cbind.h
new file mode 100644
index 00000000..314c4c3a
--- /dev/null
+++ b/cbind/mlprec/mld_cbind.h
@@ -0,0 +1,10 @@
+#ifndef MLD_CBIND_
+#define MLD_CBIND_
+
+#define MLD_VALID_PRECONDITIONER_STRINGS "NONE","DIAG","BJAC","ML","AS"
+#define MLD_VALID_PRECONDITIONER_STRING  "NONE DIAG BJAC ML AS"
+
+#include "mld_const.h"
+#include "mld_c_dprec.h"
+
+#endif
diff --git a/cbind/mlprec/mld_const.h b/cbind/mlprec/mld_const.h
new file mode 100644
index 00000000..39df6045
--- /dev/null
+++ b/cbind/mlprec/mld_const.h
@@ -0,0 +1,14 @@
+/* This file was generated by a script using the mld_base_prec_type.F90 file as a basis. */
+#ifndef MLD_CONST_H_
+#define MLD_CONST_H_
+#ifdef __cplusplus
+extern "C" { 
+#endif
+#define MLD_VERSION_STRING_ ( "2.3.0" )
+#define MLD_VERSION_MAJOR_ ( 2 )
+#define MLD_VERSION_MINOR_ ( 3 )
+#define MLD_PATCHLEVEL_ ( 0 )
+#ifdef __cplusplus
+}
+#endif
+#endif
diff --git a/cbind/mlprec/mld_dprec_cbind_mod.F90 b/cbind/mlprec/mld_dprec_cbind_mod.F90
new file mode 100644
index 00000000..062597fa
--- /dev/null
+++ b/cbind/mlprec/mld_dprec_cbind_mod.F90
@@ -0,0 +1,381 @@
+module mld_dprec_cbind_mod
+
+  use iso_c_binding
+  use mld_prec_mod
+  use psb_base_cbind_mod
+
+  type, bind(c) :: mld_c_dprec
+    type(c_ptr) :: item = c_null_ptr
+  end type mld_c_dprec
+
+contains 
+
+#if 1
+#define MLDC_DEBUG(MSG) write(*,*) __FILE__,':',__LINE__,':',MSG
+#define MLDC_ERROR(MSG) write(*,*) __FILE__,':',__LINE__,':'," ERROR: ",MSG
+#else
+#define MLDC_DEBUG(MSG)
+#define MLDC_ERROR(MSG)
+#endif
+#define mld_success_ 0
+!#define MLDC_ERR_FILTER(INFO) min(0,INFO)
+#define MLDC_ERR_FILTER(INFO) (INFO)
+#define MLDC_ERR_HANDLE(INFO) if(INFO/=mld_success_)MLDC_ERROR("ERROR!")
+    
+  function  mld_c_dprecinit(ictxt,ph,ptype) bind(c) result(res)
+    use psb_base_mod
+    use mld_prec_mod
+    implicit none 
+    
+    integer(psb_c_ipk_)  :: res   
+    type(mld_c_dprec)    :: ph
+    integer(psb_c_ipk_), value :: ictxt
+    character(c_char)     :: ptype(*)
+    integer               :: info
+    type(mld_dprec_type), pointer :: precp
+    character(len=80)     :: fptype
+
+    res = -1
+    res = -1
+    if (c_associated(ph%item)) then 
+      return 
+    end if
+
+    allocate(precp,stat=info)
+    if (info /= 0) return
+    
+    ph%item = c_loc(precp)
+
+    call stringc2f(ptype,fptype)
+    
+    call precp%init(ictxt,fptype,info) 
+     
+    res = MLDC_ERR_FILTER(info)
+    MLDC_ERR_HANDLE(res)
+    return
+  end function mld_c_dprecinit
+
+  function  mld_c_dprecseti(ph,what,val) bind(c) result(res)
+    use psb_base_mod
+    use mld_prec_mod
+    implicit none 
+    
+    integer(psb_c_ipk_) :: res   
+    type(psb_c_object_type) :: ph
+    character(c_char)     :: what(*)
+    integer(psb_c_ipk_), value :: val
+    integer               :: info
+    character(len=80)     :: fwhat
+    type(mld_dprec_type), pointer  :: precp
+
+    res = -1
+    if (c_associated(ph%item)) then 
+      call c_f_pointer(ph%item,precp)
+    else
+      return 
+    end if
+
+    call stringc2f(what,fwhat)
+
+    call mld_precset(precp,fwhat,val,info) 
+     
+    res = MLDC_ERR_FILTER(info)
+    MLDC_ERR_HANDLE(res)
+    return
+  end function mld_c_dprecseti
+    
+
+  function  mld_c_dprecsetr(ph,what,val) bind(c) result(res)
+    use psb_base_mod
+    use mld_prec_mod
+    implicit none 
+    
+    integer(psb_c_ipk_) :: res   
+    type(psb_c_object_type) :: ph
+    character(c_char)     :: what(*)
+    real(c_double), value :: val
+    integer               :: info
+    character(len=80)     :: fwhat
+    type(mld_dprec_type), pointer  :: precp
+
+    res = -1
+    if (c_associated(ph%item)) then 
+      call c_f_pointer(ph%item,precp)
+    else
+      return 
+    end if
+
+    call stringc2f(what,fwhat)
+
+    call mld_precset(precp,fwhat,val,info) 
+     
+    res = MLDC_ERR_FILTER(info)
+    MLDC_ERR_HANDLE(res)
+    return
+  end function mld_c_dprecsetr
+  
+  function  mld_c_dprecsetc(ph,what,val) bind(c) result(res)
+    use psb_base_mod
+    use mld_prec_mod
+    implicit none 
+    
+    integer(psb_c_ipk_) :: res   
+    type(psb_c_object_type) :: ph
+    character(c_char)     :: what(*), val(*)
+    integer               :: info
+    character(len=80)     :: fwhat,fval
+    type(mld_dprec_type), pointer  :: precp
+
+    res = -1
+    if (c_associated(ph%item)) then 
+      call c_f_pointer(ph%item,precp)
+    else
+      return 
+    end if
+
+    call stringc2f(what,fwhat)
+    call stringc2f(val,fval)
+
+    call mld_precset(precp,fwhat,fval,info) 
+     
+    res = MLDC_ERR_FILTER(info)
+    MLDC_ERR_HANDLE(res)
+    return
+  end function mld_c_dprecsetc
+
+  function  mld_c_dprecbld(ah,cdh,ph) bind(c) result(res)
+    use psb_base_mod
+    use mld_prec_mod
+    implicit none 
+    
+    integer(psb_c_ipk_) :: res   
+    type(psb_c_object_type)  :: ph,ah,cdh
+    integer               :: info
+    type(mld_dprec_type), pointer  :: precp
+    type(psb_dspmat_type), pointer :: ap
+    type(psb_desc_type), pointer   :: descp
+    character(len=80)     :: fptype
+
+    res = -1
+
+    if (c_associated(cdh%item)) then 
+      call c_f_pointer(cdh%item,descp)
+    else
+      return 
+    end if
+    if (c_associated(ah%item)) then 
+      call c_f_pointer(ah%item,ap)
+    else
+      return 
+    end if
+    if (c_associated(ph%item)) then 
+      call c_f_pointer(ph%item,precp)
+    else
+      return 
+    end if
+
+    call mld_precbld(ap,descp,precp,info) 
+
+    res = MLDC_ERR_FILTER(info)
+    MLDC_ERR_HANDLE(res)
+
+    return
+  end function mld_c_dprecbld
+    
+  function  mld_c_dhierarchy_build(ah,cdh,ph) bind(c) result(res)
+    use psb_base_mod
+    use mld_prec_mod
+    implicit none 
+    
+    integer(psb_c_ipk_) :: res   
+    type(psb_c_object_type)  :: ph,ah,cdh
+    integer               :: info
+    type(mld_dprec_type), pointer  :: precp
+    type(psb_dspmat_type), pointer :: ap
+    type(psb_desc_type), pointer   :: descp
+    character(len=80)     :: fptype
+
+    res = -1
+
+    if (c_associated(cdh%item)) then 
+      call c_f_pointer(cdh%item,descp)
+    else
+      return 
+    end if
+    if (c_associated(ah%item)) then 
+      call c_f_pointer(ah%item,ap)
+    else
+      return 
+    end if
+    if (c_associated(ph%item)) then 
+      call c_f_pointer(ph%item,precp)
+    else
+      return 
+    end if
+
+    call precp%hierarchy_build(ap,descp,info) 
+
+    res = MLDC_ERR_FILTER(info)
+    MLDC_ERR_HANDLE(res)
+
+    return
+  end function mld_c_dhierarchy_build
+
+  function  mld_c_dsmoothers_build(ah,cdh,ph) bind(c) result(res)
+    use psb_base_mod
+    use mld_prec_mod
+    implicit none 
+    
+    integer(psb_c_ipk_) :: res   
+    type(psb_c_object_type)  :: ph,ah,cdh
+    integer               :: info
+    type(mld_dprec_type), pointer  :: precp
+    type(psb_dspmat_type), pointer :: ap
+    type(psb_desc_type), pointer   :: descp
+    character(len=80)     :: fptype
+
+    res = -1
+
+    if (c_associated(cdh%item)) then 
+      call c_f_pointer(cdh%item,descp)
+    else
+      return 
+    end if
+    if (c_associated(ah%item)) then 
+      call c_f_pointer(ah%item,ap)
+    else
+      return 
+    end if
+    if (c_associated(ph%item)) then 
+      call c_f_pointer(ph%item,precp)
+    else
+      return 
+    end if
+
+    call precp%smoothers_build(ap,descp,info) 
+
+    res = MLDC_ERR_FILTER(info)
+    MLDC_ERR_HANDLE(res)
+
+    return
+  end function mld_c_dsmoothers_build
+    
+  function  mld_c_dkrylov(methd,&
+       & ah,ph,bh,xh,cdh,options) bind(c) result(res)
+    use psb_base_mod
+    use psb_prec_mod
+    use psb_krylov_mod
+    use psb_prec_cbind_mod
+    use psb_dkrylov_cbind_mod
+    implicit none 
+    integer(psb_c_ipk_)          :: res
+    type(psb_c_object_type) :: ah,cdh,ph,bh,xh
+    character(c_char)       :: methd(*)
+    type(solveroptions)     :: options
+    
+    res= mld_c_dkrylov_opt(methd, ah, ph, bh, xh, options%eps,cdh,  &
+         & itmax=options%itmax, iter=options%iter,&
+         & itrace=options%itrace, istop=options%istop,&
+         & irst=options%irst, err=options%err)
+    
+  end function mld_c_dkrylov
+
+
+  function  mld_c_dkrylov_opt(methd,&
+       & ah,ph,bh,xh,eps,cdh,itmax,iter,err,itrace,irst,istop) bind(c) result(res)
+    use psb_base_mod
+    use psb_prec_mod
+    use psb_krylov_mod
+    use psb_objhandle_mod
+    use psb_prec_cbind_mod
+    use psb_base_string_cbind_mod
+    implicit none 
+    integer(psb_c_ipk_)          :: res
+    type(psb_c_object_type) :: ah,cdh,ph,bh,xh
+    integer(psb_c_ipk_), value :: itmax,itrace,irst,istop
+    real(c_double), value :: eps
+    integer(psb_c_ipk_)        :: iter
+    real(c_double)        :: err
+    character(c_char)       :: methd(*)
+    type(psb_desc_type), pointer   :: descp
+    type(psb_dspmat_type), pointer :: ap
+    type(mld_dprec_type), pointer  :: precp
+    type(psb_d_vect_type), pointer :: xp, bp
+
+    integer               :: info,fitmax,fitrace,first,fistop,fiter
+    character(len=20)     :: fmethd
+    real(kind(1.d0))      :: feps,ferr
+
+    res = -1
+    if (c_associated(cdh%item)) then 
+      call c_f_pointer(cdh%item,descp)
+    else
+      return 
+    end if
+    if (c_associated(xh%item)) then 
+      call c_f_pointer(xh%item,xp)
+    else
+      return 
+    end if
+    if (c_associated(bh%item)) then 
+      call c_f_pointer(bh%item,bp)
+    else
+      return 
+    end if
+    if (c_associated(ah%item)) then 
+      call c_f_pointer(ah%item,ap)
+    else
+      return 
+    end if
+    if (c_associated(ph%item)) then 
+      call c_f_pointer(ph%item,precp)
+    else
+      return 
+    end if
+
+    
+    call stringc2f(methd,fmethd)
+    feps    = eps
+    fitmax  = itmax
+    fitrace = itrace
+    first   = irst
+    fistop  = istop
+
+    call psb_krylov(fmethd, ap, precp, bp, xp, feps, &
+         & descp, info,&
+         & itmax=fitmax,iter=fiter,itrace=fitrace,istop=fistop,&
+         & irst=first, err=ferr)
+    iter = fiter
+    err  = ferr
+    res = min(info,0)
+    
+  end function mld_c_dkrylov_opt
+
+  function  mld_c_dprecfree(ph) bind(c) result(res)
+    use psb_base_mod
+    use mld_prec_mod
+    implicit none 
+    
+    integer(psb_c_ipk_) :: res   
+    type(psb_c_object_type) :: ph
+    integer               :: info
+    type(mld_dprec_type), pointer :: precp
+    character(len=80)     :: fptype
+
+    res = -1
+    if (c_associated(ph%item)) then 
+      call c_f_pointer(ph%item,precp)
+    else
+      return 
+    end if
+
+        
+    call precp%free(info) 
+     
+    res = MLDC_ERR_FILTER(info)
+    MLDC_ERR_HANDLE(res)
+    return
+  end function mld_c_dprecfree
+
+end module mld_dprec_cbind_mod
+
diff --git a/cbind/mlprec/mld_prec_cbind_mod.F90 b/cbind/mlprec/mld_prec_cbind_mod.F90
new file mode 100644
index 00000000..e773beba
--- /dev/null
+++ b/cbind/mlprec/mld_prec_cbind_mod.F90
@@ -0,0 +1,8 @@
+module mld_prec_cbind_mod
+
+  use iso_c_binding
+  use psb_base_cbind_mod
+  use psb_prec_cbind_mod
+  use mld_dprec_cbind_mod
+
+end module mld_prec_cbind_mod
diff --git a/cbind/test/pargen/Makefile b/cbind/test/pargen/Makefile
new file mode 100644
index 00000000..36db400d
--- /dev/null
+++ b/cbind/test/pargen/Makefile
@@ -0,0 +1,53 @@
+include ../../Make.inc
+#
+# Libraries used
+#
+LIBDIR=../../lib/
+INCLUDEDIR=../../include/
+HERE=.
+
+FINCLUDES=$(FMFLAG). $(FMFLAG)$(LIBDIR) $(FMFLAG)$(PSBLAS_INCDIR)
+#PSBLAS_LIBS= -L$(PSBLAS_LIBDIR) -L$(LIBDIR)  $(CPSBLAS_LIB) $(PSBLAS_LIB)
+# -lpsb_krylov_cbind -lpsb_prec_cbind -lpsb_base_cbind 
+PSBC_LIBS= -L$(PSBLAS_CBIND_LIBDIR) \
+    -lpsb_krylov_cbind -lpsb_prec_cbind -lpsb_base_cbind 
+MLDC_LIBS=-L$(LIBDIR) -lmld_prec_cbind -lmld_prec -L$(MLDLIBDIR)
+#
+# Compilers and such
+#
+CCOPT= -g
+FINCLUDES=$(FMFLAG)$(LIBDIR) $(FMFLAG).
+CINCLUDES=-I$(LIBDIR) -I$(PSBLAS_CBIND_INCDIR) -I$(INCLUDEDIR)
+
+EXEDIR=./runs
+# 20110404 specifying UMFLIBS here is not portable!
+#UMFLIBS=-lumfpack -lamd -lcholmod -lcolamd -lcamd -lccolamd -L/usr/include/suitesparse
+#UMFFLAGS=-DHave_UMF_ -I/usr/include/suitesparse
+
+all:  mldec
+
+mldec: mldec.o
+	$(MPF90) mldec.o  -o mldec $(MLDC_LIBS)  $(PSBC_LIBS) $(PSBCLDLIBS) $(PSBLAS_LIBS) \
+	  $(UMFLIBS) $(PSBLDLIBS) -lm -lgfortran
+# \
+#	-lifcore -lifcoremt -lguide -limf -lirc -lintlc -lcxaguard -L/opt/intel/fc/10.0.023/lib/ -lm
+
+	/bin/mv mldec $(EXEDIR)
+
+.f90.o:
+	$(MPF90) $(F90COPT) $(FINCLUDES) $(FDEFINES) -c $<
+.c.o:
+	$(MPCC) $(CCOPT) $(CINCLUDES) $(CDEFINES) -c $<
+
+
+clean: 
+	/bin/rm -f mldec.o  $(EXEDIR)/mldec
+verycleanlib: 
+	(cd ../..; make veryclean)
+lib:
+	(cd ../../; make library)
+
+tests: all
+	cd runs ; ./mldec < mlde.inp
+
+
diff --git a/cbind/test/pargen/mldec.c b/cbind/test/pargen/mldec.c
new file mode 100644
index 00000000..61898540
--- /dev/null
+++ b/cbind/test/pargen/mldec.c
@@ -0,0 +1,421 @@
+/*----------------------------------------------------------------------------------*/
+/*              Parallel Sparse BLAS  v2.2					    */
+/*    (C) Copyright 2007 Salvatore Filippone    University of Rome Tor Vergata      */
+/* 										    */
+/*  Redistribution and use in source and binary forms, with or without		    */
+/*  modification, are permitted provided that the following conditions		    */
+/*  are met:									    */
+/*    1. Redistributions of source code must retain the above copyright		    */
+/*       notice, this list of conditions and the following disclaimer.		    */
+/*    2. Redistributions in binary form must reproduce the above copyright	    */
+/*       notice, this list of conditions, and the following disclaimer in the	    */
+/*       documentation and/or other materials provided with the distribution.	    */
+/*    3. The name of the PSBLAS group or the names of its contributors may	    */
+/*       not be used to endorse or promote products derived from this		    */
+/*       software without specific written permission.				    */
+/* 										    */
+/*  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS		    */
+/*  ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED	    */
+/*  TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR	    */
+/*  PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS  */
+/*  BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR	    */
+/*  CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF	    */
+/*  SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS	    */
+/*  INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN	    */
+/*  CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)	    */
+/*  ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE	    */
+/*  POSSIBILITY OF SUCH DAMAGE.							    */
+/* 										    */
+/*  										    */
+/* File: ppdec.c   								    */
+/*										    */
+/* Program: ppdec 								    */
+/* This sample program shows how to build and solve a sparse linear		    */
+/*										    */
+/* The program  solves a linear system based on the partial differential	    */
+/* equation 									    */
+/*										    */
+/* 										    */
+/*										    */
+/* The equation generated is							    */
+/*										    */
+/*   b1 d d (u)  b2  d d (u)    a1 d (u))  a2 d (u)))   			    */
+/* -   ------   -    ------  +    ----- +  ------     + a3 u = 0		    */
+/*      dx dx         dy dy         dx        dy        			    */
+/*										    */
+/* 										    */
+/* with  Dirichlet boundary conditions on the unit cube 			    */
+/*										    */
+/*    0<=x,y,z<=1								    */
+/* 										    */
+/* The equation is discretized with finite differences and uniform stepsize;	    */
+/* the resulting  discrete  equation is						    */
+/*										    */
+/* ( u(x,y,z)(2b1+2b2+a1+a2)+u(x-1,y)(-b1-a1)+u(x,y-1)(-b2-a2)+			    */
+/*  -u(x+1,y)b1-u(x,y+1)b2)*(1/h**2)						    */
+/*										    */
+/* Example adapted from: C.T.Kelley						    */
+/*    Iterative Methods for Linear and Nonlinear Equations			    */
+/*    SIAM 1995									    */
+/*										    */
+/*										    */
+/* In this sample program the index space of the discretized			    */
+/* computational domain is first numbered sequentially in a standard way, 	    */
+/* then the corresponding vector is distributed according to an HPF BLOCK	    */
+/* distribution directive.							    */
+/*										    */
+/* Boundary conditions are set in a very simple way, by adding 			    */
+/* equations of the form							    */
+/*										    */
+/*   u(x,y) = rhs(x,y)								    */
+/*                                                                                  */
+/*----------------------------------------------------------------------------------*/
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <math.h>
+
+#include "psb_cbind.h"
+#include "mld_cbind.h"
+
+#define LINEBUFSIZE 1024
+#define NBMAX       20
+
+double  a1(double x, double y, double  z)
+{
+  return(1.0);
+}
+double a2(double x, double y, double  z)
+{
+  return(20.0*y);
+}
+double a3(double x, double y, double  z)
+{
+  return(1.0);
+}
+double  a4(double x, double y, double  z)
+{
+  return(1.0);
+}
+double  b1(double x, double y, double  z)
+{
+  return(1.0);
+}
+double b2(double x, double y, double  z)
+{
+  return(1.0);
+}
+double b3(double x, double y, double  z)
+{
+  return(1.0);
+}
+
+int matgen(int ictxt, int ng,int idim,int vg[],psb_c_dspmat *ah,psb_c_descriptor *cdh,
+	   psb_c_dvector *xh, psb_c_dvector *bh, psb_c_dvector *rh)
+{
+  int iam, np;
+  int x, y, z, el,glob_row,i,info,ret;
+  double gx, gy, gz, deltah;
+  double val[10*NBMAX], zt[NBMAX];
+  int irow[10*NBMAX], icol[10*NBMAX];
+
+  info = 0;
+  psb_c_info(ictxt,&iam,&np);
+  deltah = (double) 1.0/(idim-1);
+
+  for (glob_row=1; glob_row<=ng; glob_row++) {
+
+    /* Check if I have to do something about this entry */
+    if (vg[glob_row-1] == iam) {
+      el=0;
+      if ( (glob_row%(idim*idim)) == 0) {
+	x = glob_row/(idim*idim);
+      } else {
+	x = glob_row/(idim*idim)+1;
+      }
+      if (((glob_row-(x-1)*idim*idim)%idim) == 0) {
+	y = (glob_row-(x-1)*idim*idim)/idim;
+      } else {
+	y = (glob_row-(x-1)*idim*idim)/idim+1;
+      }
+      z = glob_row-(x-1)*idim*idim-(y-1)*idim;
+      gx=x*deltah;
+      gy=y*deltah;
+      gz=z*deltah;
+      zt[0] = 0.0;
+      /*  internal point: build discretization */
+      /*  term depending on   (x-1,y,z)        */
+	
+      if (x==1) {
+	val[el] = -b1(gx,gy,gz)-a1(gx,gy,gz);
+	val[el] /= deltah*deltah;
+	zt[0] = exp(-gy*gy-gz*gz)*(-val[el]);
+      } else {
+	val[el]=-b1(gx,gy,gz) -a1(gx,gy,gz);
+	val[el] = val[el]/(deltah*deltah);
+	icol[el]=(x-2)*idim*idim+(y-1)*idim+(z);
+	el=el+1;
+      }
+      /*  term depending on     (x,y-1,z) */ 
+      if (y==1) { 
+	val[el]=-b2(gx,gy,gz)-a2(gx,gy,gz);
+	val[el] = val[el]/(deltah*deltah);
+	zt[0] = exp(-gy*gy-gz*gz)*exp(-gx)*(-val[el]);
+      } else {
+	val[el]=-b2(gx,gy,gz)-a2(gx,gy,gz);
+	val[el] = val[el]/(deltah*deltah);
+	icol[el]=(x-1)*idim*idim+(y-2)*idim+(z);
+	el=el+1;
+      }
+      /* term depending on     (x,y,z-1)*/
+      if (z==1) { 
+	val[el]=-b3(gx,gy,gz)-a3(gx,gy,gz);
+	val[el] = val[el]/(deltah*deltah);
+	zt[0] = exp(-gy*gy-gz*gz)*exp(-gx)*(-val[el]);
+      } else {
+	val[el]=-b3(gx,gy,gz)-a3(gx,gy,gz);
+	val[el] = val[el]/(deltah*deltah);
+	icol[el]=(x-1)*idim*idim+(y-1)*idim+(z-1);
+	el=el+1;
+      }
+      /* term depending on     (x,y,z)*/
+      val[el]=2*b1(gx,gy,gz)+2*b2(gx,gy,gz)+2*b3(gx,gy,gz)      
+	+a1(gx,gy,gz)+a2(gx,gy,gz)+a3(gx,gy,gz);
+      val[el] = val[el]/(deltah*deltah);
+      icol[el]=(x-1)*idim*idim+(y-1)*idim+(z);
+      el=el+1;
+      /*  term depending on     (x,y,z+1) */
+      if (z==idim) { 
+	val[el]=-b1(gx,gy,gz);
+	val[el] = val[el]/(deltah*deltah);
+	zt[0] = exp(-gy*gy-gz*gz)*exp(-gx)*(-val[el]);
+      } else {
+	val[el]=-b1(gx,gy,gz);
+	val[el] = val[el]/(deltah*deltah);
+	icol[el]=(x-1)*idim*idim+(y-1)*idim+(z+1);
+	el=el+1;
+      }
+      /* term depending on     (x,y+1,z) */
+      if (y==idim) { 
+	val[el]=-b2(gx,gy,gz);
+	val[el] = val[el]/(deltah*deltah);
+	zt[0] = exp(-gy*gy-gz*gz)*exp(-gx)*(-val[el]);
+      } else {
+	val[el]=-b2(gx,gy,gz);
+	val[el] = val[el]/(deltah*deltah);
+	icol[el]=(x-1)*idim*idim+(y)*idim+(z);
+	el=el+1;
+      }
+      /*  term depending on     (x+1,y,z) */
+      if (x<idim) { 
+	val[el]=-b3(gx,gy,gz);
+	val[el] = val[el]/(deltah*deltah);
+	icol[el]=(x)*idim*idim+(y-1)*idim+(z);
+	el=el+1;
+      }
+      for (i=0; i<el; i++) irow[i]=glob_row;
+      if ((ret=psb_c_dspins(el,irow,icol,val,ah,cdh))!=0) 
+	fprintf(stderr,"From psb_c_dspins: %d\n",ret); 
+      psb_c_dgeins(1,irow,zt,bh,cdh);
+      zt[0]=0.0;
+      psb_c_dgeins(1,irow,zt,xh,cdh);
+    }
+  }
+  
+  if ((info=psb_c_cdasb(cdh))!=0)  return(info);
+  
+  if ((info=psb_c_dspasb(ah,cdh))!=0)  return(info); 
+  
+  if ((info=psb_c_dgeasb(xh,cdh))!=0)  return(info);
+  if ((info=psb_c_dgeasb(bh,cdh))!=0)  return(info);
+  if ((info=psb_c_dgeasb(rh,cdh))!=0)  return(info);
+  return(info);
+
+}
+  
+int main(int argc, char *argv[])
+{
+  int ictxt, iam, np;
+  char methd[40], ptype[20], afmt[8], buffer[LINEBUFSIZE+1];
+  int nparms;
+  int idim,info,istop,itmax,itrace,irst,i,iter,ret;
+  mld_c_dprec *ph;
+  psb_c_dspmat *ah;
+  psb_c_dvector *bh, *xh, *rh; 
+  int *vg, ng, nb,nlr; 
+  double t1,t2,eps,err;
+  double *xv, *bv, *rv;
+  double one=1.0, zero=0.0, res2;
+  psb_c_SolverOptions options;
+  psb_c_descriptor *cdh;
+  
+
+  ictxt = psb_c_init();
+  psb_c_info(ictxt,&iam,&np);
+  fprintf(stdout,"Initialization: am %d of %d\n",iam,np);
+
+  fflush(stdout);
+  psb_c_barrier(ictxt);
+  if (iam == 0) {
+    fgets(buffer,LINEBUFSIZE,stdin);
+    sscanf(buffer,"%d ",&nparms);
+    fgets(buffer,LINEBUFSIZE,stdin);
+    sscanf(buffer,"%s",methd);
+    fgets(buffer,LINEBUFSIZE,stdin);
+    sscanf(buffer,"%s",ptype);
+    fgets(buffer,LINEBUFSIZE,stdin);
+    sscanf(buffer,"%s",afmt);
+    fgets(buffer,LINEBUFSIZE,stdin);
+    sscanf(buffer,"%d",&idim);
+    fgets(buffer,LINEBUFSIZE,stdin);
+    sscanf(buffer,"%d",&istop);
+    fgets(buffer,LINEBUFSIZE,stdin);
+    sscanf(buffer,"%d",&itmax);
+    fgets(buffer,LINEBUFSIZE,stdin);
+    sscanf(buffer,"%d",&itrace);
+    fgets(buffer,LINEBUFSIZE,stdin);
+    sscanf(buffer,"%d",&irst);
+  } 
+  /* Now broadcast the values, and check they're OK */
+  psb_c_ibcast(ictxt,1,&nparms,0);
+  psb_c_hbcast(ictxt,methd,0);
+  psb_c_hbcast(ictxt,ptype,0);
+  psb_c_hbcast(ictxt,afmt,0);
+  psb_c_ibcast(ictxt,1,&idim,0);
+  psb_c_ibcast(ictxt,1,&istop,0);
+  psb_c_ibcast(ictxt,1,&itmax,0);
+  psb_c_ibcast(ictxt,1,&itrace,0);
+  psb_c_ibcast(ictxt,1,&irst,0);
+  
+  psb_c_barrier(ictxt);
+
+  cdh=psb_c_new_descriptor();
+
+  /* Simple minded BLOCK data distribution */ 
+  ng = idim*idim*idim;
+  nb = (ng+np-1)/np;
+  if ((vg=malloc(ng*sizeof(int)))==NULL) {
+    fprintf(stderr,"On %d: malloc failure\n",iam);
+    psb_c_abort(ictxt);
+  }
+  for (i=0; i<ng; i++) {
+    vg[i] = i/nb; 
+  }
+  if ((info=psb_c_cdall_vg(ng,vg,ictxt,cdh))!=0) {
+    fprintf(stderr,"From cdall: %d\nBailing out\n",info);
+    psb_c_abort(ictxt);
+  }
+
+  bh  = psb_c_new_dvector();
+  xh  = psb_c_new_dvector();
+  rh  = psb_c_new_dvector();
+  ah  = psb_c_new_dspmat();
+
+  /* Allocate mem space for sparse matrix and vectors */ 
+  ret=psb_c_dspall(ah,cdh); 
+  psb_c_dgeall(bh,cdh);
+  psb_c_dgeall(xh,cdh);
+  psb_c_dgeall(rh,cdh);
+  /* Matrix generation  */
+  if (matgen(ictxt,ng,idim,vg,ah,cdh,xh,bh,rh) != 0) {
+    fprintf(stderr,"Error during matrix build loop\n");
+    psb_c_abort(ictxt);
+  }    
+  psb_c_barrier(ictxt);
+  /* Set up the preconditioner */ 
+  ph  = mld_c_new_dprec();
+  mld_c_dprecinit(ph,ptype,3);
+  mld_c_dprecseti(ph,"SMOOTHER_SWEEPS",2);
+  mld_c_dprecseti(ph,"SUB_FILLIN",1);
+  mld_c_dprecsetc(ph,"COARSE_SOLVE","BJAC");
+  mld_c_dprecsetc(ph,"COARSE_SUBSOLVE","ILU");
+  mld_c_dprecseti(ph,"COARSE_FILLIN",1);
+  mld_c_dprecbld(ah,cdh,ph); 
+
+  psb_c_barrier(ictxt);
+  /* Set up the solver options */ 
+  psb_c_DefaultSolverOptions(&options);
+  options.eps    = 1.e-9;
+  options.itmax  = itmax;
+  options.irst   = irst;
+  options.itrace = 1;
+  options.istop  = istop;
+
+  t1=psb_c_wtime();
+  ret=mld_c_dkrylov(methd,ah,ph,bh,xh,cdh,&options); 
+  t2=psb_c_wtime();
+  iter = options.iter;
+  err  = options.err;
+  fprintf(stderr,"From krylov: %d %lf, %d\n",iter,err,ret);
+  
+
+  /* Check 2-norm of residual on exit */ 
+  psb_c_dgeaxpby(one,bh,zero,rh,cdh); 
+  psb_c_dspmm(-one,ah,xh,one,rh,cdh); 
+  res2=psb_c_dgenrm2(rh,cdh);
+
+  if (iam==0) {
+    fprintf(stdout,"Time: %lf\n",(t2-t1));
+    fprintf(stdout,"Iter: %d\n",iter);
+    fprintf(stdout,"Err: %lg\n",err);
+    fprintf(stdout,"||r||_2: %lg\n",res2);
+    
+  }
+
+#if 0
+  bv = psb_c_dvect_get_cpy(bh);
+  nlr=psb_c_cd_get_local_rows(cdh);
+  for (i=0;i<nlr; i++)
+    fprintf(stdout,"RHS: %d %d %lf\n",iam,i,bv[i]);
+
+
+  xv = psb_c_dvect_get_cpy(xh);
+  nlr=psb_c_cd_get_local_rows(cdh);
+  for (i=0;i<nlr; i++)
+    fprintf(stdout,"SOL: %d %d %lf\n",iam,i,xv[i]);
+
+
+  rv = psb_c_dvect_get_cpy(rh);
+  nlr=psb_c_cd_get_local_rows(cdh);
+  for (i=0;i<nlr; i++)
+    fprintf(stdout,"RES: %d %d %lf\n",iam,i,rv[i]);
+
+#endif
+  
+  /* Clean up memory */ 
+  if ((info=psb_c_dgefree(xh,cdh))!=0) {
+    fprintf(stderr,"From dgefree: %d\nBailing out\n",info);
+    psb_c_abort(ictxt);
+  }
+  if ((info=psb_c_dgefree(bh,cdh))!=0) {
+    fprintf(stderr,"From dgefree: %d\nBailing out\n",info);
+    psb_c_abort(ictxt);
+  }
+  if ((info=psb_c_dgefree(rh,cdh))!=0) {
+    fprintf(stderr,"From dgefree: %d\nBailing out\n",info);
+    psb_c_abort(ictxt);
+  }
+
+  if ((info=psb_c_cdfree(cdh))!=0) {
+    fprintf(stderr,"From cdfree: %d\nBailing out\n",info);
+    psb_c_abort(ictxt);
+  }
+  fprintf(stderr,"pointer  from cdfree: %p\n",cdh->descriptor);
+ 
+  /* Clean up object handles  */   
+  free(ph);
+  free(xh);
+  free(bh);
+  free(ah);
+  free(cdh);
+
+
+  fprintf(stderr,"program completed successfully\n");
+
+  psb_c_barrier(ictxt);
+  psb_c_exit(ictxt);
+}
+
+
+
diff --git a/cbind/test/pargen/runs/mlde.inp b/cbind/test/pargen/runs/mlde.inp
new file mode 100644
index 00000000..365c689d
--- /dev/null
+++ b/cbind/test/pargen/runs/mlde.inp
@@ -0,0 +1,11 @@
+7                      Number of entries below this
+BICGSTAB                 Iterative method BICGSTAB CGS  BICG BICGSTABL RGMRES
+ML                      Preconditioner NONE  DIAG  BJAC 
+CSR                    A Storage format CSR COO JAD 
+60                     Domain size (acutal system is this**3)
+1                      Stopping criterion
+80                     MAXIT
+01                     ITRACE
+20                     IRST    restart for RGMRES  and BiCGSTABL
+
+