Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Kokkos-free einspline reference implementation #1

Open
wants to merge 2 commits into
base: kokkos
Choose a base branch
from
Open

Kokkos-free einspline reference implementation #1

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
334b592
Initial commit for Kokkos-free einspline_spo reference implementation.
jptowns Apr 25, 2019
0b4b0a7
check_spo now passes on both CPU and GPU!
jptowns Apr 25, 2019
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions src/CMakeLists.txt
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -137,7 +137,9 @@ ENDIF()
Utilities/XMLWriter.cpp
Utilities/tinyxml/tinyxml2.cpp
Utilities/qmcpack_version.cpp
Numerics/Spline2/einspline_allocator_ref.cpp
Numerics/Spline2/einspline_allocator.cpp
Numerics/Spline2/bspline_allocator_ref.cpp
Numerics/Spline2/bspline_allocator.cpp
Numerics/Spline2/MultiBsplineData.cpp
${GITREV_TMP}
Expand Down
56 changes: 32 additions & 24 deletions src/Drivers/check_spo.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -57,7 +57,6 @@ int main(int argc, char** argv)
Kokkos::initialize(argc, argv);
{ //Begin kokkos block.


// clang-format off
typedef QMCTraits::RealType RealType;
typedef ParticleSet::ParticlePos_t ParticlePos_t;
Expand Down Expand Up @@ -188,8 +187,10 @@ int main(int argc, char** argv)
app_summary() << "\nSPO coefficients size = " << SPO_coeff_size << " bytes ("
<< SPO_coeff_size_MB << " MB)" << endl;

std::cout << "Setting up the Kokkos implementation...";
spo_main.set(nx, ny, nz, norb);
spo_main.Lattice.set(lattice_b);
std::cout << "Done." << std::endl;

/*
auto& spline = spo_main.spline;
Expand Down Expand Up @@ -219,13 +220,20 @@ int main(int argc, char** argv)
cout << " psi[0] = " << psiMirror(0) << ", psi[4] = " << psiMirror(4) << endl;
*/

// REF
std::cout << "Setting up the reference implementation...";
spo_ref_main.set(nx, ny, nz, norb, nTiles);
spo_ref_main.Lattice.set(lattice_b);
/*
std::cout << "Done." << std::endl;

/*
cout << " looking at the reference version:" << endl;
cout << " coef(2,2,0,0) = " << spo_ref_main.einsplines(0).coefs_view(2,2,0,0) << endl;
cout << " coef(2,0,2,0) = " << spo_ref_main.einsplines(0).coefs_view(2,0,2,0) << endl;
cout << " coef(0,2,2,0) = " << spo_ref_main.einsplines(0).coefs_view(0,2,2,0) << endl;
intptr_t xs = spo_ref_main.einsplines(0).x_stride;
intptr_t ys = spo_ref_main.einsplines(0).y_stride;
intptr_t zs = spo_ref_main.einsplines(0).z_stride;
cout << " coef(2,2,0,0) = " << spo_ref_main.einsplines(0).coefs[xs*2+ys*2+0*zs] << endl;
cout << " coef(2,0,2,0) = " << spo_ref_main.einsplines(0).coefs[xs*2+ys*0+zs*2] << endl;
cout << " coef(0,2,2,0) = " << spo_ref_main.einsplines(0).coefs[xs*0+ys*2+zs*2] << endl;
spo_ref_main.evaluate_v(p);
cout << " doing evaluate v with position: (" << p[0] << ", " << p[1] << ", " << p[2] << ")" << endl;
cout << " psi[0] = " << spo_ref_main.psi(0).operator()(0) << ", psi[4] = " << spo_ref_main.psi(0).operator()(4) << endl;
Expand Down Expand Up @@ -300,30 +308,30 @@ int main(int argc, char** argv)
spo_ref.evaluate_vgh(pos);
// accumulate error

auto psiView = Kokkos::create_mirror_view(spo_main.psi);
auto gradView = Kokkos::create_mirror_view(spo_main.grad);
auto hessView = Kokkos::create_mirror_view(spo_main.hess);
Kokkos::deep_copy(psiView, spo_main.psi);
Kokkos::deep_copy(gradView, spo_main.grad);
Kokkos::deep_copy(hessView, spo_main.hess);
auto psi_host = Kokkos::create_mirror_view(spo_main.psi);
auto grad_host = Kokkos::create_mirror_view(spo_main.grad);
auto hess_host = Kokkos::create_mirror_view(spo_main.hess);
Kokkos::deep_copy(psi_host, spo_main.psi);
Kokkos::deep_copy(grad_host, spo_main.grad);
Kokkos::deep_copy(hess_host, spo_main.hess);

for (int ib = 0; ib < spo_ref.nBlocks; ib++)
for (int n = 0; n < spo_ref.nSplinesPerBlock; n++)
{
int spNum = ib*spo_ref.nSplinesPerBlock+n;
// value
evalVGH_v_err += std::fabs(psiView(spNum) - spo_ref.psi[ib][n]);
evalVGH_v_err += std::fabs(psi_host(spNum) - spo_ref.psi[ib][n]);
// grad
evalVGH_g_err += std::fabs(gradView(spNum, 0) - spo_ref.grad[ib](n, 0));
evalVGH_g_err += std::fabs(gradView(spNum, 1) - spo_ref.grad[ib](n, 1));
evalVGH_g_err += std::fabs(gradView(spNum, 2) - spo_ref.grad[ib](n, 2));
evalVGH_g_err += std::fabs(grad_host(spNum, 0) - spo_ref.grad[ib].data(0)[n]);
evalVGH_g_err += std::fabs(grad_host(spNum, 1) - spo_ref.grad[ib].data(1)[n]);
evalVGH_g_err += std::fabs(grad_host(spNum, 2) - spo_ref.grad[ib].data(2)[n]);
// hess
evalVGH_h_err += std::fabs(hessView(spNum, 0) - spo_ref.hess[ib](n, 0));
evalVGH_h_err += std::fabs(hessView(spNum, 1) - spo_ref.hess[ib](n, 1));
evalVGH_h_err += std::fabs(hessView(spNum, 2) - spo_ref.hess[ib](n, 2));
evalVGH_h_err += std::fabs(hessView(spNum, 3) - spo_ref.hess[ib](n, 3));
evalVGH_h_err += std::fabs(hessView(spNum, 4) - spo_ref.hess[ib](n, 4));
evalVGH_h_err += std::fabs(hessView(spNum, 5) - spo_ref.hess[ib](n, 5));
evalVGH_h_err += std::fabs(hess_host(spNum, 0) - spo_ref.hess[ib].data(0)[n]);
evalVGH_h_err += std::fabs(hess_host(spNum, 1) - spo_ref.hess[ib].data(1)[n]);
evalVGH_h_err += std::fabs(hess_host(spNum, 2) - spo_ref.hess[ib].data(2)[n]);
evalVGH_h_err += std::fabs(hess_host(spNum, 3) - spo_ref.hess[ib].data(3)[n]);
evalVGH_h_err += std::fabs(hess_host(spNum, 4) - spo_ref.hess[ib].data(4)[n]);
evalVGH_h_err += std::fabs(hess_host(spNum, 5) - spo_ref.hess[ib].data(5)[n]);
}
if (ur[iel] > accept)
{
Expand All @@ -350,12 +358,12 @@ int main(int argc, char** argv)
//spo.evaluate_v(pos);
spo_ref.evaluate_v(pos);
// accumulate error
auto psiView = Kokkos::create_mirror_view(spo_main.psi);
Kokkos::deep_copy(psiView, spo_main.psi);
auto psi_host = Kokkos::create_mirror_view(spo_main.psi);
Kokkos::deep_copy(psi_host, spo_main.psi);

for (int ib = 0; ib < spo_ref.nBlocks; ib++)
for (int n = 0; n < spo_ref.nSplinesPerBlock; n++)
evalV_v_err += std::fabs(psiView(ib*spo_ref.nSplinesPerBlock+n) - spo_ref.psi[ib][n]);
evalV_v_err += std::fabs(psi_host(ib*spo_ref.nSplinesPerBlock+n) - spo_ref.psi[ib][n]);
}
} // els
} // ions
Expand Down
5 changes: 4 additions & 1 deletion src/Drivers/miniqmc.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -403,13 +403,16 @@ int main(int argc, char** argv)
const int teamID = partition_id;
Timers[Timer_Init]->start();
// create and initialize movers
std::cerr << "Making a mover...\n";
Mover thiswalker(myPrimes[teamID], ions);
// create a spo view in each Mover
std::cerr << "Making the spo for this mover...\n";
thiswalker.spo = build_SPOSet_view(useRef, spo_main, team_size, teamID);

// create wavefunction per mover
std::cerr << "Making the wave function for this mover...\n";
build_WaveFunction(useRef, thiswalker.wavefunction, ions, thiswalker.els, thiswalker.rng, enableJ3);

std::cerr << "Finished making stuff for this mover...\n";
// initial computing
thiswalker.els.update();
thiswalker.wavefunction.evaluateLog(thiswalker.els);
Expand Down
33 changes: 33 additions & 0 deletions src/Numerics/Einspline/bspline_ref.h
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,33 @@
/////////////////////////////////////////////////////////////////////////////
// einspline: a library for creating and evaluating B-splines //
// Copyright (C) 2007 Kenneth P. Esler, Jr. //
// //
// This program is free software; you can redistribute it and/or modify //
// it under the terms of the GNU General Public License as published by //
// the Free Software Foundation; either version 2 of the License, or //
// (at your option) any later version. //
// //
// This program is distributed in the hope that it will be useful, //
// but WITHOUT ANY WARRANTY; without even the implied warranty of //
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the //
// GNU General Public License for more details. //
// //
// You should have received a copy of the GNU General Public License //
// along with this program; if not, write to the Free Software //
// Foundation, Inc., 51 Franklin Street, Fifth Floor, //
// Boston, MA 02110-1301 USA //
/////////////////////////////////////////////////////////////////////////////

#ifndef BSPLINE_H
#define BSPLINE_H

#include "bspline_base.h"
////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////
//// Bspline structure definitions ////
////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////
#include "bspline_structs.h"
#include "multi_bspline_structs_ref.h"

#endif
91 changes: 91 additions & 0 deletions src/Numerics/Einspline/multi_bspline_structs_ref.h
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,91 @@
/////////////////////////////////////////////////////////////////////////////
// einspline: a library for creating and evaluating B-splines //
// Copyright (C) 2007 Kenneth P. Esler, Jr. //
// //
// This program is free software; you can redistribute it and/or modify //
// it under the terms of the GNU General Public License as published by //
// the Free Software Foundation; either version 2 of the License, or //
// (at your option) any later version. //
// //
// This program is distributed in the hope that it will be useful, //
// but WITHOUT ANY WARRANTY; without even the implied warranty of //
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the //
// GNU General Public License for more details. //
// //
// You should have received a copy of the GNU General Public License //
// along with this program; if not, write to the Free Software //
// Foundation, Inc., 51 Franklin Street, Fifth Floor, //
// Boston, MA 02110-1301 USA //
/////////////////////////////////////////////////////////////////////////////

#ifndef MULTI_BSPLINE_STRUCTS_STD_H
#define MULTI_BSPLINE_STRUCTS_STD_H

#include <inttypes.h>
#include <stdlib.h>

///////////////////////////
// Single precision real //
///////////////////////////
struct multi_UBspline_3d_s
{
spline_code spcode;
type_code tcode;
float* restrict coefs;
intptr_t x_stride, y_stride, z_stride;
Ugrid x_grid, y_grid, z_grid;
BCtype_s xBC, yBC, zBC;
int num_splines;
size_t coefs_size;
};

///////////////////////////
// Double precision real //
///////////////////////////
struct multi_UBspline_3d_d
{
spline_code spcode;
type_code tcode;
double* restrict coefs;
intptr_t x_stride, y_stride, z_stride;
Ugrid x_grid, y_grid, z_grid;
BCtype_d xBC, yBC, zBC;
int num_splines;
size_t coefs_size;
};

//////////////////////////////
// Single precision complex //
//////////////////////////////
struct multi_UBspline_3d_c
{
spline_code spcode;
type_code tcode;
complex_float* restrict coefs;
intptr_t x_stride, y_stride, z_stride;
Ugrid x_grid, y_grid, z_grid;
BCtype_c xBC, yBC, zBC;
int num_splines;
size_t coefs_size;
// temporary storage for laplacian components
complex_float* restrict lapl3;
};

//////////////////////////////
// Double precision complex //
//////////////////////////////
struct multi_UBspline_3d_z
{
spline_code spcode;
type_code tcode;
complex_double* restrict coefs;
intptr_t x_stride, y_stride, z_stride;
Ugrid x_grid, y_grid, z_grid;
BCtype_z xBC, yBC, zBC;
int num_splines;
size_t coefs_size;
// temporary storage for laplacian components
complex_double* restrict lapl3;
};

#endif
48 changes: 20 additions & 28 deletions src/Numerics/Spline2/MultiBsplineRef.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -36,7 +36,7 @@ struct MultiBsplineRef
using spliner_type = typename bspline_traits<T, 3>::SplineType;

MultiBsplineRef() {}
MultiBsplineRef(const MultiBsplineRef& in) = default;
MultiBsplineRef(const MultiBsplineRef& in) = delete;
MultiBsplineRef& operator=(const MultiBsplineRef& in) = delete;

/** compute values vals[0,num_splines)
Expand All @@ -45,12 +45,8 @@ struct MultiBsplineRef
* evaluate_vgh(r,psi,grad,hess,ip).
*/

void evaluate_v(const spliner_type* restrict spline_m,
T x,
T y,
T z,
T* restrict vals,
size_t num_splines) const;
void evaluate_v(
const spliner_type* restrict spline_m, T x, T y, T z, T* restrict vals, size_t num_splines) const;

void evaluate_vgl(const spliner_type* restrict spline_m,
T x,
Expand All @@ -72,12 +68,8 @@ struct MultiBsplineRef
};

template<typename T>
inline void MultiBsplineRef<T>::evaluate_v(const spliner_type* restrict spline_m,
T x,
T y,
T z,
T* restrict vals,
size_t num_splines) const
inline void MultiBsplineRef<T>::evaluate_v(
const spliner_type* restrict spline_m, T x, T y, T z, T* restrict vals, size_t num_splines) const
{
x -= spline_m->x_grid.start;
y -= spline_m->y_grid.start;
Expand Down Expand Up @@ -183,12 +175,12 @@ inline void MultiBsplineRef<T>::evaluate_vgl(const spliner_type* restrict spline
T sum0 = c[0] * coefsv + c[1] * coefsvzs + c[2] * coefsv2zs + c[3] * coefsv3zs;
T sum1 = dc[0] * coefsv + dc[1] * coefsvzs + dc[2] * coefsv2zs + dc[3] * coefsv3zs;
T sum2 = d2c[0] * coefsv + d2c[1] * coefsvzs + d2c[2] * coefsv2zs + d2c[3] * coefsv3zs;
gx[n] += pre10 * sum0;
gy[n] += pre01 * sum0;
gz[n] += pre00 * sum1;
lx[n] += pre20 * sum0;
ly[n] += pre02 * sum0;
lz[n] += pre00 * sum2;
gx[n] += pre10 * sum0;
gy[n] += pre01 * sum0;
gz[n] += pre00 * sum1;
lx[n] += pre20 * sum0;
ly[n] += pre02 * sum0;
lz[n] += pre00 * sum2;
vals[n] += pre00 * sum0;
}
}
Expand Down Expand Up @@ -290,15 +282,15 @@ inline void MultiBsplineRef<T>::evaluate_vgh(const spliner_type* restrict spline
T sum1 = dc[0] * coefsv + dc[1] * coefsvzs + dc[2] * coefsv2zs + dc[3] * coefsv3zs;
T sum2 = d2c[0] * coefsv + d2c[1] * coefsvzs + d2c[2] * coefsv2zs + d2c[3] * coefsv3zs;

hxx[n] += pre20 * sum0;
hxy[n] += pre11 * sum0;
hxz[n] += pre10 * sum1;
hyy[n] += pre02 * sum0;
hyz[n] += pre01 * sum1;
hzz[n] += pre00 * sum2;
gx[n] += pre10 * sum0;
gy[n] += pre01 * sum0;
gz[n] += pre00 * sum1;
hxx[n] += pre20 * sum0;
hxy[n] += pre11 * sum0;
hxz[n] += pre10 * sum1;
hyy[n] += pre02 * sum0;
hyz[n] += pre01 * sum1;
hzz[n] += pre00 * sum2;
gx[n] += pre10 * sum0;
gy[n] += pre01 * sum0;
gz[n] += pre00 * sum1;
vals[n] += pre00 * sum0;
}
}
Expand Down
Loading