Fixed large bug, small optimizations

1 files changed, 5 insertions, 13 deletions

diff --git a/fourier_integrator.cpp b/fourier_integrator.cpp
index 7fc8ab4..829a074 100644
--- a/fourier_integrator.cpp
+++ b/fourier_integrator.cpp
@@ -152,21 +152,18 @@ int main(int argc, char* argv[]) {
       }
       std::vector<Complex> dfCt = fft.fourier(dfC);
 
-      std::vector<Complex> Ctnew(Ct.size());
-      std::vector<Complex> Rtnew(Rt.size());
-
       for (unsigned i = 0; i < Rt.size(); i++) {
         Real ω = i * Δω;
-        Rtnew[i] = (1.0 + pow(y, 2) * RddfCt[i] * Rt[i]) / (z + 1i * ω);
+        Rt[i] = (1.0 + pow(y, 2) * RddfCt[i] * Rt[i]) / (z + 1i * ω);
       }
 
-      for (unsigned i = 0; i < Rt.size(); i++) {
+      for (unsigned i = 0; i < Ct.size(); i++) {
         Real ω = i * Δω;
-        Ctnew[i] = (Γ₀ * std::conj(Rt[i]) / (1 + pow(τ₀ * ω, 2)) + pow(y, 2) * (RddfCt[i] * Ct[i] + dfCt[i] * std::conj(Rt[i]))) / (z + 1i * ω);
+        Ct[i] = (Γ₀ * std::conj(Rt[i]) / (1 + pow(τ₀ * ω, 2)) + pow(y, 2) * (RddfCt[i] * Ct[i] + dfCt[i] * std::conj(Rt[i]))) / (z + 1i * ω);
       }
 
-      std::vector<Real> Cnew = fft.inverse(Ctnew);
-      std::vector<Real> Rnew = fft.inverse(Rtnew);
+      std::vector<Real> Cnew = fft.inverse(Ct);
+      std::vector<Real> Rnew = fft.inverse(Rt);
 
       ΔC = 0;
       for (unsigned i = 0; i < Cnew.size(); i++) {
@@ -174,7 +171,6 @@ int main(int argc, char* argv[]) {
         ΔC += pow(Rnew[i] - R[i], 2);
       }
 
-      /*
       for (unsigned i = 0; i < Cnew.size(); i++) {
         C[i] += γ * (Cnew[i] - C[i]);
       }
@@ -182,10 +178,6 @@ int main(int argc, char* argv[]) {
       for (unsigned i = 0; i < Rnew.size(); i++) {
         R[i] += γ * (Rnew[i] - R[i]);
       }
-      */
-
-      C = Cnew;
-      R = Rnew;
 
       z *= Cnew[0];