diff --git a/vlaplex.f90 b/vlaplex.f90
index 633924f..a804b59 100644
--- a/vlaplex.f90
+++ b/vlaplex.f90
@@ -67,6 +67,7 @@ program VlaPlEx
   real(dp), allocatable, dimension(:,:):: T_i
   real(dp), allocatable, dimension(:):: n_e
   real(dp), allocatable, dimension(:):: Zave
+  real(dp), allocatable, dimension(:):: Z_list
   real(dp), allocatable, dimension(:):: diag, diag_low, diag_high
   real(dp), allocatable, dimension(:,:):: A
   real(dp), allocatable, dimension(:):: Res
@@ -157,13 +158,14 @@ program VlaPlEx
 
   write(*, '(A,ES0.4e3)') 'CFL: ', dt*vf/dr
 
-  nz = 2
+  nz = 16
   ! Allocate vectors
   allocate(f_i(1:nz,1:nr,1:nv), f_i_old(1:nz,1:nr,1:nv))
   allocate(n_i(1:nz,1:nr))
   allocate(sum_ni(1:nr))
   allocate(u_i(1:nz,1:nr), E_i(1:nr), T_i(1:nz,1:nr))
   allocate(Zave(1:nr))
+  allocate(Z_list(1:nz))
   allocate(n_e(1:nr))
   allocate(phi(1:nr), phi_old(1:nr), E(1:nr))
   allocate(fCum_i(1:nz,1:nv))
@@ -177,6 +179,7 @@ program VlaPlEx
   n_e     = 0.0_dp
   T_e     = 0.0_dp
   Zave    = 0.0_dp
+  Z_list  = (/(I, I = 1, nz) /)
   phi     = 0.0_dp
   phi_old = 0.0_dp
   E       = 0.0_dp
@@ -268,7 +271,6 @@ program VlaPlEx
     f_i_old(:,:,nv) = 0.0_dp
     sum_ni = 0.0_dp
     ! Advect in the r direction
-    
     do iz = 1, nz
       !$omp parallel do
       do i = 1, nr
@@ -300,7 +302,6 @@ program VlaPlEx
       !$omp end parallel do
       sum_ni = sum_ni + Zave * n_i(iz,:)
     end do
-    
 
     ! Assume quasi-neutrality to start iterating
     n_e = 1.0_dp/nz * sum_ni