From 6614482a641f4aaf5ccf53641429134b7c8298c1 Mon Sep 17 00:00:00 2001
From: Alexander Schmidt <alex@treefish.org>
Date: Mon, 18 Jul 2016 22:25:57 +0200
Subject: [PATCH] Fixed them all

---
 CMakeLists.txt |  2 +-
 complexpair.h  |  3 ---
 sea.cpp        | 42 +++++++++++++++++++++++-------------------
 sea.h          |  6 ++----
 4 files changed, 26 insertions(+), 27 deletions(-)
 delete mode 100644 complexpair.h

diff --git a/CMakeLists.txt b/CMakeLists.txt
index 986edd4..b6b7e6b 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -18,4 +18,4 @@ add_executable(seamulator seamulator.cpp
 			  seaview.cpp)
 
 target_link_libraries(seamulator ${OPENGL_LIBRARIES} ${GLUT_LIBRARY}
-				 ${FFTW_LIBRARIES} ${FFTW_LIBRARIES_THREADS})
+				 ${FFTW_LIBRARIES})
diff --git a/complexpair.h b/complexpair.h
deleted file mode 100644
index 9eea2a4..0000000
--- a/complexpair.h
+++ /dev/null
@@ -1,3 +0,0 @@
-#pragma once
-
-using ComplexPair = std::pair<std::complex<double>, std::complex<double>>;
diff --git a/sea.cpp b/sea.cpp
index b76a9dd..98c5c4b 100644
--- a/sea.cpp
+++ b/sea.cpp
@@ -6,7 +6,7 @@
 
 #include "watersurface.h"
 
-const double Sea::PHILLIPS_CONSTANT{0.0000003};
+const double Sea::PHILLIPS_CONSTANT{0.0000002};
 const double Sea::GRAVITATIONAL_CONSTANT{9.8};
 
 Sea::Sea(WaterSurfacePtr surface) :
@@ -16,12 +16,9 @@ Sea::Sea(WaterSurfacePtr surface) :
   m_randomGenerator{m_randomDevice()},
   m_normalDistribution{0.0, 1.0}
 {
-  m_fourierAmplitudes.resize(pow(m_surface->size(), 2));
+  m_fourierAmplitudes.resize(pow(m_surface->size() + 1, 2));
   generateFourierAmplitudes();
 
-  fftw_init_threads();
-  fftw_plan_with_nthreads(4);
-
   m_fftwIn = (fftw_complex*)
     fftw_malloc(sizeof(fftw_complex) * pow(m_surface->size(), 2));
   m_fftwOut = (fftw_complex*)
@@ -63,8 +60,8 @@ void Sea::update()
       const double omega = sqrt(GRAVITATIONAL_CONSTANT * k);
 
       std::complex<double> amplitude =
-	fourierAmplitudeAt(n, m).first * exp(1i * omega * runtime) +
-	fourierAmplitudeAt(n, m).second * exp(-1i * omega * runtime);
+       	fourierAmplitudeAt(n, m) * exp(1i * omega * runtime) +
+       	std::conj(fourierAmplitudeAt(-n, -m)) * exp(-1i * omega * runtime);
 
       const int positiveN = (n + m_surface->size()) % m_surface->size();
       int fftwIndex = positiveM + positiveN * m_surface->size();
@@ -96,16 +93,7 @@ double Sea::phillipsSpectrum(double k_x, double k_y) const
     cosineFactor;
 }
 
-ComplexPair Sea::generateFourierAmplitude(double k_x, double k_y)
-{
-  std::complex<double> cDist(m_normalDistribution(m_randomGenerator),
-			     m_normalDistribution(m_randomGenerator));
-
-  return {cDist * sqrt(phillipsSpectrum(k_x, k_y)) / sqrt(2),
-      std::conj(cDist) * sqrt(phillipsSpectrum(-k_x, -k_y)) / sqrt(2)};
-}
-
-ComplexPair& Sea::fourierAmplitudeAt(int n, int m)
+std::complex<double>& Sea::fourierAmplitudeAt(int n, int m)
 {
   return m_fourierAmplitudes.at
     (n + m_surface->size()/2 +
@@ -120,12 +108,28 @@ double Sea::spatialFrequencyForIndex(int n) const
 void Sea::generateFourierAmplitudes()
 {
   for (int m = -m_surface->size()/2; m < m_surface->size()/2; ++m) {
+    const double k_y = spatialFrequencyForIndex(m);
+
     for (int n = -m_surface->size()/2; n < m_surface->size()/2; ++n) {
+      const double k_x = spatialFrequencyForIndex(n);
+
+      std::complex<double> cDist(m_normalDistribution(m_randomGenerator),
+				 m_normalDistribution(m_randomGenerator));
+
       fourierAmplitudeAt(n, m) =
-	generateFourierAmplitude(spatialFrequencyForIndex(n),
-				 spatialFrequencyForIndex(m));
+	cDist * sqrt(phillipsSpectrum(k_x, k_y)) / sqrt(2);
     }
   }
 
+  for (int n = -m_surface->size()/2; n < m_surface->size()/2; ++n) {
+    fourierAmplitudeAt(n, m_surface->size()/2) =
+      fourierAmplitudeAt(n, -m_surface->size()/2);
+  }
+
+  for (int m = -m_surface->size()/2; m < m_surface->size()/2; ++m) {
+    fourierAmplitudeAt(m_surface->size()/2, m) =
+      fourierAmplitudeAt(-m_surface->size()/2, m);
+  }
+
   fourierAmplitudeAt(0, 0) = {0, 0};
 }
diff --git a/sea.h b/sea.h
index 20f5d63..798b92d 100644
--- a/sea.h
+++ b/sea.h
@@ -8,7 +8,6 @@
 
 #include "seafwd.h"
 
-#include "complexpair.h"
 #include "watersurfacefwd.h"
 
 class Sea
@@ -30,14 +29,13 @@ class Sea
   std::random_device m_randomDevice;
   std::mt19937 m_randomGenerator;
   std::normal_distribution<> m_normalDistribution;
-  std::vector<ComplexPair> m_fourierAmplitudes;
+  std::vector<std::complex<double>> m_fourierAmplitudes;
   fftw_complex *m_fftwIn, *m_fftwOut;
   fftw_plan m_fftwPlan;
   std::chrono::time_point<std::chrono::system_clock> m_startTime;
 
   double phillipsSpectrum(double k_x, double k_y) const;
-  ComplexPair generateFourierAmplitude(double k_x, double k_y);
-  ComplexPair& fourierAmplitudeAt(int n, int m);
+  std::complex<double>& fourierAmplitudeAt(int n, int m);
   void generateFourierAmplitudes();
   double spatialFrequencyForIndex(int n) const;
   double getRuntime() const;
-- 
2.39.2