diff --git a/src/modules/BundleIOModule.f90 b/src/modules/BundleIOModule.f90
index 541970a..e0b6b58 100644
--- a/src/modules/BundleIOModule.f90
+++ b/src/modules/BundleIOModule.f90
@@ -57,6 +57,7 @@ subroutine FMS_WriteFBranching(B1)
       call FMS_Branching(B1, Prob)
 
       write (IUnit, '(1x,f11.2,50(1x,f11.9))') B1%CurrentTime, (Prob(i), i=1, nstate), FMS_Norm(B1)
+      flush (IUnit)
 
    end subroutine FMS_WriteFBranching
 
@@ -100,6 +101,7 @@ subroutine FMS_WriteFEnergy(B1)
 
       ! write energies
       write (IUnit, '(f10.2,6(1x,f14.9))') B1%CurrentTime, PotQM, KinQM, PotQM + KinQM, PotCl, KinCl, PotCl + KinCl
+      flush (IUnit)
 
    end subroutine FMS_WriteFEnergy
 
@@ -182,6 +184,7 @@ subroutine FMS_WriteFCorr(B1)
 
       Corr = overlap_bundle(B1, B_init)
       write (iUnit, '(f10.2,4(1x,f12.9))') B1%CurrentTime, abs(Corr)**2, Corr
+      flush (IUnit)
 
    contains
 
diff --git a/src/modules/TrajectoryIOModule.f90 b/src/modules/TrajectoryIOModule.f90
index 8e4c701..441a2d3 100644
--- a/src/modules/TrajectoryIOModule.f90
+++ b/src/modules/TrajectoryIOModule.f90
@@ -209,6 +209,7 @@ subroutine FMS_WriteFTrajDump(T)
       write (iunit, 2) T%get_time(), (T%Particle(i)%get_pos(), i=1, T%NumParticles), &
          ((T%Particle(i)%get_mom(j), j=1, T%Particle(i)%NumDimensions), i=1, T%NumParticles), &
          T%Phase, real(T%Amplitude), aimag(T%Amplitude), FMS_Weight(T), dble(T%StateID)
+      flush (iunit)
 
    end subroutine FMS_WriteFTrajDump
 !>
@@ -235,6 +236,7 @@ subroutine FMS_WriteFAmp(T)
       end if
 
       write (iunit, '(f10.2,4(f10.4))') T%get_time(), FMS_Weight(T), T%Amplitude
+      flush (iunit)
 
    end subroutine FMS_WriteFAmp
 !>
@@ -548,6 +550,7 @@ subroutine FMS_WriteFCouple(T)
       ! 4. Write the couplings
       write (IUnit, 2) T%get_time(), (sqrt(sum(FMS_Coupling(T, i, j)**2)), j=1, nstate), &
          (FMS_CoupDotVel(T, j), j=1, nstate)
+      flush (IUnit)
 
    end subroutine FMS_WriteFCouple
 !>
@@ -613,6 +616,7 @@ subroutine FMS_WriteFSOCouple(T)
       end if
 
       write (IUnit, 2) T%get_time(), Coup(1:k)
+      flush (IUnit)
 
    end subroutine FMS_WriteFSOCouple
 !>
@@ -916,6 +920,7 @@ subroutine FMS_WriteFTDipole(T)
       ! 4. Write transition dipole
       write (iunit, 2) T%get_time(), (sqrt(sum(FMS_TransDipole(T, j)**2)), j=2, nstate), &
          (FMS_TransDipole(T, j), j=2, nstate)
+      flush (iunit)
 
    end subroutine FMS_WriteFTDipole
 !>
@@ -966,6 +971,7 @@ subroutine FMS_WriteFPotEn(T)
       else
          write (iunit, 2) time, (Potential(T, j), j=1, nstate), Kinetic(T) + Potential(T, i)
       end if
+      flush (iunit)
 
    end subroutine FMS_WriteFPotEn
 !>