LIBINT  2.6.0
GenericGaussDeriv.impl.h
1 /*
2  * Copyright (C) 2004-2019 Edward F. Valeev
3  *
4  * This file is part of Libint.
5  *
6  * Libint is free software: you can redistribute it and/or modify
7  * it under the terms of the GNU Lesser General Public License as published by
8  * the Free Software Foundation, either version 3 of the License, or
9  * (at your option) any later version.
10  *
11  * Libint is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  * GNU Lesser General Public License for more details.
15  *
16  * You should have received a copy of the GNU Lesser General Public License
17  * along with Libint. If not, see <http://www.gnu.org/licenses/>.
18  *
19  */
20 
21 #ifndef _libint2_src_lib_libint_genericgaussderivimpl_h_
22 #define _libint2_src_lib_libint_genericgaussderivimpl_h_
23 
24 #include <GenericGaussDeriv.h>
25 
26 namespace libint2 {
27 
32  template <int L,
33  bool vectorize>
34  void
35  GenericGaussDeriv<L,vectorize>::compute(const Libint_t* inteval,
36  LIBINT2_REALTYPE* target,
37  const LIBINT2_REALTYPE* src0,
38  const LIBINT2_REALTYPE* src1,
39  unsigned int highdim,
40  unsigned int lowdim,
41  unsigned int dir,
42  const LIBINT2_REALTYPE (&two_alpha)[LIBINT2_MAX_VECLEN]
43  )
44  {
45 
46  const unsigned int veclen = vectorize ? inteval->veclen : 1;
47  const unsigned int lveclen = lowdim * veclen;
48 
49  const unsigned int N = INT_NCART(L);
50  const unsigned int Np1 = INT_NCART(L+1);
51  const unsigned int Nm1 = L > 0 ? INT_NCART(L-1) : 0; // L == 0 case will not use Nm1
52 
53  for(unsigned int h=0, target_idx=0; h<highdim; ++h) {
54 
55  const unsigned int hNp1 = h * Np1;
56  const unsigned int hNm1 = h * Nm1;
57 
58  int ax, ay, az;
59  FOR_CART(ax, ay, az, L)
60 
61  int a[3]; a[0] = ax; a[1] = ay; a[2] = az;
62 
63  // d |a) / dRi = 2 alpha a+1_i
64  ++a[dir];
65  const unsigned int iap1 = INT_CARTINDEX(L+1,a[0],a[1]);
66  const unsigned int ap1_offset = (hNp1 + iap1) * lveclen;
67  const LIBINT2_REALTYPE* src0_ptr = src0 + ap1_offset;
68  --a[dir];
69 
70  const bool have_am1 = (a[dir] > 0);
71  // d |a) / dRi -= a_i a-1_i
72  unsigned int iam1;
73  unsigned int am1_offset;
74  const LIBINT2_REALTYPE* src1_ptr;
75  if (have_am1) {
76  --a[dir];
77  iam1 = INT_CARTINDEX(L-1,a[0],a[1]);
78  am1_offset = (hNm1 + iam1) * lveclen;
79  src1_ptr = src1 + am1_offset;
80  ++a[dir];
81  }
82  LIBINT2_REALTYPE adir_real = LIBINT2_REALTYPE(a[dir]);
83 
84  if (have_am1)
85  for(unsigned int l = 0, lv=0; l < lowdim; ++l) {
86  for(unsigned int v=0; v<veclen; ++v, ++lv, ++target_idx) {
87  target[target_idx] = two_alpha[v] * src0_ptr[lv] - adir_real * src1_ptr[lv];
88  }
89  }
90  else
91  for(unsigned int l = 0, lv=0; l < lowdim; ++l) {
92  for(unsigned int v=0; v<veclen; ++v, ++lv, ++target_idx) {
93  target[target_idx] = two_alpha[v] * src0_ptr[lv];
94  }
95  }
96 
97 #if LIBINT2_FLOP_COUNT
98  inteval->nflops[0] += (have_am1 ? 3 : 1) * lveclen;
99 #endif
100 
101  END_FOR_CART // end of loop over a
102 
103  }
104  }
105 
106 }; // namespace libint2
107 
108 #endif // header guard
109 
static void compute(const Libint_t *inteval, LIBINT2_REALTYPE *target, const LIBINT2_REALTYPE *src0, const LIBINT2_REALTYPE *src1, unsigned int highdim, unsigned int lowdim, unsigned int dir, const LIBINT2_REALTYPE(&two_alpha)[LIBINT2_MAX_VECLEN])
builds ( ...
Definition: GenericGaussDeriv.impl.h:35
Defaults definitions for various parameters assumed by Libint.
Definition: algebra.cc:24