Commit 34dddd25 authored by Martin Reinecke's avatar Martin Reinecke
Browse files

tweaks

parent bddd29f1
......@@ -3,6 +3,12 @@
#include <iostream>
#include <algorithm>
#ifdef __GNUC__
#define RESTRICT __restrict__
#else
#define RESTRICT
#endif
using namespace std;
namespace py = pybind11;
......@@ -403,13 +409,17 @@ a_c_c to_grid (const a_d_c &uv_, const a_c_c &vis_,
for (int ipart=0; ipart<nvis; ++ipart)
{
hlp.prep_write(uv[2*ipart], uv[2*ipart+1]);
auto ptr = hlp.p0;
for (int cu=0; cu<hlp.w; ++cu)
auto * RESTRICT ptr = hlp.p0;
int w = hlp.w;
auto v = vis[ipart];
const double * RESTRICT ku = hlp.kernel.data();
const double * RESTRICT kv = hlp.kernel.data()+hlp.w;
for (int cu=0; cu<w; ++cu)
{
complex<double> tmp = vis[ipart]*hlp.kernel[cu];
for (int cv=0; cv<hlp.w; ++cv)
*ptr++ += tmp*hlp.kernel[hlp.w+cv];
ptr+=hlp.sv-hlp.w;
complex<double> tmp = v*ku[cu];
for (int cv=0; cv<w; ++cv)
ptr[cv] += tmp*kv[cv];
ptr+=hlp.sv;
}
}
} // end of parallel region
......@@ -467,14 +477,17 @@ a_c_c from_grid (const a_d_c &uv_, const a_c_c &grid_,
{
hlp.prep_read(uv[2*ipart], uv[2*ipart+1]);
complex<double> r = 0.;
auto ptr = hlp.p0;
for (int cu=0; cu<hlp.w; ++cu)
auto * RESTRICT ptr = hlp.p0;
int w = hlp.w;
const double * RESTRICT ku = hlp.kernel.data();
const double * RESTRICT kv = hlp.kernel.data()+hlp.w;
for (int cu=0; cu<w; ++cu)
{
complex<double> tmp = 0.;
for (int cv=0; cv<hlp.w; ++cv)
tmp += (*ptr++) * hlp.kernel[hlp.w+cv];
r+=tmp*hlp.kernel[cu];
ptr += hlp.sv-hlp.w;
for (int cv=0; cv<w; ++cv)
tmp += ptr[cv] * kv[cv];
r += tmp*ku[cu];
ptr += hlp.sv;
}
vis[ipart] = r;
}
......@@ -512,14 +525,13 @@ a_d_c correction_factors (size_t n, size_t nval, double epsilon)
auto w = int(-log10(epsilon)+1.9999);
auto beta = 2.3*w;
auto p = int(1.5*w+2);
double h = 2*pi/n;
double alpha = pi*w/n;
vector<double> x, wgt;
legendre_prep(2*p,x,wgt);
auto psi = x;
for (auto &v:psi)
v = exp(beta*(sqrt(1-v*v)-1.));
int odim[] = {nval};
int odim[] = {int(nval)};
a_d_c res(odim);
auto val = res.mutable_data();
for (size_t k=0; k<nval; ++k)
......
Supports Markdown
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