(git:374b731)
Loading...
Searching...
No Matches
grid_hip_prepare_pab.h
Go to the documentation of this file.
1/*----------------------------------------------------------------------------*/
2/* CP2K: A general program to perform molecular dynamics simulations */
3/* Copyright 2000-2024 CP2K developers group <https://cp2k.org> */
4/* */
5/* SPDX-License-Identifier: BSD-3-Clause */
6/*----------------------------------------------------------------------------*/
7
8/*
9 * NOTE : derived from the reference and GPU backends
10 * Authors :
11 - Dr Mathieu Taillefumier (ETH Zurich / CSCS)
12 - Advanced Micro Devices, Inc.
13*/
14
15#include "../common/grid_constants.h"
16#include <stdbool.h>
17#include <stdio.h>
18#include <stdlib.h>
19
20namespace rocm_backend {
21/*******************************************************************************
22 * \brief Implementation of function GRID_FUNC_AB, ie. identity transformation.
23 ******************************************************************************/
24template <typename T>
25__device__ __inline__ void prepare_pab_AB(const orbital a, const orbital b,
26 const T pab_val, const int n,
27 T *cab) {
28
29 // simply copy pab to cab
30 prep_term<T>(a, b, pab_val, n, cab);
31}
32
33/*******************************************************************************
34 * \brief Implementation of function GRID_FUNC_DADB.
35 ******************************************************************************/
36template <typename T>
37__device__ __inline__ void
38prepare_pab_DADB(const orbital a, const orbital b, const T zeta, const T zetb,
39 const T pab_val, const int n, T *cab) {
40
41 // creates cab such that mapping it with pgf_a pgf_b
42 // is equivalent to mapping pab with 0.5 * (nabla pgf_a) . (nabla pgf_b)
43 // (ddx pgf_a ) (ddx pgf_b) = (lax pgf_{a-1x} - 2*zeta*pgf_{a+1x})*(lbx
44 // pgf_{b-1x} - 2*zetb*pgf_{b+1x})
45
46 for (int i = 0; i < 3; i++) {
47 prep_term<T>(down(i, a), down(i, b), 0.5 * a.l[i] * b.l[i] * pab_val, n,
48 cab);
49 prep_term<T>(down(i, a), up(i, b), -1.0 * a.l[i] * zetb * pab_val, n, cab);
50 prep_term<T>(up(i, a), down(i, b), -1.0 * zeta * b.l[i] * pab_val, n, cab);
51 prep_term<T>(up(i, a), up(i, b), 2.0 * zeta * zetb * pab_val, n, cab);
52 }
53}
54
55/*******************************************************************************
56 * \brief Implementation of function GRID_FUNC_ADBmDAB_{X,Y,Z}.
57 ******************************************************************************/
58template <typename T>
59__device__ void prepare_pab_ADBmDAB(const int idir, const orbital a,
60 const orbital b, const T zeta, const T zetb,
61 const T pab_val, const int n, T *cab) {
62
63 // creates cab such that mapping it with pgf_a pgf_b
64 // is equivalent to mapping pab with
65 // pgf_a (nabla_{idir} pgf_b) - (nabla_{idir} pgf_a) pgf_b
66 // ( pgf_a ) (ddx pgf_b) - (ddx pgf_a)( pgf_b ) =
67 // pgf_a *(lbx pgf_{b-1x} - 2*zetb*pgf_{b+1x}) -
68 // (lax pgf_{a-1x} - 2*zeta*pgf_{a+1x}) pgf_b
69
70 prep_term<T>(a, down(idir, b), +b.l[idir] * pab_val, n, cab);
71 prep_term<T>(a, up(idir, b), -2.0 * zetb * pab_val, n, cab);
72 prep_term<T>(down(idir, a), b, -a.l[idir] * pab_val, n, cab);
73 prep_term<T>(up(idir, a), b, +2.0 * zeta * pab_val, n, cab);
74}
75
76/*******************************************************************************
77 * \brief Implementation of function GRID_FUNC_ARDBmDARB_{X,Y,Z}{X,Y,Z}.
78 ******************************************************************************/
79template <typename T>
80__device__ void prepare_pab_ARDBmDARB(const int idir, const int ir,
81 const orbital a, const orbital b,
82 const T zeta, const T zetb,
83 const T pab_val, const int n, T *cab) {
84
85 // creates cab such that mapping it with pgf_a pgf_b
86 // is equivalent to mapping pab with
87 // pgf_a (r-Rb)_{ir} (nabla_{idir} pgf_b) - (nabla_{idir} pgf_a) (r-Rb)_{ir}
88 // pgf_b ( pgf_a )(r-Rb)_{ir} (ddx pgf_b) - (ddx pgf_a) (r-Rb)_{ir} ( pgf_b )
89 // =
90 // pgf_a *(lbx pgf_{b-1x+1ir} - 2*zetb*pgf_{b+1x+1ir})
91 // -
92 // (lax pgf_{a-1x} - 2*zeta*pgf_{a+1x}) pgf_{b+1ir}
93
94 prep_term<T>(a, down(idir, up(ir, b)), b.l[idir] * pab_val, n, cab);
95 prep_term<T>(a, up(idir, up(ir, b)), -2.0 * zetb * pab_val, n, cab);
96 prep_term<T>(down(idir, a), up(ir, b), -a.l[idir] * pab_val, n, cab);
97 prep_term<T>(up(idir, a), up(ir, b), +2.0 * zeta * pab_val, n, cab);
98}
99
100/*******************************************************************************
101 * \brief Implementation of function GRID_FUNC_DABpADB_{X,Y,Z}.
102 * \author Ole Schuett
103 ******************************************************************************/
104template <typename T>
105__device__ void prepare_pab_DABpADB(const int idir, const orbital a,
106 const orbital b, const T zeta, const T zetb,
107 const T pab_val, const int n, T *cab) {
108
109 // creates cab such that mapping it with pgf_a pgf_b
110 // is equivalent to mapping pab with
111 // pgf_a (nabla_{idir} pgf_b) + (nabla_{idir} pgf_a) pgf_b
112 // ( pgf_a ) (ddx pgf_b) + (ddx pgf_a)( pgf_b ) =
113 // pgf_a *(lbx pgf_{b-1x} + 2*zetb*pgf_{b+1x}) +
114 // (lax pgf_{a-1x} + 2*zeta*pgf_{a+1x}) pgf_b
115
116 prep_term<T>(a, down(idir, b), b.l[idir] * pab_val, n, cab);
117 prep_term<T>(a, up(idir, b), -2.0 * zetb * pab_val, n, cab);
118 prep_term<T>(down(idir, a), b, a.l[idir] * pab_val, n, cab);
119 prep_term<T>(up(idir, a), b, -2.0 * zeta * pab_val, n, cab);
120}
121
122/*******************************************************************************
123 * \brief Implementation of function GRID_FUNC_{DX,DY,DZ}.
124 * \author Ole Schuett
125 ******************************************************************************/
126template <typename T>
127__device__ __inline__ void
128prepare_pab_Di(const int ider, const orbital a, const orbital b, const T zeta,
129 const T zetb, const T pab_val, const int n, T *cab) {
130
131 // creates cab such that mapping it with pgf_a pgf_b
132 // is equivalent to mapping pab with
133 // d_{ider1} pgf_a d_{ider1} pgf_b
134 // dx pgf_a dx pgf_b =
135 // (lax pgf_{a-1x})*(lbx pgf_{b-1x}) - 2*zetb*lax*pgf_{a-1x}*pgf_{b+1x}
136 // -
137 // lbx pgf_{b-1x}*2*zeta*pgf_{a+1x}+ 4*zeta*zetab*pgf_{a+1x}pgf_{b+1x}
138
139 prep_term<T>(down(ider, a), down(ider, b), a.l[ider] * b.l[ider] * pab_val, n,
140 cab);
141 prep_term<T>(down(ider, a), up(ider, b), -2.0 * a.l[ider] * zetb * pab_val, n,
142 cab);
143 prep_term<T>(up(ider, a), down(ider, b), -2.0 * zeta * b.l[ider] * pab_val, n,
144 cab);
145 prep_term<T>(up(ider, a), up(ider, b), +4.0 * zeta * zetb * pab_val, n, cab);
146}
147
148/*******************************************************************************
149 * \brief Helper for grid_prepare_pab_DiDj.
150 ******************************************************************************/
151template <typename T>
152__device__ __inline__ void oneterm_dijdij(const int idir, const T func_a,
153 const orbital a, const orbital b,
154 const T zetb, const int n, T *cab) {
155
156 int i1, i2;
157 if (idir == 0) {
158 i1 = 0;
159 i2 = 1;
160 } else if (idir == 1) {
161 i1 = 1;
162 i2 = 2;
163 } else if (idir == 2) {
164 i1 = 2;
165 i2 = 0;
166 } else {
167 return; // error
168 }
169
170 prep_term<T>(a, down(i1, down(i2, b)), b.l[i1] * b.l[i2] * func_a, n, cab);
171 prep_term<T>(a, up(i1, down(i2, b)), -2.0 * zetb * b.l[i2] * func_a, n, cab);
172 prep_term<T>(a, down(i1, up(i2, b)), -2.0 * zetb * b.l[i1] * func_a, n, cab);
173 prep_term<T>(a, up(i1, up(i2, b)), +4.0 * zetb * zetb * func_a, n, cab);
174}
175
176/*******************************************************************************
177 * \brief Implementation of function GRID_FUNC_{DXDY,DYDZ,DZDX}
178 * \author Ole Schuett
179 ******************************************************************************/
180template <typename T>
181__device__ __inline__ void
182prepare_pab_DiDj(const int ider1, const int ider2, const orbital a,
183 const orbital b, const T zeta, const T zetb, const T pab_val,
184 const int n, T *cab) {
185
186 // creates cab such that mapping it with pgf_a pgf_b
187 // is equivalent to mapping pab with
188 // d_{ider1} pgf_a d_{ider1} pgf_b
189
190 const T func_a1 = a.l[ider1] * a.l[ider2] * pab_val;
191 oneterm_dijdij<T>(ider1, func_a1, down(ider1, down(ider2, a)), b, zetb, n,
192 cab);
193
194 const T func_a2 = -2.0 * zeta * a.l[ider2] * pab_val;
195 oneterm_dijdij<T>(ider1, func_a2, up(ider1, down(ider2, a)), b, zetb, n, cab);
196
197 const T func_a3 = -2.0 * zeta * a.l[ider1] * pab_val;
198 oneterm_dijdij<T>(ider1, func_a3, down(ider1, up(ider2, a)), b, zetb, n, cab);
199
200 const T func_a4 = 4.0 * zeta * zeta * pab_val;
201 oneterm_dijdij<T>(ider1, func_a4, up(ider1, up(ider2, a)), b, zetb, n, cab);
202}
203
204/*******************************************************************************
205 * \brief Helper for grid_prepare_pab_Di2.
206 * \author Ole Schuett
207 ******************************************************************************/
208template <typename T>
209__device__ void oneterm_diidii(const int idir, const T func_a, const orbital a,
210 const orbital b, const T zetb, const int n,
211 T *cab) {
212
213 prep_term<T>(a, down(idir, down(idir, b)),
214 b.l[idir] * (b.l[idir] - 1) * func_a, n, cab);
215 prep_term<T>(a, b, -2.0 * zetb * (2 * b.l[idir] + 1) * func_a, n, cab);
216 prep_term<T>(a, up(idir, up(idir, b)), +4.0 * zetb * zetb * func_a, n, cab);
217}
218
219/*******************************************************************************
220 * \brief Implementation of function GRID_FUNC_{DXDX,DYDY,DZDZ}
221 ******************************************************************************/
222template <typename T>
223__device__ __inline__ void
224prepare_pab_Di2(const int ider, const orbital a, const orbital b, const T zeta,
225 const T zetb, const T pab_val, const int n, T *cab) {
226
227 // creates cab such that mapping it with pgf_a pgf_b
228 // is equivalent to mapping pab with
229 // dd_{ider1} pgf_a dd_{ider1} pgf_b
230
231 const T func_a1 = a.l[ider] * (a.l[ider] - 1) * pab_val;
232 oneterm_diidii<T>(ider, func_a1, down(ider, down(ider, a)), b, zetb, n, cab);
233
234 const T func_a2 = -2.0 * zeta * (2 * a.l[ider] + 1) * pab_val;
235 oneterm_diidii<T>(ider, func_a2, a, b, zetb, n, cab);
236
237 const T func_a3 = 4.0 * zeta * zeta * pab_val;
238 oneterm_diidii<T>(ider, func_a3, up(ider, up(ider, a)), b, zetb, n, cab);
239}
240
241/*******************************************************************************
242 * \brief Transforms a given element of the density matrix according to func.
243 * \param func Transformation function to apply, one of GRID_FUNC_*.
244 * \param {a,b} Orbital angular momenta.
245 * \param zet_{a,b} Gaussian exponents.
246 * \param pab_val Input matrix element of pab.
247 * \param n Leading dimensions of output matrix cab.
248 * \param cab Output matrix.
249 ******************************************************************************/
250template <typename T>
251__device__ __inline__ void
252prepare_pab(const enum grid_func func, const orbital a, const orbital b,
253 const T zeta, const T zetb, const T pab_val, const int n, T *cab) {
254
255 // This switch statment will be in an inner loop but only with few iterations.
256 switch (func) {
257 case GRID_FUNC_AB:
258 prepare_pab_AB<T>(a, b, pab_val, n, cab);
259 break;
260 case GRID_FUNC_DADB:
261 prepare_pab_DADB<T>(a, b, zeta, zetb, pab_val, n, cab);
262 break;
264 prepare_pab_ADBmDAB<T>(0, a, b, zeta, zetb, pab_val, n, cab);
265 break;
267 prepare_pab_ADBmDAB<T>(1, a, b, zeta, zetb, pab_val, n, cab);
268 break;
270 prepare_pab_ADBmDAB<T>(2, a, b, zeta, zetb, pab_val, n, cab);
271 break;
273 prepare_pab_ARDBmDARB<T>(0, 0, a, b, zeta, zetb, pab_val, n, cab);
274 break;
276 prepare_pab_ARDBmDARB<T>(0, 1, a, b, zeta, zetb, pab_val, n, cab);
277 break;
279 prepare_pab_ARDBmDARB<T>(0, 2, a, b, zeta, zetb, pab_val, n, cab);
280 break;
282 prepare_pab_ARDBmDARB<T>(1, 0, a, b, zeta, zetb, pab_val, n, cab);
283 break;
285 prepare_pab_ARDBmDARB<T>(1, 1, a, b, zeta, zetb, pab_val, n, cab);
286 break;
288 prepare_pab_ARDBmDARB<T>(1, 2, a, b, zeta, zetb, pab_val, n, cab);
289 break;
291 prepare_pab_ARDBmDARB<T>(2, 0, a, b, zeta, zetb, pab_val, n, cab);
292 break;
294 prepare_pab_ARDBmDARB<T>(2, 1, a, b, zeta, zetb, pab_val, n, cab);
295 break;
297 prepare_pab_ARDBmDARB<T>(2, 2, a, b, zeta, zetb, pab_val, n, cab);
298 break;
300 prepare_pab_DABpADB<T>(0, a, b, zeta, zetb, pab_val, n, cab);
301 break;
303 prepare_pab_DABpADB<T>(1, a, b, zeta, zetb, pab_val, n, cab);
304 break;
306 prepare_pab_DABpADB<T>(2, a, b, zeta, zetb, pab_val, n, cab);
307 break;
308 case GRID_FUNC_DX:
309 prepare_pab_Di<T>(0, a, b, zeta, zetb, pab_val, n, cab);
310 break;
311 case GRID_FUNC_DY:
312 prepare_pab_Di<T>(1, a, b, zeta, zetb, pab_val, n, cab);
313 break;
314 case GRID_FUNC_DZ:
315 prepare_pab_Di<T>(2, a, b, zeta, zetb, pab_val, n, cab);
316 break;
317 case GRID_FUNC_DXDY:
318 prepare_pab_DiDj<T>(0, 1, a, b, zeta, zetb, pab_val, n, cab);
319 break;
320 case GRID_FUNC_DYDZ:
321 prepare_pab_DiDj<T>(1, 2, a, b, zeta, zetb, pab_val, n, cab);
322 break;
323 case GRID_FUNC_DZDX:
324 prepare_pab_DiDj<T>(2, 0, a, b, zeta, zetb, pab_val, n, cab);
325 break;
326 case GRID_FUNC_DXDX:
327 prepare_pab_Di2<T>(0, a, b, zeta, zetb, pab_val, n, cab);
328 break;
329 case GRID_FUNC_DYDY:
330 prepare_pab_Di2<T>(1, a, b, zeta, zetb, pab_val, n, cab);
331 break;
332 case GRID_FUNC_DZDZ:
333 prepare_pab_Di2<T>(2, a, b, zeta, zetb, pab_val, n, cab);
334 break;
335 default:
336 break; // Error: Unknown ga_gb_function - do nothing.
337 }
338}
339
340/*******************************************************************************
341 * \brief Returns difference in angular momentum range for given func.
342 * \author Ole Schuett
343 ******************************************************************************/
345 ldiffs_value ldiffs;
346
347 switch (func) {
348 case GRID_FUNC_AB:
349 ldiffs.la_max_diff = 0;
350 ldiffs.la_min_diff = 0;
351 ldiffs.lb_max_diff = 0;
352 ldiffs.lb_min_diff = 0;
353 break;
354 case GRID_FUNC_DADB:
361 ldiffs.la_max_diff = +1;
362 ldiffs.la_min_diff = -1;
363 ldiffs.lb_max_diff = +1;
364 ldiffs.lb_min_diff = -1;
365 break;
375 ldiffs.la_max_diff = +1;
376 ldiffs.la_min_diff = -1;
377 ldiffs.lb_max_diff = +2; // this is legit
378 ldiffs.lb_min_diff = -1;
379 break;
380 case GRID_FUNC_DX:
381 case GRID_FUNC_DY:
382 case GRID_FUNC_DZ:
383 ldiffs.la_max_diff = +1;
384 ldiffs.la_min_diff = -1;
385 ldiffs.lb_max_diff = +1;
386 ldiffs.lb_min_diff = -1;
387 break;
388 case GRID_FUNC_DXDY:
389 case GRID_FUNC_DYDZ:
390 case GRID_FUNC_DZDX:
391 case GRID_FUNC_DXDX:
392 case GRID_FUNC_DYDY:
393 case GRID_FUNC_DZDZ:
394 ldiffs.la_max_diff = +2;
395 ldiffs.la_min_diff = -2;
396 ldiffs.lb_max_diff = +2;
397 ldiffs.lb_min_diff = -2;
398 break;
399 default:
400 fprintf(stderr, "Error: Unknown ga_gb_function %i.\n", func);
401 abort();
402 }
403
404 return ldiffs;
405}
406} // namespace rocm_backend
grid_func
@ GRID_FUNC_ADBmDAB_X
@ GRID_FUNC_DZ
@ GRID_FUNC_DABpADB_X
@ GRID_FUNC_DY
@ GRID_FUNC_ARDBmDARB_YY
@ GRID_FUNC_DX
@ GRID_FUNC_ARDBmDARB_XY
@ GRID_FUNC_ARDBmDARB_ZY
@ GRID_FUNC_DYDY
@ GRID_FUNC_ARDBmDARB_XZ
@ GRID_FUNC_DZDX
@ GRID_FUNC_ARDBmDARB_YZ
@ GRID_FUNC_AB
@ GRID_FUNC_DYDZ
@ GRID_FUNC_ADBmDAB_Z
@ GRID_FUNC_DXDY
@ GRID_FUNC_DXDX
@ GRID_FUNC_ADBmDAB_Y
@ GRID_FUNC_DABpADB_Z
@ GRID_FUNC_ARDBmDARB_ZX
@ GRID_FUNC_DZDZ
@ GRID_FUNC_ARDBmDARB_YX
@ GRID_FUNC_ARDBmDARB_ZZ
@ GRID_FUNC_ARDBmDARB_XX
@ GRID_FUNC_DABpADB_Y
@ GRID_FUNC_DADB
static void const int const int i
__device__ void prepare_pab_ADBmDAB(const int idir, const orbital a, const orbital b, const T zeta, const T zetb, const T pab_val, const int n, T *cab)
Implementation of function GRID_FUNC_ADBmDAB_{X,Y,Z}.
__device__ __inline__ void oneterm_dijdij(const int idir, const T func_a, const orbital a, const orbital b, const T zetb, const int n, T *cab)
Helper for grid_prepare_pab_DiDj.
__device__ __inline__ void prepare_pab_Di2(const int ider, const orbital a, const orbital b, const T zeta, const T zetb, const T pab_val, const int n, T *cab)
Implementation of function GRID_FUNC_{DXDX,DYDY,DZDZ}.
__device__ __inline__ void prepare_pab_DiDj(const int ider1, const int ider2, const orbital a, const orbital b, const T zeta, const T zetb, const T pab_val, const int n, T *cab)
Implementation of function GRID_FUNC_{DXDY,DYDZ,DZDX}.
__device__ __inline__ void prepare_pab_AB(const orbital a, const orbital b, const T pab_val, const int n, T *cab)
Implementation of function GRID_FUNC_AB, ie. identity transformation.
ldiffs_value prepare_get_ldiffs(const enum grid_func func)
Returns difference in angular momentum range for given func.
__device__ void prepare_pab_DABpADB(const int idir, const orbital a, const orbital b, const T zeta, const T zetb, const T pab_val, const int n, T *cab)
Implementation of function GRID_FUNC_DABpADB_{X,Y,Z}.
__device__ __inline__ void prepare_pab_Di(const int ider, const orbital a, const orbital b, const T zeta, const T zetb, const T pab_val, const int n, T *cab)
Implementation of function GRID_FUNC_{DX,DY,DZ}.
__device__ void prepare_pab_ARDBmDARB(const int idir, const int ir, const orbital a, const orbital b, const T zeta, const T zetb, const T pab_val, const int n, T *cab)
Implementation of function GRID_FUNC_ARDBmDARB_{X,Y,Z}{X,Y,Z}.
__inline__ __device__ orbital down(const int i, const orbital &a)
Decrease i'th component of given orbital angular momentum.
__device__ __inline__ orbital up(const int i, const orbital &a)
Increase i'th component of given orbital angular momentum.
__device__ void oneterm_diidii(const int idir, const T func_a, const orbital a, const orbital b, const T zetb, const int n, T *cab)
Helper for grid_prepare_pab_Di2.
__device__ __inline__ void prepare_pab(const enum grid_func func, const orbital a, const orbital b, const T zeta, const T zetb, const T pab_val, const int n, T *cab)
Transforms a given element of the density matrix according to func.
__device__ __inline__ void prepare_pab_DADB(const orbital a, const orbital b, const T zeta, const T zetb, const T pab_val, const int n, T *cab)
Implementation of function GRID_FUNC_DADB.
Differences in angular momentum.
Orbital angular momentum.