diff --git a/vlaplex.f90 b/vlaplex.f90
index 6bc3423..bef0e04 100644
--- a/vlaplex.f90
+++ b/vlaplex.f90
@@ -95,7 +95,7 @@ program VlaPlEx
 
   ! Set input parameters (remember these have to be in non-dimensional units)
   c_s      = sqrt(11.0_dp * gamma_i * 1.0_dp)
-  bc_file = 'bc_fa_T30Z11.csv'
+  bc_file = 'bc.csv'
   call boundaryConditions%init(bc_file)
 
   ! Set domain boundaries (non-dimensional units)
@@ -104,7 +104,7 @@ program VlaPlEx
   dr =   1.0e-6_dp / L_ref
   nr = nint((rf - r0) / dr) + 1
   dr =      (rf - r0) / float(nr-1)
-  print *, '#R gridpoints: ', nr 
+
   allocate(r(1:nr))
   do i = 1, nr
     r(i) = dr * float(i-1) + r0
@@ -121,7 +121,7 @@ program VlaPlEx
   dv = 1.0e-1_dp
   nv = nint((vf - v0) / dv) + 1
   dv =      (vf - v0) / float(nv-1)
-  print *, '#V gridpoints: ', nv
+
   allocate(v(1:nv))
   do j = 1, nv
     v(j) = dv * float(j-1) + v0
@@ -139,7 +139,6 @@ program VlaPlEx
   dt = 1.0e-2_dp*dr/c_s
   nt = nint((tf - t0) / dt) 
   dt =      (tf - t0) / float(nt)
-  print *, '#timesteps: ', nt 
 
   everyOutput = nint(1.0e-9_dp/t_ref/dt)
   if (everyOutput  == 0) then
@@ -237,9 +236,6 @@ program VlaPlEx
   do t = 1, nt
     time = t * dt + t0
     call boundaryConditions%get(time, n_bc, u_bc, Temp_bc, Zave_bc)
-    print *, 'Time: ', time * t_ref
-    print *, 'Temp_bc: ', Temp_bc
-    print *, 'Zave_bc: ', Zave_bc
     z_inj = minloc(abs(Z_list - T_to_Z(Temp_bc)),1)
     Zave_bc = Z_list(z_inj)
     u_bc = sqrt(Zave_bc * Temp_bc)
@@ -254,9 +250,6 @@ program VlaPlEx
     f_i(:,1,j0:nv) = f_i_old(:,1,j0:nv)
 
     T_e = Temp_bc
-    print *, 'TtoZ: ', T_to_Z(Temp_bc)
-    print *, 'Injection species: ', Zave_bc
-    print *, '-------------------------'
 
     ! r = rf, v<0
     f_i_old(:,nr,1:j0-1) = 0.0_dp