diff --git a/tests/pdegen/ppde.f90 b/tests/pdegen/ppde.f90
index 8db81503..a99a1002 100644
--- a/tests/pdegen/ppde.f90
+++ b/tests/pdegen/ppde.f90
@@ -432,9 +432,9 @@ contains
     integer                        :: ictxt, info
     character                      :: afmt*5
     type(psb_d_sparse_mat)         :: a
-    real(psb_dpk_)           :: zt(nb),glob_x,glob_y,glob_z
+    real(psb_dpk_)           :: zt(nb),x,y,z
     integer                  :: m,n,nnz,glob_row,nlr,i,ii,ib,k
-    integer                  :: x,y,z,ia,indx_owner
+    integer                  :: ix,iy,iz,ia,indx_owner
     integer                  :: np, iam, nr, nt
     integer                  :: element
     integer, allocatable     :: irow(:),icol(:),myidx(:)
@@ -526,20 +526,20 @@ contains
         glob_row=myidx(i)
         ! compute gridpoint coordinates
         if (mod(glob_row,(idim*idim)) == 0) then
-          x = glob_row/(idim*idim)
+          ix = glob_row/(idim*idim)
         else
-          x = glob_row/(idim*idim)+1
+          ix = glob_row/(idim*idim)+1
         endif
-        if (mod((glob_row-(x-1)*idim*idim),idim) == 0) then
-          y = (glob_row-(x-1)*idim*idim)/idim
+        if (mod((glob_row-(ix-1)*idim*idim),idim) == 0) then
+          iy = (glob_row-(ix-1)*idim*idim)/idim
         else
-          y = (glob_row-(x-1)*idim*idim)/idim+1
+          iy = (glob_row-(ix-1)*idim*idim)/idim+1
         endif
-        z = glob_row-(x-1)*idim*idim-(y-1)*idim
-        ! glob_x, glob_y, glob_x coordinates
-        glob_x=x*deltah
-        glob_y=y*deltah
-        glob_z=z*deltah
+        iz = glob_row-(ix-1)*idim*idim-(iy-1)*idim
+        ! x, y, x coordinates
+        x = ix*deltah
+        y = iy*deltah
+        z = iz*deltah
 
         ! check on boundary points 
         zt(k) = 0.d0
@@ -547,99 +547,81 @@ contains
         !   
         !  term depending on   (x-1,y,z)
         !
-        if (x==1) then 
-          val(element)=-b1(glob_x,glob_y,glob_z)&
-               & -a1(glob_x,glob_y,glob_z)
+        if (ix==1) then 
+          val(element)=-b1(x,y,z)-a1(x,y,z)
           val(element) = val(element)/(deltah*&
                & deltah)
-          zt(k) = exp(-glob_y**2-glob_z**2)*(-val(element))
+          zt(k) = exp(-y**2-z**2)*(-val(element))
         else
-          val(element)=-b1(glob_x,glob_y,glob_z)&
-               & -a1(glob_x,glob_y,glob_z)
+          val(element)=-b1(x,y,z)-a1(x,y,z)
           val(element) = val(element)/(deltah*&
                & deltah)
-          icol(element) = (x-2)*idim*idim+(y-1)*idim+(z)
+          icol(element) = (ix-2)*idim*idim+(iy-1)*idim+(iz)
           irow(element) = glob_row
           element       = element+1
         endif
         !  term depending on     (x,y-1,z)
-        if (y==1) then 
-          val(element)=-b2(glob_x,glob_y,glob_z)&
-               & -a2(glob_x,glob_y,glob_z)
+        if (iy==1) then 
+          val(element)=-b2(x,y,z)-a2(x,y,z)
           val(element) = val(element)/(deltah*&
                & deltah)
-          zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element))  
+          zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))  
         else
-          val(element)=-b2(glob_x,glob_y,glob_z)&
-               & -a2(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          icol(element) = (x-1)*idim*idim+(y-2)*idim+(z)
+          val(element)=-b2(x,y,z)-a2(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          icol(element) = (ix-1)*idim*idim+(iy-2)*idim+(iz)
           irow(element) = glob_row
           element       = element+1
         endif
         !  term depending on     (x,y,z-1)
-        if (z==1) then 
-          val(element)=-b3(glob_x,glob_y,glob_z)&
-               & -a3(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element))  
+        if (iz==1) then 
+          val(element)=-b3(x,y,z)-a3(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))  
         else
-          val(element)=-b3(glob_x,glob_y,glob_z)&
-               & -a3(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          icol(element) = (x-1)*idim*idim+(y-1)*idim+(z-1)
+          val(element)=-b3(x,y,z)-a3(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          icol(element) = (ix-1)*idim*idim+(iy-1)*idim+(iz-1)
           irow(element) = glob_row
           element       = element+1
         endif
         !  term depending on     (x,y,z)
-        val(element)=2*b1(glob_x,glob_y,glob_z)&
-             & +2*b2(glob_x,glob_y,glob_z)&
-             & +2*b3(glob_x,glob_y,glob_z)&
-             & +a1(glob_x,glob_y,glob_z)&
-             & +a2(glob_x,glob_y,glob_z)&
-             & +a3(glob_x,glob_y,glob_z)
-        val(element) = val(element)/(deltah*&
-             & deltah)
-        icol(element) = (x-1)*idim*idim+(y-1)*idim+(z)
+        val(element)=2*b1(x,y,z) + 2*b2(x,y,z)&
+             & + 2*b3(x,y,z) + a1(x,y,z)&
+             & + a2(x,y,z) + a3(x,y,z)
+        val(element) = val(element)/(deltah*deltah)
+        icol(element) = (ix-1)*idim*idim+(iy-1)*idim+(iz)
         irow(element) = glob_row
         element       = element+1                  
         !  term depending on     (x,y,z+1)
-        if (z==idim) then 
-          val(element)=-b1(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element))  
+        if (iz==idim) then 
+          val(element)=-b1(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))  
         else
-          val(element)=-b1(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          icol(element) = (x-1)*idim*idim+(y-1)*idim+(z+1)
+          val(element)=-b1(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          icol(element) = (ix-1)*idim*idim+(iy-1)*idim+(iz+1)
           irow(element) = glob_row
           element       = element+1
         endif
         !  term depending on     (x,y+1,z)
-        if (y==idim) then 
-          val(element)=-b2(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element))  
+        if (iy==idim) then 
+          val(element)=-b2(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))  
         else
-          val(element)=-b2(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          icol(element) = (x-1)*idim*idim+(y)*idim+(z)
+          val(element)=-b2(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          icol(element) = (ix-1)*idim*idim+(iy)*idim+(iz)
           irow(element) = glob_row
           element       = element+1
         endif
         !  term depending on     (x+1,y,z)
-        if (x<idim) then 
-          val(element)=-b3(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          icol(element) = (x)*idim*idim+(y-1)*idim+(z)
+        if (ix<idim) then 
+          val(element)=-b3(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          icol(element) = (ix)*idim*idim+(iy-1)*idim+(iz)
           irow(element) = glob_row
           element       = element+1
         endif
diff --git a/tests/pdegen/spde.f90 b/tests/pdegen/spde.f90
index 6d4a2c68..7b7f2d26 100644
--- a/tests/pdegen/spde.f90
+++ b/tests/pdegen/spde.f90
@@ -419,9 +419,9 @@ contains
     integer                        :: ictxt, info
     character                      :: afmt*5
     type(psb_s_sparse_mat)    :: a
-    real(psb_spk_)           :: zt(nb),glob_x,glob_y,glob_z
+    real(psb_spk_)           :: zt(nb),x,y,z
     integer                  :: m,n,nnz,glob_row,nlr,i,ii,ib,k
-    integer                  :: x,y,z,ia,indx_owner
+    integer                  :: ix,iy,iz,ia,indx_owner
     integer                  :: np, iam, nr, nt
     integer                  :: element
     integer, allocatable     :: irow(:),icol(:),myidx(:)
@@ -513,20 +513,20 @@ contains
         glob_row=myidx(i)
         ! compute gridpoint coordinates
         if (mod(glob_row,(idim*idim)) == 0) then
-          x = glob_row/(idim*idim)
+          ix = glob_row/(idim*idim)
         else
-          x = glob_row/(idim*idim)+1
+          ix = glob_row/(idim*idim)+1
         endif
-        if (mod((glob_row-(x-1)*idim*idim),idim) == 0) then
-          y = (glob_row-(x-1)*idim*idim)/idim
+        if (mod((glob_row-(ix-1)*idim*idim),idim) == 0) then
+          iy = (glob_row-(ix-1)*idim*idim)/idim
         else
-          y = (glob_row-(x-1)*idim*idim)/idim+1
+          iy = (glob_row-(ix-1)*idim*idim)/idim+1
         endif
-        z = glob_row-(x-1)*idim*idim-(y-1)*idim
-        ! glob_x, glob_y, glob_x coordinates
-        glob_x=x*deltah
-        glob_y=y*deltah
-        glob_z=z*deltah
+        iz = glob_row-(ix-1)*idim*idim-(iy-1)*idim
+        ! x, y, x coordinates
+        x = ix*deltah
+        y = iy*deltah
+        z = iz*deltah
 
         ! check on boundary points 
         zt(k) = 0.d0
@@ -534,99 +534,81 @@ contains
         !   
         !  term depending on   (x-1,y,z)
         !
-        if (x==1) then 
-          val(element)=-b1(glob_x,glob_y,glob_z)&
-               & -a1(glob_x,glob_y,glob_z)
+        if (ix==1) then 
+          val(element)=-b1(x,y,z)-a1(x,y,z)
           val(element) = val(element)/(deltah*&
                & deltah)
-          zt(k) = exp(-glob_y**2-glob_z**2)*(-val(element))
+          zt(k) = exp(-y**2-z**2)*(-val(element))
         else
-          val(element)=-b1(glob_x,glob_y,glob_z)&
-               & -a1(glob_x,glob_y,glob_z)
+          val(element)=-b1(x,y,z)-a1(x,y,z)
           val(element) = val(element)/(deltah*&
                & deltah)
-          icol(element) = (x-2)*idim*idim+(y-1)*idim+(z)
+          icol(element) = (ix-2)*idim*idim+(iy-1)*idim+(iz)
           irow(element) = glob_row
           element       = element+1
         endif
         !  term depending on     (x,y-1,z)
-        if (y==1) then 
-          val(element)=-b2(glob_x,glob_y,glob_z)&
-               & -a2(glob_x,glob_y,glob_z)
+        if (iy==1) then 
+          val(element)=-b2(x,y,z)-a2(x,y,z)
           val(element) = val(element)/(deltah*&
                & deltah)
-          zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element))  
+          zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))  
         else
-          val(element)=-b2(glob_x,glob_y,glob_z)&
-               & -a2(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          icol(element) = (x-1)*idim*idim+(y-2)*idim+(z)
+          val(element)=-b2(x,y,z)-a2(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          icol(element) = (ix-1)*idim*idim+(iy-2)*idim+(iz)
           irow(element) = glob_row
           element       = element+1
         endif
         !  term depending on     (x,y,z-1)
-        if (z==1) then 
-          val(element)=-b3(glob_x,glob_y,glob_z)&
-               & -a3(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element))  
+        if (iz==1) then 
+          val(element)=-b3(x,y,z)-a3(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))  
         else
-          val(element)=-b3(glob_x,glob_y,glob_z)&
-               & -a3(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          icol(element) = (x-1)*idim*idim+(y-1)*idim+(z-1)
+          val(element)=-b3(x,y,z)-a3(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          icol(element) = (ix-1)*idim*idim+(iy-1)*idim+(iz-1)
           irow(element) = glob_row
           element       = element+1
         endif
         !  term depending on     (x,y,z)
-        val(element)=2*b1(glob_x,glob_y,glob_z)&
-             & +2*b2(glob_x,glob_y,glob_z)&
-             & +2*b3(glob_x,glob_y,glob_z)&
-             & +a1(glob_x,glob_y,glob_z)&
-             & +a2(glob_x,glob_y,glob_z)&
-             & +a3(glob_x,glob_y,glob_z)
-        val(element) = val(element)/(deltah*&
-             & deltah)
-        icol(element) = (x-1)*idim*idim+(y-1)*idim+(z)
+        val(element)=2*b1(x,y,z) + 2*b2(x,y,z)&
+             & + 2*b3(x,y,z) + a1(x,y,z)&
+             & + a2(x,y,z) + a3(x,y,z)
+        val(element) = val(element)/(deltah*deltah)
+        icol(element) = (ix-1)*idim*idim+(iy-1)*idim+(iz)
         irow(element) = glob_row
         element       = element+1                  
         !  term depending on     (x,y,z+1)
-        if (z==idim) then 
-          val(element)=-b1(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element))  
+        if (iz==idim) then 
+          val(element)=-b1(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))  
         else
-          val(element)=-b1(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          icol(element) = (x-1)*idim*idim+(y-1)*idim+(z+1)
+          val(element)=-b1(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          icol(element) = (ix-1)*idim*idim+(iy-1)*idim+(iz+1)
           irow(element) = glob_row
           element       = element+1
         endif
         !  term depending on     (x,y+1,z)
-        if (y==idim) then 
-          val(element)=-b2(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          zt(k) = exp(-glob_y**2-glob_z**2)*exp(-glob_x)*(-val(element))  
+        if (iy==idim) then 
+          val(element)=-b2(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          zt(k) = exp(-y**2-z**2)*exp(-x)*(-val(element))  
         else
-          val(element)=-b2(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          icol(element) = (x-1)*idim*idim+(y)*idim+(z)
+          val(element)=-b2(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          icol(element) = (ix-1)*idim*idim+(iy)*idim+(iz)
           irow(element) = glob_row
           element       = element+1
         endif
         !  term depending on     (x+1,y,z)
-        if (x<idim) then 
-          val(element)=-b3(glob_x,glob_y,glob_z)
-          val(element) = val(element)/(deltah*&
-               & deltah)
-          icol(element) = (x)*idim*idim+(y-1)*idim+(z)
+        if (ix<idim) then 
+          val(element)=-b3(x,y,z)
+          val(element) = val(element)/(deltah*deltah)
+          icol(element) = (ix)*idim*idim+(iy-1)*idim+(iz)
           irow(element) = glob_row
           element       = element+1
         endif