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