Skip to content

GitLab

P
pypocketfft
Commit 68a9c113 authored by Peter Bell's avatar Peter Bell Committed by Martin Reinecke
Browse files
Options

Simplify general_ functions handling of inplace transforms

parent 50476533
Hide whitespace changes
Inline Side-by-side
Showing with 30 additions and 69 deletions
pocketfft_hdronly.h
View file @ 68a9c113
...@@ -2902,25 +2902,16 @@ template<typename T> POCKETFFT_NOINLINE void general_c( ...@@ -2902,25 +2902,16 @@ template<typename T> POCKETFFT_NOINLINE void general_c(
while (it.remaining()>0) while (it.remaining()>0)
{ {
it.advance(1); it.advance(1);
auto tdata = reinterpret_cast<cmplx<T> *>(storage.data()); auto buf = it.stride_out() == sizeof(cmplx<T>) ?
if ((&tin[0]==&out[0]) && (it.stride_out()==sizeof(cmplx<T>))) // fully in-place &out[it.oofs(0)] : reinterpret_cast<cmplx<T> *>(storage.data());
forward ? plan->forward (&out[it.oofs(0)], fct)
: plan->backward(&out[it.oofs(0)], fct); if (buf != &tin[it.iofs(0)])
else if (it.stride_out()==sizeof(cmplx<T>)) // compute FFT in output location
{
for (size_t i=0; i<len; ++i)
out[it.oofs(i)] = tin[it.iofs(i)];
forward ? plan->forward (&out[it.oofs(0)], fct)
: plan->backward(&out[it.oofs(0)], fct);
}
else
{
for (size_t i=0; i<len; ++i) for (size_t i=0; i<len; ++i)
tdata[i] = tin[it.iofs(i)]; buf[i] = tin[it.iofs(i)];
forward ? plan->forward (tdata, fct) : plan->backward(tdata, fct); forward ? plan->forward (buf, fct) : plan->backward(buf, fct);
if (buf != &out[it.oofs(0)])
for (size_t i=0; i<len; ++i) for (size_t i=0; i<len; ++i)
out[it.oofs(i)] = tdata[i]; out[it.oofs(i)] = buf[i];
}
} }
} // end of parallel region } // end of parallel region
fct = T(1); // factor has been applied, use 1 for remaining axes fct = T(1); // factor has been applied, use 1 for remaining axes
...@@ -3034,23 +3025,16 @@ template<typename Trafo, typename T> POCKETFFT_NOINLINE void general_dcst( ...@@ -3034,23 +3025,16 @@ template<typename Trafo, typename T> POCKETFFT_NOINLINE void general_dcst(
while (it.remaining()>0) while (it.remaining()>0)
{ {
it.advance(1); it.advance(1);
auto tdata = reinterpret_cast<T *>(storage.data()); auto buf = it.stride_out() == sizeof(T) ? &out[it.oofs(0)]
if ((&tin[0]==&out[0]) && (it.stride_out()==sizeof(T))) // fully in-place : reinterpret_cast<T *>(storage.data());
plan->exec(&out[it.oofs(0)], fct, ortho, type, cosine);
else if (it.stride_out()==sizeof(T)) // compute FFT in output location if (buf != &tin[it.iofs(0)])
{
for (size_t i=0; i<len; ++i)
out[it.oofs(i)] = tin[it.iofs(i)];
plan->exec(&out[it.oofs(0)], fct, ortho, type, cosine);
}
else
{
for (size_t i=0; i<len; ++i) for (size_t i=0; i<len; ++i)
tdata[i] = tin[it.iofs(i)]; buf[i] = tin[it.iofs(i)];
plan->exec(tdata, fct, ortho, type, cosine); plan->exec(buf, fct, ortho, type, cosine);
if (buf != &out[it.oofs(0)])
for (size_t i=0; i<len; ++i) for (size_t i=0; i<len; ++i)
out[it.oofs(i)] = tdata[i]; out[it.oofs(i)] = buf[i];
}
} }
} // end of parallel region } // end of parallel region
fct = T(1); // factor has been applied, use 1 for remaining axes fct = T(1); // factor has been applied, use 1 for remaining axes
...@@ -3242,45 +3226,22 @@ template<typename T> POCKETFFT_NOINLINE void general_r( ...@@ -3242,45 +3226,22 @@ template<typename T> POCKETFFT_NOINLINE void general_r(
while (it.remaining()>0) while (it.remaining()>0)
{ {
it.advance(1); it.advance(1);
auto tdata = reinterpret_cast<T *>(storage.data()); auto buf = it.stride_out() == sizeof(T) ?
if ((&tin[0]==&out[0]) && (it.stride_out()==sizeof(T))) // fully in-place &out[it.oofs(0)] : reinterpret_cast<T *>(storage.data());
{
if ((!r2c) && forward) if (buf != &tin[it.iofs(0)])
for (size_t i=2; i<len; i+=2)
out[it.oofs(i)] = -out[it.oofs(i)];
forward ? plan->forward (&out[it.oofs(0)], fct)
: plan->backward(&out[it.oofs(0)], fct);
if (r2c && (!forward))
for (size_t i=2; i<len; i+=2)
out[it.oofs(i)] = -out[it.oofs(i)];
}
else if (it.stride_out()==sizeof(T)) // compute FFT in output location
{
for (size_t i=0; i<len; ++i)
out[it.oofs(i)] = tin[it.iofs(i)];
if ((!r2c) && forward)
for (size_t i=2; i<len; i+=2)
out[it.oofs(i)] = -out[it.oofs(i)];
forward ? plan->forward (&out[it.oofs(0)], fct)
: plan->backward(&out[it.oofs(0)], fct);
if (r2c && (!forward))
for (size_t i=2; i<len; i+=2)
out[it.oofs(i)] = -out[it.oofs(i)];
}
else
{
for (size_t i=0; i<len; ++i) for (size_t i=0; i<len; ++i)
tdata[i] = tin[it.iofs(i)]; buf[i] = tin[it.iofs(i)];
if ((!r2c) && forward) if ((!r2c) && forward)
for (size_t i=2; i<len; i+=2) for (size_t i=2; i<len; i+=2)
tdata[i] = -tdata[i]; buf[i] = -buf[i];
forward ? plan->forward (tdata, fct) : plan->backward(tdata, fct); forward ? plan->forward(buf, fct) : plan->backward(buf, fct);
if (r2c && (!forward)) if (r2c && (!forward))
for (size_t i=2; i<len; i+=2) for (size_t i=2; i<len; i+=2)
tdata[i] = -tdata[i]; buf[i] = -buf[i];
if (buf != &out[it.oofs(0)])
for (size_t i=0; i<len; ++i) for (size_t i=0; i<len; ++i)
out[it.oofs(i)] = tdata[i]; out[it.oofs(i)] = buf[i];
}
} }
} // end of parallel region } // end of parallel region
fct = T(1); // factor has been applied, use 1 for remaining axes fct = T(1); // factor has been applied, use 1 for remaining axes
......
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