Recovered the ouflowBoundary code

src/modules/init/moduleInput.f90

@@ -900,6 +900,9 @@ MODULE moduleInput
 
               call initQuasiNeutrality(bound, s_incident)
 
+            CASE('outflowAdaptive')
+              ALLOCATE(boundaryOutflowAdaptive:: bound)
+
             CASE DEFAULT
               CALL criticalError('Boundary type ' // bType // ' undefined', 'readBoundary')
 

src/modules/mesh/moduleMesh.f90

@@ -789,7 +789,15 @@ MODULE moduleMesh
 
   end type boundaryQuasiNeutrality
 
-  !Wrapper for boundary types (one per species)
+  !Boundary for quasi-neutral outflow adjusting reflection coefficient
+  type, public, extends(boundaryParticleGeneric):: boundaryOutflowAdaptive
+    real(8):: outflowCurrent
+    real(8):: reflectionFraction
+    contains
+      procedure, pass:: apply => outflowAdaptive
+
+  end type boundaryOutflowAdaptive
+
   interface
     module subroutine reflection(self, edge, part)
       use moduleSpecies
@@ -854,6 +862,15 @@ MODULE moduleMesh
 
     end subroutine quasiNeutrality
 
+    module subroutine outflowAdaptive(self, edge, part)
+      use moduleSpecies
+
+      class(boundaryOutflowAdaptive), intent(inout):: self
+      class(meshEdge), intent(inout):: edge
+      class(particle), intent(inout):: part
+
+    end subroutine outflowAdaptive
+
     ! Generic basic boundary conditions to use internally in the code
     module subroutine genericReflection(edge, part) 
       use moduleSpecies

src/modules/mesh/moduleMesh@boundaryParticle.f90

@@ -487,6 +487,36 @@ submodule(moduleMesh) boundaryParticle
 
     end subroutine quasiNeutrality_print
 
+    ! outflowAdaptive
+    ! Adjust the reflection coefficient of the boundary to maintain a quasi-neutral outflow
+    module subroutine outflowAdaptive(self, edge, part)
+      use moduleSpecies
+      use moduleRefParam, only: v_ref
+      implicit none
+
+      class(boundaryOutflowAdaptive), intent(inout):: self
+      class(meshEdge), intent(inout):: edge
+      class(particle), intent(inout):: part
+      real(8):: v_cut
+
+      v_cut = 2.d0 * 40.d3/v_ref !This will be the drag velocity of the ions in the future
+
+      if (dot_product(part%v,-edge%normal) > v_cut) then
+        ! print *,part%v(1), v_cut
+        ! part%v = part%v + v_cut*edge%normal
+        part%v(1) = part%v(1) - v_cut
+        ! print *,part%v(1)
+        call genericReflection(edge, part)
+        ! print *,part%v(1)
+        ! print *
+
+      else
+        call genericTransparent(edge, part)
+
+      end if
+
+    end subroutine outflowAdaptive
+
     ! Generic boundary conditions for internal use
     ! reflection
     module subroutine genericReflection(edge, part)

src/modules/moduleInject.f90

@@ -281,7 +281,7 @@ MODULE moduleInject
                                                      self%v(3)%obj%randomVel() /)
             !If injecting a no-drift distribution and velocity is negative, reflect
             if ((self%vMod == 0.D0) .and. &
-                (dot_product(direction, partInj(n)%v) < 0.D0)) then
+                (dot_product(partInj(n)%v, direction) <= 0.D0)) then
               partInj(n)%v = - partInj(n)%v
 
             end if