From 2ca7dcb853acadba95b77f0173ba6e40a38c077f Mon Sep 17 00:00:00 2001 From: Stan Gerald Moore Date: Wed, 27 Oct 2021 13:41:29 -0600 Subject: [PATCH] Use lighter construct than powint --- src/KOKKOS/compute_orientorder_atom_kokkos.cpp | 14 ++++++++------ src/compute_orientorder_atom.cpp | 10 ++++++---- 2 files changed, 14 insertions(+), 10 deletions(-) diff --git a/src/KOKKOS/compute_orientorder_atom_kokkos.cpp b/src/KOKKOS/compute_orientorder_atom_kokkos.cpp index 548c503628..42d18a3df8 100644 --- a/src/KOKKOS/compute_orientorder_atom_kokkos.cpp +++ b/src/KOKKOS/compute_orientorder_atom_kokkos.cpp @@ -26,7 +26,7 @@ #include "atom_masks.h" #include "kokkos.h" #include "math_const.h" -#include "math_special_kokkos.h" +#include "math_special.h" #include "memory_kokkos.h" #include "neigh_list.h" #include "neigh_request.h" @@ -38,7 +38,7 @@ using namespace LAMMPS_NS; using namespace MathConst; -using namespace MathSpecialKokkos; +using namespace MathSpecial; using namespace std; #ifdef DBL_EPSILON @@ -546,6 +546,7 @@ void ComputeOrientOrderAtomKokkos::calc_boop2(int ncount, int ii) co for (int il = 0; il < nqlist; il++) { int l = d_qlist[il]; double wlsum = 0.0; + int sgn = 1; for (int m1 = -l; m1 <= 0; m1++) { for (int m2 = 0; m2 <= ((-m1)>>1); m2++) { const int m3 = -(m1 + m2); @@ -556,8 +557,6 @@ void ComputeOrientOrderAtomKokkos::calc_boop2(int ncount, int ii) co // structure enforces visiting only one member of each // such symmetry (invariance) group. - //const int sgn = 1 - 2*(m1&1); - const int sgn = powint(-1,m1); // sgn = (-1)^m // m1 <= 0, and Qlm[-m] = (-1)^m*conjg(Qlm[m]) SNAcomplex Q1Q2; Q1Q2.re = (d_qnm(ii,il,-m1).re*d_qnm(ii,il,m2).re + d_qnm(ii,il,-m1).im*d_qnm(ii,il,m2).im)*sgn; @@ -565,6 +564,8 @@ void ComputeOrientOrderAtomKokkos::calc_boop2(int ncount, int ii) co const double Q1Q2Q3 = Q1Q2.re*d_qnm(ii,il,m3).re - Q1Q2.im*d_qnm(ii,il,m3).im; const double c = d_w3jlist[widx_count++]; wlsum += Q1Q2Q3*c; + + sgn = -sgn; // sgn = (-1)^m } } d_qnarray(i,jj++) = wlsum/d_qnormfac2(il); @@ -599,13 +600,14 @@ void ComputeOrientOrderAtomKokkos::calc_boop2(int ncount, int ii) co } else { const double qnfac = d_qnormfac(il)/d_qnarray(i,il); + int sgn = 1; for (int m = -l; m < 0; m++) { // Computed only qnm for m>=0. // qnm[-m] = (-1)^m * conjg(qnm[m]) - //const int sgn = 1 - 2*(m&1); // sgn = (-1)^m - const int sgn = powint(-1,m); // sgn = (-1)^m d_qnarray(i,jj++) = d_qnm(ii,il,-m).re * qnfac * sgn; d_qnarray(i,jj++) = -d_qnm(ii,il,-m).im * qnfac * sgn; + + sgn = -sgn; // sgn = (-1)^m } for (int m = 0; m < l+1; m++) { d_qnarray(i,jj++) = d_qnm(ii,il,m).re * qnfac; diff --git a/src/compute_orientorder_atom.cpp b/src/compute_orientorder_atom.cpp index 6f4f86504d..1bc11ff259 100644 --- a/src/compute_orientorder_atom.cpp +++ b/src/compute_orientorder_atom.cpp @@ -524,6 +524,7 @@ void ComputeOrientOrderAtom::calc_boop(double **rlist, for (int il = 0; il < nqlist; il++) { int l = qlist[il]; double wlsum = 0.0; + int sgn = 1; for (int m1 = -l; m1 <= 0; m1++) { for (int m2 = 0; m2 <= ((-m1)>>1); m2++) { const int m3 = -(m1 + m2); @@ -534,14 +535,14 @@ void ComputeOrientOrderAtom::calc_boop(double **rlist, // structure enforces visiting only one member of each // such symmetry (invariance) group. - //const int sgn = 1 - 2*(m1&1); - const int sgn = powint(-1,m1); // sgn = (-1)^m // m1 <= 0, and Qlm[-m] = (-1)^m*conjg(Qlm[m]) const double Q1Q2_r = (qnm_r[il][-m1]*qnm_r[il][m2] + qnm_i[il][-m1]*qnm_i[il][m2])*sgn; const double Q1Q2_i = (qnm_r[il][-m1]*qnm_i[il][m2] - qnm_i[il][-m1]*qnm_r[il][m2])*sgn; const double Q1Q2Q3 = Q1Q2_r*qnm_r[il][m3] - Q1Q2_i*qnm_i[il][m3]; const double c = w3jlist[widx_count++]; wlsum += Q1Q2Q3*c; + + sgn = -sgn; // sgn = (-1)^m } } qn[jj++] = wlsum/qnormfac2[il]; @@ -576,13 +577,14 @@ void ComputeOrientOrderAtom::calc_boop(double **rlist, } else { double qnfac = qnormfac[il]/qn[il]; + int sgn = 1; for (int m = -l; m < 0; m++) { // Computed only qnm for m>=0. // qnm[-m] = (-1)^m * conjg(qnm[m]) - //const int sgn = 1 - 2*(m&1); // sgn = (-1)^m - const int sgn = powint(-1,m); // sgn = (-1)^m qn[jj++] = qnm_r[il][-m] * qnfac * sgn; qn[jj++] = -qnm_i[il][-m] * qnfac * sgn; + + sgn = -sgn; // sgn = (-1)^m } for (int m = 0; m < l+1; m++) { qn[jj++] = qnm_r[il][m] * qnfac;