1!!  gen1int: compute one-electron integrals using rotational London atomic-orbitals
2!!  Copyright 2009-2012 Bin Gao, and Andreas Thorvaldsen
3!!
4!!  gen1int is free software: you can redistribute it and/or modify
5!!  it under the terms of the GNU Lesser General Public License as published by
6!!  the Free Software Foundation, either version 3 of the License, or
7!!  (at your option) any later version.
8!!
9!!  gen1int 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 Lesser General Public License for more details.
13!!
14!!  You should have received a copy of the GNU Lesser General Public License
15!!  along with gen1int. If not, see <http://www.gnu.org/licenses/>.
16!!
17!!  This file reorders the integrals. You may not need it if your basis
18!!  functions are in the same sequence as those in Gen1Int.
19!!
20!!  2011-08-03, Bin Gao
21!!  * revised for contracted GTOs, also implemented reordering GTOs either on bra or ket center
22!!
23!!  2010-10-07, Bin Gao
24!!  * first version
25
26#include "stdout.h"
27
28  !> \brief reorders the contracted real solid-harmonic Gaussians on bra or ket center
29  !> \author Bin Gao
30  !> \date 2011-08-03
31  !> \param ang_ket is the orbital quantum number (or angular number) of ket center
32  !> \param num_sgto_ket is the number of basis functions on ket center,
33  !>        (equals to \f$2\var(ang_ket)+1\f$)
34  !> \param mag_ket contains the magnetic numbers of basis functions on ket center
35  !> \param dim_sgto_bra is the dimension of SGTOs on bra center
36  !> \param num_contr_ket is the number of contractions of ket center
37  !> \param num_opt is the number of operators
38  !> \param gen_ints contains the contracted integrals from Gen1Int
39  !> \return ro_ints contains the reordered integrals according to given \var(mag_ket)
40  subroutine reorder_sgtos(ang_ket, num_sgto_ket, mag_ket, dim_sgto_bra, &
41                           num_contr_ket, num_opt, gen_ints, ro_ints)
42    use xkind
43    implicit none
44    integer, intent(in) :: ang_ket
45    integer, intent(in) :: num_sgto_ket
46    integer, intent(in) :: mag_ket(num_sgto_ket)
47    integer, intent(in) :: dim_sgto_bra
48    integer, intent(in) :: num_contr_ket
49    integer, intent(in) :: num_opt
50    real(REALK), intent(in) :: gen_ints(dim_sgto_bra,num_sgto_ket,num_contr_ket,num_opt)
51    real(REALK), intent(out) :: ro_ints(dim_sgto_bra,num_sgto_ket,num_contr_ket,num_opt)
52!f2py intent(in) :: ang_ket
53!f2py intent(hide) :: num_sgto_ket
54!f2py intent(in) :: mag_ket
55!f2py intent(hide) :: dim_sgto_bra
56!f2py intent(hide) :: num_contr_ket
57!f2py intent(hide) :: num_opt
58!f2py intent(in) :: gen_ints
59!f2py depend(num_sgto_ket) :: gen_ints
60!f2py intent(out) :: ro_ints
61!f2py depend(dim_sgto_bra) :: ro_ints
62!f2py depend(num_sgto_ket) :: ro_ints
63!f2py depend(num_contr_ket) :: ro_ints
64!f2py depend(num_opt) :: ro_ints
65    integer ang_ket_one               !orbital quantum number of ket +1
66    integer ibra, iket, jcontr, iopt  !incremental recorders
67    integer addr_ket                  !address of SGTOs on ket center in Gen1Int
68#if defined(XTIME)
69    real(REALK) curr_time             !current CPU time
70    ! sets current CPU time
71    call xtimer_set(curr_time)
72#endif
73    ang_ket_one = ang_ket+1
74    ! loops over operators
75    do iopt = 1, num_opt
76      ! loops over contractions on ket center
77      do jcontr = 1, num_contr_ket
78        ! loops over SGTOs on ket center
79        do iket = 1, num_sgto_ket
80          addr_ket = mag_ket(iket)+ang_ket_one
81          ! loops over SGTOs on bra center
82          do ibra = 1, dim_sgto_bra
83            ro_ints(ibra,iket,jcontr,iopt) = gen_ints(ibra,addr_ket,jcontr,iopt)
84          end do
85        end do
86      end do
87    end do
88#if defined(XTIME)
89    ! prints the CPU elapsed time
90    call xtimer_view(curr_time, "reorder_sgtos", STDOUT)
91#endif
92    return
93  end subroutine reorder_sgtos
94
95  !> \brief reorders the integrals of contracted real solid-harmonic Gaussians
96  !> \author Bin Gao
97  !> \date 2010-10-07
98  !> \param ang_bra is the orbital quantum number (or angular number) of bra center
99  !> \param num_sgto_bra is the number of basis functions on bra center,
100  !>        (equals to \f$2\var(ang_bra)+1\f$)
101  !> \param mag_bra contains the magnetic numbers of basis functions on bra center
102  !> \param ang_ket is the orbital quantum number (or angular number) of ket center
103  !> \param num_sgto_ket is the number of basis functions on ket center,
104  !>        (equals to \f$2\var(ang_ket)+1\f$)
105  !> \param mag_ket contains the magnetic numbers of basis functions on ket center
106  !> \param num_contr_bra is the number of contractions of bra center
107  !> \param num_contr_ket is the number of contractions of ket center
108  !> \param num_opt is the number of operators
109  !> \param gen_ints contains the contracted integrals from Gen1Int
110  !> \return ro_ints contains the reordered integrals according to given \var(mag_bra)
111  !>         and \var(mag_ket)
112  subroutine reorder_sgto_ints(ang_bra, num_sgto_bra, mag_bra, &
113                               ang_ket, num_sgto_ket, mag_ket, &
114                               num_contr_bra, num_contr_ket,   &
115                               num_opt, gen_ints, ro_ints)
116    use xkind
117    implicit none
118    integer, intent(in) :: ang_bra
119    integer, intent(in) :: num_sgto_bra
120    integer, intent(in) :: mag_bra(num_sgto_bra)
121    integer, intent(in) :: ang_ket
122    integer, intent(in) :: num_sgto_ket
123    integer, intent(in) :: mag_ket(num_sgto_ket)
124    integer, intent(in) :: num_contr_bra
125    integer, intent(in) :: num_contr_ket
126    integer, intent(in) :: num_opt
127    real(REALK), intent(in) :: gen_ints(num_sgto_bra,num_contr_bra, &
128                                        num_sgto_ket,num_contr_ket,num_opt)
129    real(REALK), intent(out) :: ro_ints(num_sgto_bra,num_contr_bra, &
130                                        num_sgto_ket,num_contr_ket,num_opt)
131!f2py intent(in) :: ang_bra
132!f2py intent(hide) :: num_sgto_bra
133!f2py intent(in) :: mag_bra
134!f2py intent(in) :: ang_ket
135!f2py intent(hide) :: num_sgto_ket
136!f2py intent(in) :: mag_ket
137!f2py intent(hide) :: num_contr_bra
138!f2py intent(hide) :: num_contr_ket
139!f2py intent(hide) :: num_opt
140!f2py intent(in) :: gen_ints
141!f2py depend(num_sgto_bra) :: gen_ints
142!f2py depend(num_sgto_ket) :: gen_ints
143!f2py intent(out) :: ro_ints
144!f2py depend(num_sgto_bra) :: ro_ints
145!f2py depend(num_contr_bra) :: ro_ints
146!f2py depend(num_sgto_ket) :: ro_ints
147!f2py depend(num_contr_ket) :: ro_ints
148!f2py depend(num_opt) :: ro_ints
149    integer ang_bra_one, ang_ket_one          !orbital quantum numbers of bra and ket +1
150    integer ibra, iket, icontr, jcontr, iopt  !incremental recorders
151    integer addr_bra, addr_ket                !addresses of SGTOs on bra and ket centers in Gen1Int
152#if defined(XTIME)
153    real(REALK) curr_time                     !current CPU time
154    ! sets current CPU time
155    call xtimer_set(curr_time)
156#endif
157    ang_bra_one = ang_bra+1
158    ang_ket_one = ang_ket+1
159    ! loops over operators
160    do iopt = 1, num_opt
161      ! loops over contractions on ket center
162      do jcontr = 1, num_contr_ket
163        ! loops over SGTOs on ket center
164        do iket = 1, num_sgto_ket
165          addr_ket = mag_ket(iket)+ang_ket_one
166          ! loops over contractions on bra center
167          do icontr = 1, num_contr_bra
168            ! loops over SGTOs on bra center
169            do ibra = 1, num_sgto_bra
170              addr_bra = mag_bra(ibra)+ang_bra_one
171              ro_ints(ibra,icontr,iket,jcontr,iopt) &
172                = gen_ints(addr_bra,icontr,addr_ket,jcontr,iopt)
173            end do
174          end do
175        end do
176      end do
177    end do
178#if defined(XTIME)
179    ! prints the CPU elapsed time
180    call xtimer_view(curr_time, "reorder_sgto_ints", STDOUT)
181#endif
182    return
183  end subroutine reorder_sgto_ints
184
185  !> \brief reorders the integrals of contracted Cartesian Gaussians on bra or ket center
186  !> \author Bin Gao
187  !> \date 2010-10-07
188  !> \param ang_ket is the orbital quantum number (or angular number) of ket center
189  !> \param num_cgto_ket is the number of basis functions on ket center,
190  !>        (equals to $(\var(ang_ket)+1)(\var(ang_ket)+2)/2\f$)
191  !> \param power_ket contains the Cartesian powers of basis functions on ket center
192  !> \param dim_cgto_bra is the dimension of CGTOs on bra center
193  !> \param num_contr_ket is the number of contractions of ket center
194  !> \param num_opt is the number of operators
195  !> \param gen_ints contains the contracted integrals from Gen1Int
196  !> \return ro_ints contains the reordered integrals according to given \var(power_ket)
197  subroutine reorder_cgtos(ang_ket, num_cgto_ket, power_ket, dim_cgto_bra, &
198                           num_contr_ket, num_opt, gen_ints, ro_ints)
199    use xkind
200    implicit none
201    integer, intent(in) :: ang_ket
202    integer, intent(in) :: num_cgto_ket
203    integer, intent(in) :: power_ket(3,num_cgto_ket)
204    integer, intent(in) :: dim_cgto_bra
205    integer, intent(in) :: num_contr_ket
206    integer, intent(in) :: num_opt
207    real(REALK), intent(in) :: gen_ints(dim_cgto_bra,num_cgto_ket,num_contr_ket,num_opt)
208    real(REALK), intent(out) :: ro_ints(dim_cgto_bra,num_cgto_ket,num_contr_ket,num_opt)
209!f2py intent(in) :: ang_ket
210!f2py intent(hide) :: num_cgto_ket
211!f2py intent(in) :: power_ket
212!f2py intent(hide) :: dim_cgto_bra
213!f2py intent(hide) :: num_contr_ket
214!f2py intent(hide) :: num_opt
215!f2py intent(in) :: gen_ints
216!f2py depend(num_cgto_ket) :: gen_ints
217!f2py intent(out) :: ro_ints
218!f2py depend(dim_cgto_bra) :: ro_ints
219!f2py depend(num_cgto_ket) :: ro_ints
220!f2py depend(num_contr_ket) :: ro_ints
221!f2py depend(num_opt) :: ro_ints
222    integer dang_ket_tri              !twice of the angular number on ket +3
223    integer ibra, iket, jcontr, iopt  !incremental recorders
224    integer addr_ket                  !address of CGTOs on ket center in Gen1Int
225#if defined(XTIME)
226    real(REALK) curr_time             !current CPU time
227    ! sets current CPU time
228    call xtimer_set(curr_time)
229#endif
230    ! index of x^{l}y^{m}z^{n} with l+m+n=angm is 1+m+(2*angm+3-n)*n/2
231    dang_ket_tri = ang_ket+ang_ket+3
232    ! loops over operators
233    do iopt = 1, num_opt
234      ! loops over contractions on ket center
235      do jcontr = 1, num_contr_ket
236        ! loops over CGTOs on ket center
237        do iket = 1, num_cgto_ket
238          addr_ket = 1+power_ket(2,iket)+(dang_ket_tri-power_ket(3,iket))*power_ket(3,iket)/2
239          ! loops over CGTOs on bra center
240          do ibra = 1, dim_cgto_bra
241            ro_ints(ibra,iket,jcontr,iopt) = gen_ints(ibra,addr_ket,jcontr,iopt)
242          end do
243        end do
244      end do
245    end do
246#if defined(XTIME)
247    ! prints the CPU elapsed time
248    call xtimer_view(curr_time, "reorder_cgtos", STDOUT)
249#endif
250    return
251  end subroutine reorder_cgtos
252
253  !> \brief reorders the integrals of contracted Cartesian Gaussians
254  !> \author Bin Gao
255  !> \date 2010-10-07
256  !> \param ang_bra is the orbital quantum number (or angular number) of bra center
257  !> \param num_cgto_bra is the number of basis functions on bra center,
258  !>        (equals to $(\var(ang_bra)+1)(\var(ang_bra)+2)/2\f$)
259  !> \param power_bra contains the Cartesian powers of basis functions on bra center
260  !> \param ang_ket is the orbital quantum number (or angular number) of ket center
261  !> \param num_cgto_ket is the number of basis functions on ket center,
262  !>        (equals to $(\var(ang_ket)+1)(\var(ang_ket)+2)/2\f$)
263  !> \param power_ket contains the Cartesian powers of basis functions on ket center
264  !> \param num_contr_bra is the number of contractions of bra center
265  !> \param num_contr_ket is the number of contractions of ket center
266  !> \param num_opt is the number of operators
267  !> \param gen_ints contains the contracted integrals from Gen1Int
268  !> \return ro_ints contains the reordered integrals according to given \var(power_bra)
269  !>         and \var(powre_ket)
270  subroutine reorder_cgto_ints(ang_bra, num_cgto_bra, power_bra, &
271                               ang_ket, num_cgto_ket, power_ket, &
272                               num_contr_bra, num_contr_ket,     &
273                               num_opt, gen_ints, ro_ints)
274    use xkind
275    implicit none
276    integer, intent(in) :: ang_bra
277    integer, intent(in) :: num_cgto_bra
278    integer, intent(in) :: power_bra(3,num_cgto_bra)
279    integer, intent(in) :: ang_ket
280    integer, intent(in) :: num_cgto_ket
281    integer, intent(in) :: power_ket(3,num_cgto_ket)
282    integer, intent(in) :: num_contr_bra
283    integer, intent(in) :: num_contr_ket
284    integer, intent(in) :: num_opt
285    real(REALK), intent(in) :: gen_ints(num_cgto_bra,num_contr_bra, &
286                                        num_cgto_ket,num_contr_ket,num_opt)
287    real(REALK), intent(out) :: ro_ints(num_cgto_bra,num_contr_bra, &
288                                        num_cgto_ket,num_contr_ket,num_opt)
289!f2py intent(in) :: ang_bra
290!f2py intent(hide) :: num_cgto_bra
291!f2py intent(in) :: power_bra
292!f2py intent(in) :: ang_ket
293!f2py intent(hide) :: num_cgto_ket
294!f2py intent(in) :: power_ket
295!f2py intent(hide) :: num_contr_bra
296!f2py intent(hide) :: num_contr_ket
297!f2py intent(hide) :: num_opt
298!f2py intent(in) :: gen_ints
299!f2py depend(num_cgto_bra) :: gen_ints
300!f2py depend(num_cgto_ket) :: gen_ints
301!f2py intent(out) :: ro_ints
302!f2py depend(num_cgto_bra) :: ro_ints
303!f2py depend(num_contr_bra) :: ro_ints
304!f2py depend(num_cgto_ket) :: ro_ints
305!f2py depend(num_contr_ket) :: ro_ints
306!f2py depend(num_opt) :: ro_ints
307    integer dang_bra_tri, dang_ket_tri        !twice of the angular numbers on bra and ket +3
308    integer ibra, iket, icontr, jcontr, iopt  !incremental recorders
309    integer addr_bra, addr_ket                !addresses of CGTOs on bra and ket centers in Gen1Int
310#if defined(XTIME)
311    real(REALK) curr_time                     !current CPU time
312    ! sets current CPU time
313    call xtimer_set(curr_time)
314#endif
315    ! index of x^{l}y^{m}z^{n} with l+m+n=angm is 1+m+(2*angm+3-n)*n/2
316    dang_bra_tri = ang_bra+ang_bra+3
317    dang_ket_tri = ang_ket+ang_ket+3
318    ! loops over operators
319    do iopt = 1, num_opt
320      ! loops over contractions on ket center
321      do jcontr = 1, num_contr_ket
322        ! loops over CGTOs on ket center
323        do iket = 1, num_cgto_ket
324          addr_ket = 1+power_ket(2,iket)+(dang_ket_tri-power_ket(3,iket))*power_ket(3,iket)/2
325          ! loops over contractions on bra center
326          do icontr = 1, num_contr_bra
327            ! loops over CGTOs on bra center
328            do ibra = 1, num_cgto_bra
329              addr_bra = 1+power_bra(2,ibra)+(dang_bra_tri-power_bra(3,ibra))*power_bra(3,ibra)/2
330              ro_ints(ibra,icontr,iket,jcontr,iopt) &
331                = gen_ints(addr_bra,icontr,addr_ket,jcontr,iopt)
332            end do
333          end do
334        end do
335      end do
336    end do
337#if defined(XTIME)
338    ! prints the CPU elapsed time
339    call xtimer_view(curr_time, "reorder_cgto_ints", STDOUT)
340#endif
341    return
342  end subroutine reorder_cgto_ints
343