Commit 307544f5 authored by Martin Reinecke's avatar Martin Reinecke
Browse files

Merge branch 'incremental_updates' into 'master'

synchronize with pocketfft C++

See merge request !43
parents 61b41fca 9f9338fc
/* /*
This file is part of pocketfft. This file is part of pocketfft.
Copyright (C) 2010-2020 Max-Planck-Society Copyright (C) 2010-2021 Max-Planck-Society
Copyright (C) 2019-2020 Peter Bell Copyright (C) 2019-2020 Peter Bell
For the odd-sized DCT-IV transforms: For the odd-sized DCT-IV transforms:
...@@ -52,12 +52,12 @@ SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. ...@@ -52,12 +52,12 @@ SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
#endif #endif
#include <cmath> #include <cmath>
#include <cstring>
#include <cstdlib> #include <cstdlib>
#include <stdexcept> #include <stdexcept>
#include <memory> #include <memory>
#include <vector> #include <vector>
#include <complex> #include <complex>
#include <algorithm>
#if POCKETFFT_CACHE_SIZE!=0 #if POCKETFFT_CACHE_SIZE!=0
#include <array> #include <array>
#include <mutex> #include <mutex>
...@@ -1302,9 +1302,9 @@ template<bool fwd, typename T> void pass11 (size_t ido, size_t l1, ...@@ -1302,9 +1302,9 @@ template<bool fwd, typename T> void pass11 (size_t ido, size_t l1,
} }
} }
#undef PARTSTEP11 #undef POCKETFFT_PARTSTEP11
#undef PARTSTEP11a0 #undef POCKETFFT_PARTSTEP11a0
#undef PARTSTEP11a #undef POCKETFFT_PARTSTEP11a
#undef POCKETFFT_PREP11 #undef POCKETFFT_PREP11
template<bool fwd, typename T> void passg (size_t ido, size_t ip, template<bool fwd, typename T> void passg (size_t ido, size_t ip,
...@@ -1456,7 +1456,7 @@ template<bool fwd, typename T> void pass_all(T c[], T0 fct) const ...@@ -1456,7 +1456,7 @@ template<bool fwd, typename T> void pass_all(T c[], T0 fct) const
for (size_t i=0; i<length; ++i) for (size_t i=0; i<length; ++i)
c[i] = ch[i]*fct; c[i] = ch[i]*fct;
else else
memcpy (c,p1,length*sizeof(T)); std::copy_n (p1, length, c);
} }
else else
if (fct!=1.) if (fct!=1.)
...@@ -2204,19 +2204,19 @@ template<typename T> void radbg(size_t ido, size_t ip, size_t l1, ...@@ -2204,19 +2204,19 @@ template<typename T> void radbg(size_t ido, size_t ip, size_t l1,
} }
} }
template<typename T> void copy_and_norm(T *c, T *p1, size_t n, T0 fct) const template<typename T> void copy_and_norm(T *c, T *p1, T0 fct) const
{ {
if (p1!=c) if (p1!=c)
{ {
if (fct!=1.) if (fct!=1.)
for (size_t i=0; i<n; ++i) for (size_t i=0; i<length; ++i)
c[i] = fct*p1[i]; c[i] = fct*p1[i];
else else
memcpy (c,p1,n*sizeof(T)); std::copy_n (p1, length, c);
} }
else else
if (fct!=1.) if (fct!=1.)
for (size_t i=0; i<n; ++i) for (size_t i=0; i<length; ++i)
c[i] *= fct; c[i] *= fct;
} }
...@@ -2224,16 +2224,16 @@ template<typename T> void radbg(size_t ido, size_t ip, size_t l1, ...@@ -2224,16 +2224,16 @@ template<typename T> void radbg(size_t ido, size_t ip, size_t l1,
template<typename T> void exec(T c[], T0 fct, bool r2hc) const template<typename T> void exec(T c[], T0 fct, bool r2hc) const
{ {
if (length==1) { c[0]*=fct; return; } if (length==1) { c[0]*=fct; return; }
size_t n=length, nf=fact.size(); size_t nf=fact.size();
arr<T> ch(n); arr<T> ch(length);
T *p1=c, *p2=ch.data(); T *p1=c, *p2=ch.data();
if (r2hc) if (r2hc)
for(size_t k1=0, l1=n; k1<nf;++k1) for(size_t k1=0, l1=length; k1<nf;++k1)
{ {
size_t k=nf-k1-1; size_t k=nf-k1-1;
size_t ip=fact[k].fct; size_t ip=fact[k].fct;
size_t ido=n / l1; size_t ido=length / l1;
l1 /= ip; l1 /= ip;
if(ip==4) if(ip==4)
radf4(ido, l1, p1, p2, fact[k].tw); radf4(ido, l1, p1, p2, fact[k].tw);
...@@ -2251,7 +2251,7 @@ template<typename T> void radbg(size_t ido, size_t ip, size_t l1, ...@@ -2251,7 +2251,7 @@ template<typename T> void radbg(size_t ido, size_t ip, size_t l1,
for(size_t k=0, l1=1; k<nf; k++) for(size_t k=0, l1=1; k<nf; k++)
{ {
size_t ip = fact[k].fct, size_t ip = fact[k].fct,
ido= n/(ip*l1); ido= length/(ip*l1);
if(ip==4) if(ip==4)
radb4(ido, l1, p1, p2, fact[k].tw); radb4(ido, l1, p1, p2, fact[k].tw);
else if(ip==2) else if(ip==2)
...@@ -2266,7 +2266,7 @@ template<typename T> void radbg(size_t ido, size_t ip, size_t l1, ...@@ -2266,7 +2266,7 @@ template<typename T> void radbg(size_t ido, size_t ip, size_t l1,
l1*=ip; l1*=ip;
} }
copy_and_norm(c,p1,n,fct); copy_and_norm(c,p1,fct);
} }
private: private:
...@@ -2437,13 +2437,12 @@ template<typename T0> class fftblue ...@@ -2437,13 +2437,12 @@ template<typename T0> class fftblue
tmp[m].Set(c[m], zero); tmp[m].Set(c[m], zero);
fft<true>(tmp.data(),fct); fft<true>(tmp.data(),fct);
c[0] = tmp[0].r; c[0] = tmp[0].r;
memcpy (c+1, tmp.data()+1, (n-1)*sizeof(T)); std::copy_n (&tmp[1].r, n-1, &c[1]);
} }
else else
{ {
tmp[0].Set(c[0],c[0]*0); tmp[0].Set(c[0],c[0]*0);
memcpy (reinterpret_cast<void *>(tmp.data()+1), std::copy_n (c+1, n-1, &tmp[1].r);
reinterpret_cast<void *>(c+1), (n-1)*sizeof(T));
if ((n&1)==0) tmp[n/2].i=T0(0)*c[0]; if ((n&1)==0) tmp[n/2].i=T0(0)*c[0];
for (size_t m=1; 2*m<n; ++m) for (size_t m=1; 2*m<n; ++m)
tmp[n-m].Set(tmp[m].r, -tmp[m].i); tmp[n-m].Set(tmp[m].r, -tmp[m].i);
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment