diff options
-rw-r--r-- | fourier.cpp | 5 | ||||
-rw-r--r-- | fourier.hpp | 2 | ||||
-rw-r--r-- | fourier_integrator.cpp | 8 | ||||
-rw-r--r-- | get_energy.cpp | 2 |
4 files changed, 9 insertions, 8 deletions
diff --git a/fourier.cpp b/fourier.cpp index 8943521..d2f9ad8 100644 --- a/fourier.cpp +++ b/fourier.cpp @@ -85,8 +85,9 @@ std::tuple<std::vector<Complex>, std::vector<Complex>> RddfCtdfCt(FourierTransfo return {RddfCt, dfCt}; } -Real estimateZ(FourierTransform& fft, const std::vector<Real>& C, const std::vector<Complex>& Ct, const std::vector<Real>& R, const std::vector<Complex>& Rt, unsigned p, unsigned s, Real λ, Real y) { +Real estimateZ(FourierTransform& fft, const std::vector<Real>& C, const std::vector<Complex>& Ct, const std::vector<Real>& R, const std::vector<Complex>& Rt, unsigned p, unsigned s, Real λ, Real τ₀, Real y) { auto [RddfCt, dfCt] = RddfCtdfCt(fft, C, R, p, s, λ); + Real Γ₀ = 1 + τ₀ / 2; - return ((std::conj(Rt[0]) + pow(y, 2) * (RddfCt[0] * Ct[0] + dfCt[0] * std::conj(Rt[0]))) / Ct[0]).real(); + return ((Γ₀ * std::conj(Rt[0]) + pow(y, 2) * (RddfCt[0] * Ct[0] + dfCt[0] * std::conj(Rt[0]))) / Ct[0]).real(); } diff --git a/fourier.hpp b/fourier.hpp index 086660a..c4c1cf2 100644 --- a/fourier.hpp +++ b/fourier.hpp @@ -31,4 +31,4 @@ public: std::string fourierFile(std::string prefix, unsigned p, unsigned s, Real λ, Real τ₀, Real y, unsigned log2n, Real τₘₐₓ); Real energy(const std::vector<Real>& C, const std::vector<Real>& R, unsigned p, unsigned s, Real λ, Real y, Real Δτ); std::tuple<std::vector<Complex>, std::vector<Complex>> RddfCtdfCt(FourierTransform& fft, const std::vector<Real>& C, const std::vector<Real>& R, unsigned p, unsigned s, Real λ); -Real estimateZ(FourierTransform& fft, const std::vector<Real>& C, const std::vector<Complex>& Ct, const std::vector<Real>& R, const std::vector<Complex>& Rt, unsigned p, unsigned s, Real λ, Real y); +Real estimateZ(FourierTransform& fft, const std::vector<Real>& C, const std::vector<Complex>& Ct, const std::vector<Real>& R, const std::vector<Complex>& Rt, unsigned p, unsigned s, Real λ, Real τ₀, Real y); diff --git a/fourier_integrator.cpp b/fourier_integrator.cpp index d83fb9d..24a7f1d 100644 --- a/fourier_integrator.cpp +++ b/fourier_integrator.cpp @@ -17,7 +17,7 @@ int main(int argc, char* argv[]) { unsigned maxIterations = 1000; Real ε = 1e-14; - Real γ = 0; + Real γ = 1; bool loadData = false; @@ -68,8 +68,8 @@ int main(int argc, char* argv[]) { Real Δτ = τₘₐₓ / M_PI / n; Real Δω = M_PI / τₘₐₓ; - Real z = (-1+sqrt(1+2*τ₀)) / (2 * τ₀); - Real Γ₀ = 1; + Real z = 0.5; + Real Γ₀ = 1 + τ₀ / 2; std::vector<Real> C(2 * n); std::vector<Real> R(2 * n); @@ -103,7 +103,7 @@ int main(int argc, char* argv[]) { Ct = fft.fourier(C); Rt = fft.fourier(R); - z = estimateZ(fft, C, Ct, R, Rt, p, s, λ, y); + z = estimateZ(fft, C, Ct, R, Rt, p, s, λ, τ₀, y); } while (y += Δy, y <= yₘₐₓ) { diff --git a/get_energy.cpp b/get_energy.cpp index b6e835f..6e14d9b 100644 --- a/get_energy.cpp +++ b/get_energy.cpp @@ -75,7 +75,7 @@ int main(int argc, char* argv[]) { std::vector<Complex> Ct = fft.fourier(C); std::vector<Complex> Rt = fft.fourier(R); - Real z = estimateZ(fft, C, Ct, R, Rt, p, s, λ, y); + Real z = estimateZ(fft, C, Ct, R, Rt, p, s, λ, τ₀, y); std::cout << y << " " << e << " " << Ct[0].real() << " " << z << std::endl; } |