diff --git a/doc/user-manual/fpakc_UserManual.pdf b/doc/user-manual/fpakc_UserManual.pdf
index 3f86031..8af12d5 100644
Binary files a/doc/user-manual/fpakc_UserManual.pdf and b/doc/user-manual/fpakc_UserManual.pdf differ
diff --git a/doc/user-manual/fpakc_UserManual.tex b/doc/user-manual/fpakc_UserManual.tex
index c2dfd9e..7a17942 100644
--- a/doc/user-manual/fpakc_UserManual.tex
+++ b/doc/user-manual/fpakc_UserManual.tex
@@ -48,9 +48,8 @@
 \newglossaryentry{openblas}{name={OpenBLAS},description={Open-source implementation of BLAS and LAPACK APIs}}
 \newglossaryentry{git}{name={Git},description={Git is a distributed version-control system for tracking changes in a set of files}}
 \newglossaryentry{gitlab}{name={GitLab},description={GitLab is a web-based lifecycle tool that provides a \Gls{git}-repository manager}}
-\newglossaryentry{gmsh}{name={Gmsh},description={A three-dimensional finite element mesh generator with built-in pre- and post-processing facilities}}
-\newglossaryentry{gnuplot}{name={Gnuplot},description={A portable command-line driven graphing utility for Linux, OS/2, MS Windows, OSX, VMS, and many other platforms}}
-\newglossaryentry{paraview}{name={ParaView},description={Tool for post-processing and visualizing data from various formats like vtu}}
+\newglossaryentry{gmsh}{name={Gmsh},description={A three-dimensional finite element mesh generator with built-in pre- and post-processing facilities.}}
+\newglossaryentry{gnuplot}{name={Gnuplot},description={A portable command-line driven graphing utility for Linux, OS/2, MS Windows, OSX, VMS, and many other platforms.}}
 
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 \bibliography{bibliography}
@@ -847,14 +846,14 @@ make
   \chapter{Example runs}\label{ch:exampleRuns}
     This chapter presents a description of the different example files distributed with \acrshort{fpakc}.
     All examples in the repository have a \lstinline|README.txt| file and a reference output.
-    Plotting of the output is done with \Gls{gnuplot} or \Gls{paraview}.
+    Plotting of the output is done with \Gls{gnuplot} or \Gls{gmsh}.
     
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     \section{1D Emissive Cathode (1D\_Cathode)}
       Emission from a 1D cathode in both, Cartesian and radial coordinates.
       Both cases insert the same number of electrons from the minimum coordinate and have the same boundary conditions for particles and the electrostatic field.
       This case is useful to illustrate how \acrshort{fpakc} can deal with different geometries by just modifying some parameters in the input file.
-      The same mesh file (\lstinline|mesh.csv|) is used for both cases.
+      The same mesh file (\lstinline|mesh.msh|) is used for both cases.
       
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     \section{0D \ce{Ar}-\ce{Ar+} Elastic Collision (0D\_Argon)}
diff --git a/src/modules/init/moduleInput.f90 b/src/modules/init/moduleInput.f90
index 9494517..077e9d6 100644
--- a/src/modules/init/moduleInput.f90
+++ b/src/modules/init/moduleInput.f90
@@ -509,7 +509,6 @@ MODULE moduleInput
       USE moduleRefParam
       USE moduleList
       USE json_module
-      use moduleMesh, only: qSpecies
       IMPLICIT NONE
 
       TYPE(json_file), INTENT(inout):: config
@@ -518,7 +517,7 @@ MODULE moduleInput
       CHARACTER(:), ALLOCATABLE:: speciesType
       REAL(8):: mass, charge
       LOGICAL:: found
-      INTEGER:: s
+      INTEGER:: i
       CHARACTER(:), ALLOCATABLE:: linkName
       INTEGER:: linkID
 
@@ -530,8 +529,8 @@ MODULE moduleInput
       ALLOCATE(species(1:nSpecies))
 
       !Reads information of individual species
-      DO s = 1, nSpecies
-        WRITE(iString, '(I2)') s
+      DO i = 1, nSpecies
+        WRITE(iString, '(I2)') i
         object = 'species(' // TRIM(iString) // ')'
         CALL config%get(object // '.type', speciesType, found)
         CALL config%get(object // '.mass', mass,        found)
@@ -540,12 +539,12 @@ MODULE moduleInput
         !Allocate species depending on type and assign specific parameters
         SELECT CASE(speciesType)
         CASE ("neutral")
-          ALLOCATE(species(s)%obj, source=speciesNeutral())
+          ALLOCATE(species(i)%obj, source=speciesNeutral())
 
         CASE ("charged")
           CALL config%get(object // '.charge', charge,      found)
           IF (.NOT. found) CALL criticalError("Required parameter charge not found for species " // object, 'readSpecies')
-          ALLOCATE(species(s)%obj, source=speciesCharged(q  = charge, &
+          ALLOCATE(species(i)%obj, source=speciesCharged(q  = charge, &
                                                          qm = charge/mass))
 
         CASE DEFAULT
@@ -553,32 +552,18 @@ MODULE moduleInput
 
         END SELECT
         !Assign shared parameters for all species
-        CALL config%get(object // '.name',   species(s)%obj%name,   found)
-        CALL config%get(object // '.weight', species(s)%obj%weight, found)
-        species(s)%obj%n = s
-        species(s)%obj%m = mass
-
-      END DO
-
-      ! Allocate the vector with the species charges for calculating the EM field
-      ALLOCATE(qSpecies(1:nSpecies))
-      DO s = 1, nSpecies
-        SELECT TYPE(sp => species(s)%obj)
-        TYPE IS (speciesCharged)
-          qSpecies(s) = sp%q
-
-        CLASS DEFAULT
-          qSpecies(s) = 0.D0
-
-        END SELECT
+        CALL config%get(object // '.name',   species(i)%obj%name,   found)
+        CALL config%get(object // '.weight', species(i)%obj%weight, found)
+        species(i)%obj%n = i
+        species(i)%obj%m = mass
 
       END DO
 
       !Read relations between species
-      DO s = 1, nSpecies
-        WRITE(iString, '(I2)') s
+      DO i = 1, nSpecies
+        WRITE(iString, '(I2)') i
         object = 'species(' // TRIM(iString) // ')'
-        SELECT TYPE(sp => species(s)%obj)
+        SELECT TYPE(sp => species(i)%obj)
         TYPE IS (speciesNeutral)
           !Get species linked ion
           CALL config%get(object // '.ion', linkName, found)
@@ -936,7 +921,7 @@ MODULE moduleInput
       TYPE(json_file), INTENT(inout):: config
       CHARACTER(:), ALLOCATABLE:: object
       LOGICAL:: found
-      INTEGER:: b
+      INTEGER:: b, s
       CHARACTER(2):: bString
       character(:), allocatable:: bType
 
@@ -991,6 +976,20 @@ MODULE moduleInput
 
       end do
 
+      ! TODO: Move this to the init of species
+      ALLOCATE(qSpecies(1:nSpecies))
+      DO s = 1, nSpecies
+        SELECT TYPE(sp => species(s)%obj)
+        TYPE IS (speciesCharged)
+          qSpecies(s) = sp%q
+
+        CLASS DEFAULT
+          qSpecies(s) = 0.D0
+
+        END SELECT
+
+      END DO
+
     END SUBROUTINE readBoundaryEM
 
     subroutine readPhysicalSurfaces(config)
@@ -1082,6 +1081,7 @@ MODULE moduleInput
       integer:: b, ps, s
       integer:: e
       integer:: nVolColl
+      integer:: boundaryIndex
 
       object = 'geometry'
 
diff --git a/src/modules/mesh/0D/moduleMesh0D.f90 b/src/modules/mesh/0D/moduleMesh0D.f90
index 21bfde8..133ae8e 100644
--- a/src/modules/mesh/0D/moduleMesh0D.f90
+++ b/src/modules/mesh/0D/moduleMesh0D.f90
@@ -43,7 +43,7 @@ MODULE moduleMesh0D
 
       CLASS(meshNode0D), INTENT(out):: self
       INTEGER, INTENT(in):: n
-      REAL(8), INTENT(in):: r(1:3) ! NOTE: Required by interface but unused
+      REAL(8), INTENT(in):: r(1:3) !Unused variable
 
       self%n = n
 
@@ -117,7 +117,7 @@ MODULE moduleMesh0D
     PURE FUNCTION fPsi0D(Xi, nNodes) RESULT(fPsi)
       IMPLICIT NONE
 
-      REAL(8), INTENT(in):: Xi(1:3) ! NOTE: Required by interface but unused
+      REAL(8), INTENT(in):: Xi(1:3)
       INTEGER, INTENT(in):: nNodes
       REAL(8):: fPsi(1:nNodes)
 
@@ -128,7 +128,7 @@ MODULE moduleMesh0D
     PURE FUNCTION dPsi0D(Xi, nNodes) RESULT(dPsi)
       IMPLICIT NONE
 
-      REAL(8), INTENT(in):: Xi(1:3) ! NOTE: Required by interface but unused
+      REAL(8), INTENT(in):: Xi(1:3)
       INTEGER, INTENT(in):: nNodes
       REAL(8):: dPsi(1:3,1:nNodes)
 
@@ -142,7 +142,7 @@ MODULE moduleMesh0D
       CLASS(meshCell0D), INTENT(in):: self
       INTEGER, INTENT(in):: nNodes
       REAL(8), INTENT(in):: dPsi(1:3,1:nNodes)
-      REAL(8):: pDer(1:3, 1:3) ! NOTE: Required by interface but unused
+      REAL(8):: pDer(1:3, 1:3)
 
       pDer = 0.D0
 
@@ -205,7 +205,7 @@ MODULE moduleMesh0D
       IMPLICIT NONE
 
       CLASS(meshCell0D), INTENT(in):: self
-      REAL(8), INTENT(in):: r(1:3) ! NOTE: Required by interface but unused
+      REAL(8), INTENT(in):: r(1:3)
       REAL(8):: xN(1:3)
 
       xN = 0.D0
@@ -215,7 +215,7 @@ MODULE moduleMesh0D
     PURE FUNCTION inside0D(Xi) RESULT(ins)
       IMPLICIT NONE
 
-      REAL(8), INTENT(in):: Xi(1:3) ! NOTE: Required by interface but unused
+      REAL(8), INTENT(in):: Xi(1:3)
       LOGICAL:: ins
 
       ins = .TRUE.
@@ -226,7 +226,7 @@ MODULE moduleMesh0D
       IMPLICIT NONE
 
       CLASS(meshCell0D), INTENT(in):: self
-      REAL(8), INTENT(in):: Xi(1:3) ! NOTE: Required by interface but unused
+      REAL(8), INTENT(in):: Xi(1:3)
       CLASS(meshElement), POINTER, INTENT(out):: neighbourElement
 
       neighbourElement => NULL()
@@ -237,7 +237,7 @@ MODULE moduleMesh0D
     PURE FUNCTION detJ0D(pDer) RESULT(dJ)
       IMPLICIT NONE
 
-      REAL(8), INTENT(in):: pDer(1:3, 1:3) ! NOTE: Required by interface but unused
+      REAL(8), INTENT(in):: pDer(1:3, 1:3)
       REAL(8):: dJ
 
       dJ = 0.D0
@@ -247,7 +247,7 @@ MODULE moduleMesh0D
     PURE FUNCTION invJ0D(pDer) RESULT(invJ)
       IMPLICIT NONE
 
-      REAL(8), INTENT(in):: pDer(1:3, 1:3) ! NOTE: Required by interface but unused
+      REAL(8), INTENT(in):: pDer(1:3, 1:3)
       REAL(8):: invJ(1:3,1:3)
 
       invJ = 0.D0
diff --git a/src/modules/mesh/1DCart/moduleMesh1DCart.f90 b/src/modules/mesh/1DCart/moduleMesh1DCart.f90
index 117b77d..25a3fd7 100644
--- a/src/modules/mesh/1DCart/moduleMesh1DCart.f90
+++ b/src/modules/mesh/1DCart/moduleMesh1DCart.f90
@@ -156,7 +156,7 @@ MODULE moduleMesh1DCart
       IMPLICIT NONE
 
       CLASS(meshEdge1DCart), INTENT(in):: self
-      REAL(8), DIMENSION(1:3), INTENT(in):: r0 ! NOTE: Required by interface but unused
+      REAL(8), DIMENSION(1:3), INTENT(in):: r0
       REAL(8), DIMENSION(1:3):: r
 
       r = (/ self%x, 0.D0, 0.D0 /)
diff --git a/src/modules/mesh/1DRad/moduleMesh1DRad.f90 b/src/modules/mesh/1DRad/moduleMesh1DRad.f90
index 3f6df2b..ca09561 100644
--- a/src/modules/mesh/1DRad/moduleMesh1DRad.f90
+++ b/src/modules/mesh/1DRad/moduleMesh1DRad.f90
@@ -156,7 +156,7 @@ MODULE moduleMesh1DRad
       IMPLICIT NONE
 
       CLASS(meshEdge1DRad), INTENT(in):: self
-      REAL(8), DIMENSION(1:3), INTENT(in):: r0 ! NOTE: Required by interface but unused
+      REAL(8), DIMENSION(1:3), INTENT(in):: r0
       REAL(8), DIMENSION(1:3):: r
 
       r = (/ self%r, 0.D0, 0.D0 /)
diff --git a/src/modules/mesh/inout/0D/moduleMeshInput0D.f90 b/src/modules/mesh/inout/0D/moduleMeshInput0D.f90
index 798d2b2..d968ba3 100644
--- a/src/modules/mesh/inout/0D/moduleMeshInput0D.f90
+++ b/src/modules/mesh/inout/0D/moduleMeshInput0D.f90
@@ -34,7 +34,7 @@ MODULE moduleMeshInput0D
       IMPLICIT NONE
 
       CLASS(meshGeneric), INTENT(inout):: self
-      CHARACTER(:), ALLOCATABLE, INTENT(in):: filename ! NOTE: Required by interface but unused
+      CHARACTER(:), ALLOCATABLE, INTENT(in):: filename !Dummy file, not used
       REAL(8):: r(1:3)
 
       !Allocates one node
diff --git a/src/modules/mesh/inout/text/moduleMeshInputText.f90 b/src/modules/mesh/inout/text/moduleMeshInputText.f90
index c5fe9ef..9946bba 100644
--- a/src/modules/mesh/inout/text/moduleMeshInputText.f90
+++ b/src/modules/mesh/inout/text/moduleMeshInputText.f90
@@ -197,7 +197,6 @@ module moduleMeshInputText
       fileID = 10
 
       open(fileID, file=trim(filename))
-      nNodes = 0
 
       do 
         read(fileID, *, iostat=reason) line
diff --git a/src/modules/mesh/inout/text/moduleMeshOutputText.f90 b/src/modules/mesh/inout/text/moduleMeshOutputText.f90
index e243c80..d69ddda 100644
--- a/src/modules/mesh/inout/text/moduleMeshOutputText.f90
+++ b/src/modules/mesh/inout/text/moduleMeshOutputText.f90
@@ -27,6 +27,7 @@ module moduleMeshOutputText
     subroutine writeCollOutput(self, fileID)
       use moduleMesh
       use moduleCollisions
+      use moduleRefParam, only: L_ref
       implicit none
 
       class(meshGeneric), intent(in):: self
diff --git a/src/modules/mesh/moduleMesh@boundaryParticle.f90 b/src/modules/mesh/moduleMesh@boundaryParticle.f90
index de4edee..9e79835 100644
--- a/src/modules/mesh/moduleMesh@boundaryParticle.f90
+++ b/src/modules/mesh/moduleMesh@boundaryParticle.f90
@@ -99,6 +99,7 @@ submodule(moduleMesh) boundaryParticle
       implicit none
 
       class(boundaryParticleGeneric), allocatable, intent(inout):: boundary
+      integer:: e, et
 
       allocate(boundaryQuasiNeutrality:: boundary)
 
@@ -335,7 +336,6 @@ submodule(moduleMesh) boundaryParticle
       real(8):: EF_dir
       real(8):: alpha
 
-      alpha = 0.85d0
 
       if (associated(edge%e1)) then
         cell => edge%e1
diff --git a/src/modules/mesh/moduleMeshCommon.f90 b/src/modules/mesh/moduleMeshCommon.f90
index e4413aa..ee5f13b 100644
--- a/src/modules/mesh/moduleMeshCommon.f90
+++ b/src/modules/mesh/moduleMeshCommon.f90
@@ -28,7 +28,7 @@ module moduleMeshCommon
     pure function fPsiPoint(Xi, nNodes) RESULT(fPsi)
       implicit none
 
-      real(8), intent(in):: Xi(1:3) ! NOTE: Required by interface but unused
+      real(8), intent(in):: Xi(1:3)
       integer, intent(in):: nNodes
       real(8):: fPsi(1:nNodes)
 
@@ -104,7 +104,7 @@ module moduleMeshCommon
     pure function dPsiSegm(Xi, nNodes) result(dPsi)
       implicit none
 
-      real(8), intent(in):: Xi(1:3) ! NOTE: Required by interface but unused
+      real(8), intent(in):: Xi(1:3)
       integer, intent(in):: nNodes
       real(8):: dPsi(1:3,1:nNodes)
 
@@ -142,7 +142,7 @@ module moduleMeshCommon
     pure function dPsiTria(Xi, nNodes) result(dPsi)
       implicit none
 
-      real(8), intent(in):: Xi(1:3) ! NOTE: Required by interface but unused
+      real(8), intent(in):: Xi(1:3)
       integer, intent(in):: nNodes
       real(8):: dPsi(1:3,1:nNodes)
 
@@ -157,7 +157,7 @@ module moduleMeshCommon
     pure function dPsiTetra(Xi, nNodes) result(dPsi)
       implicit none
 
-      real(8), intent(in):: Xi(1:3) ! NOTE: Required by interface but unused
+      real(8), intent(in):: Xi(1:3)
       integer, intent(in):: nNodes
       real(8):: dPsi(1:3, 1:nNodes)
 
diff --git a/src/modules/moduleCollisions.f90 b/src/modules/moduleCollisions.f90
index dfe2338..36eca2e 100644
--- a/src/modules/moduleCollisions.f90
+++ b/src/modules/moduleCollisions.f90
@@ -553,7 +553,7 @@ MODULE moduleCollisions
       IMPLICIT NONE
 
       CLASS(collisionBinaryChargeExchange), INTENT(in):: self
-      REAL(8), INTENT(in):: vRel ! NOTE: Required by itnerface but unused
+      REAL(8), INTENT(in):: vRel
       TYPE(particle), INTENT(inout), TARGET:: part_i, part_j
 
       SELECT TYPE(sp => part_i%species)
diff --git a/src/modules/moduleInject.f90 b/src/modules/moduleInject.f90
index 324d83b..44c0f80 100644
--- a/src/modules/moduleInject.f90
+++ b/src/modules/moduleInject.f90
@@ -48,7 +48,9 @@ MODULE moduleInject
       REAL(8):: tauInject
       REAL(8), INTENT(in):: flow
       CHARACTER(:), ALLOCATABLE, INTENT(in):: units
-      INTEGER:: e
+      INTEGER:: e, et
+      INTEGER:: phSurface(1:mesh%numEdges)
+      INTEGER:: nVolColl
       REAL(8):: fluxPerStep = 0.D0
 
       self%id          =  i
@@ -68,8 +70,8 @@ MODULE moduleInject
 
       !Calculates total area
       self%surface = 0.D0
-      DO e = 1, self%nEdges
-        self%surface = self%surface + self%edges(e)%obj%surface
+      DO et = 1, self%nEdges
+        self%surface = self%surface + self%edges(et)%obj%surface
 
       END DO
 
@@ -119,8 +121,8 @@ MODULE moduleInject
       IF (particlesPerEdge > 0) THEN
         ! Particles per edge defined by the user
         self%particlesPerEdge = particlesPerEdge
-        DO e = 1, self%nEdges
-          self%weightPerEdge(e) = fluxPerStep * self%edges(e)%obj%surface / real(particlesPerEdge)
+        DO et = 1, self%nEdges
+          self%weightPerEdge(et) = fluxPerStep * self%edges(et)%obj%surface / real(particlesPerEdge)
 
         END DO
 
@@ -129,8 +131,8 @@ MODULE moduleInject
       ELSE
         ! No particles assigned per edge, use the species weight
         self%weightPerEdge = self%species%weight
-        DO e = 1, self%nEdges
-          self%particlesPerEdge(e) = max(1,FLOOR(fluxPerStep*self%edges(e)%obj%surface / self%species%weight))
+        DO et = 1, self%nEdges
+          self%particlesPerEdge(et) = max(1,FLOOR(fluxPerStep*self%edges(et)%obj%surface / self%species%weight))
         END DO
 
         self%nParticles = SUM(self%particlesPerEdge)
diff --git a/src/modules/moduleProbe.f90 b/src/modules/moduleProbe.f90
index 1e1339e..754d56a 100644
--- a/src/modules/moduleProbe.f90
+++ b/src/modules/moduleProbe.f90
@@ -253,6 +253,7 @@ MODULE moduleProbe
     END SUBROUTINE doProbes
 
     SUBROUTINE outputProbes()
+      USE moduleCaseParam, ONLY: timeStep
       IMPLICIT NONE
 
       INTEGER:: i
diff --git a/src/modules/solver/electromagnetic/moduleEM.f90 b/src/modules/solver/electromagnetic/moduleEM.f90
index 4b115aa..b8d9b55 100644
--- a/src/modules/solver/electromagnetic/moduleEM.f90
+++ b/src/modules/solver/electromagnetic/moduleEM.f90
@@ -121,6 +121,7 @@ MODULE moduleEM
       REAL(8), INTENT(in):: phi(1:n)
       REAL(8):: n_e(1:n)
       REAL(8):: n_e0 = 1.0D16, phi_0 = -500.0D0, T_e = 11604.0
+      INTEGER:: i
 
       n_e = n_e0 / n_ref * exp(qe * (phi*Volt_ref - phi_0) / (kb * T_e))
 
diff --git a/src/modules/solver/pusher/modulePusher.f90 b/src/modules/solver/pusher/modulePusher.f90
index a7dac57..f75c733 100644
--- a/src/modules/solver/pusher/modulePusher.f90
+++ b/src/modules/solver/pusher/modulePusher.f90
@@ -231,8 +231,8 @@ MODULE modulePusher
       USE moduleSpecies
       IMPLICIT NONE
 
-      TYPE(particle), INTENT(inout):: part ! NOTE: Required by interface but unused
-      REAL(8), INTENT(in):: tauIn ! NOTE: Required by interface but unused
+      TYPE(particle), INTENT(inout):: part
+      REAL(8), INTENT(in):: tauIn
 
     END SUBROUTINE push0D