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 dumps the information of contracted GTOs, derivatives and operators.
18!!
19!!  2011-03-24, Bin Gao:
20!!  * first version
21
22  !> \brief dumps the information of contracted GTOs and partial derivatives
23  !> \author Bin Gao
24  !> \date 2011-03-24
25  !> \param type_gto specifies the type of GTOs
26  !> \param idx_gto is the atomic index of GTO center
27  !> \param coord_gto contains the coordinates of GTO center
28  !> \param angular_num is the angular number (s=0, p=1, d=2, ...)
29  !> \param num_prim is the number of primitive Gaussians
30  !> \param exponents contains the exponents of primitive Gaussians
31  !> \param num_contr is the number of contractions
32  !> \param contr_coef contains the contraction coefficients
33  !> \param order_part_mag is the order of partial magnetic derivatives
34  !> \param order_part_ram is the order of partial derivatives with respect to
35  !>        total rotational angular momentum
36  !> \param order_part_geo is the order of partial geometric derivatives
37  !> \param io_viewer is the logical unit number of the viewer
38  subroutine dump_gto_pd(type_gto, idx_gto, coord_gto, angular_num,      &
39                         num_prim, exponents, num_contr, contr_coef,     &
40                         order_part_mag, order_part_ram, order_part_geo, &
41                         io_viewer)
42    use xkind
43    implicit none
44    character*(*), intent(in) :: type_gto
45    integer, intent(in) :: idx_gto
46    real(REALK), intent(in) :: coord_gto(3)
47    integer, intent(in) :: angular_num
48    integer, intent(in) :: num_prim
49    real(REALK), intent(in) :: exponents(num_prim)
50    integer, intent(in) :: num_contr
51    real(REALK), intent(in) :: contr_coef(num_contr,num_prim)
52    integer, intent(in) :: order_part_mag
53    integer, intent(in) :: order_part_ram
54    integer, intent(in) :: order_part_geo
55    integer, intent(in) :: io_viewer
56!f2py intent(in) :: type_gto
57!f2py intent(in) :: idx_gto
58!f2py intent(in) :: coord_gto
59!f2py intent(in) :: angular_num
60!f2py intent(hide) :: num_prim
61!f2py intent(in) :: exponents
62!f2py intent(hide) :: num_contr
63!f2py intent(in) :: contr_coef
64!f2py depend(num_prim) :: contr_coef
65!f2py intent(in) :: order_part_mag
66!f2py intent(in) :: order_part_ram
67!f2py intent(in) :: order_part_geo
68!f2py intent(in) :: io_viewer
69    integer iprim  !incremental recorders over primitives
70    ! information of GTOs
71    write(io_viewer,100) "contracted "//type_gto
72    write(io_viewer,100) "atom center:", idx_gto
73    write(io_viewer,110) "coordinates:", coord_gto
74    write(io_viewer,100) "angular number:", angular_num
75    write(io_viewer,100) "number of primitives:", num_prim
76    write(io_viewer,100) "number of contractions:", num_contr
77    write(io_viewer,100) "    exponents      coefficients"
78    do iprim = 1, num_prim
79      write(io_viewer,120) exponents(iprim), contr_coef(1:min(4,num_contr),iprim)
80      write(io_viewer,130) contr_coef(5:num_contr,iprim)
81    end do
82    ! information of partial derivatives
83    if (order_part_mag>0) &
84      write(io_viewer,100) "order of partial magnetic derivatives:", order_part_mag
85    if (order_part_ram>0)                                                       &
86      write(io_viewer,100)                                                      &
87      "order of partial derivatives w.r.t. total rotational angular momentum ", &
88      order_part_ram
89    if (order_part_geo>0) &
90      write(io_viewer,100) "order of partial geometric derivatives", order_part_geo
91    return
92100 format("dump_gto_pd>> ",A,I8)
93110 format("dump_gto_pd>> ",A,3F16.8)
94120 format("dump_gto_pd>> ",5Es16.8)
95130 format("dump_gto_pd>> ",16X,4Es16.8)
96  end subroutine dump_gto_pd
97
98  !> \brief dumps the information of total derivatives
99  !> \author Bin Gao
100  !> \date 2011-03-24
101  !> \param order_mag is the order of total magnetic derivatives
102  !> \param order_ram is the order of total derivatives with respect to
103  !>        total rotational angular momentum
104  !> \param num_cents is the number of differentiated centers
105  !> \param idx_cent contains the indices of differentiated centers
106  !> \param order_cent contains the order of total geometric derivatives of
107  !>        the corresponding differentiated centers
108  !> \param io_viewer is the logical unit number of the viewer
109  subroutine dump_total_derv(order_mag, order_ram, num_cents, &
110                             idx_cent, order_cent, io_viewer)
111    use xkind
112    implicit none
113    integer, intent(in) :: order_mag
114    integer, intent(in) :: order_ram
115    integer, intent(in) :: num_cents
116    integer, intent(in) :: idx_cent(num_cents)
117    integer, intent(in) :: order_cent(num_cents)
118    integer, intent(in) :: io_viewer
119!f2py intent(in) :: order_mag
120!f2py intent(in) :: order_ram
121!f2py intent(hide) :: num_cents
122!f2py intent(in) :: idx_cent
123!f2py intent(in) :: order_cent
124!f2py depend(num_cents) :: order_cent
125!f2py intent(in) :: io_viewer
126    integer icent  !incremental recorders over centers
127    if (order_mag>0) &
128      write(io_viewer,100) "order of total magnetic derivatives:", order_mag
129    if (order_ram>0) write(io_viewer,100)                                     &
130      "order of total derivatives w.r.t. total rotational angular momentum:", &
131      order_ram
132    if (num_cents>0) then
133      write(io_viewer,100) "number of differentiated centers:", num_cents
134      write(io_viewer,110) (idx_cent(icent), order_cent(icent), icent=1,num_cents)
135    end if
136    return
137100 format("dump_total_derv>> ",A,I8)
138110 format("dump_total_derv>> center(order):",5(I6,"(",I3,")"))
139  end subroutine dump_total_derv
140
141  !> \brief dumps the information of Cartesian multipole moments
142  !> \author Bin Gao
143  !> \date 2011-03-24
144  !> \param order_elec is the order of electronic derivatives
145  !> \param idx_diporg is the atomic center of dipole origin (<1 for non-atomic center)
146  !> \param dipole_origin contains the coordinates of dipole origin
147  !> \param scal_const is the scale constant for Cartesian multipole moments
148  !> \param order_mom is the order of Cartesian multipole moments
149  !> \param order_geo_mom is the order of geometric derivatives on dipole origin
150  !> \param io_viewer is the logical unit number of the viewer
151  subroutine dump_carmom(order_elec, idx_diporg, dipole_origin, scal_const, &
152                         order_mom, order_geo_mom, io_viewer)
153    use xkind
154    implicit none
155    integer, intent(in) :: order_elec
156    integer, intent(in) :: idx_diporg
157    real(REALK), intent(in) :: dipole_origin(3)
158    real(REALK), intent(in) :: scal_const
159    integer, intent(in) :: order_mom
160    integer, intent(in) :: order_geo_mom
161    integer, intent(in) :: io_viewer
162!f2py intent(in) :: order_elec
163!f2py intent(in) :: idx_diporg
164!f2py intent(in) :: dipole_origin
165!f2py intent(in) :: scal_const
166!f2py intent(in) :: order_mom
167!f2py intent(in) :: order_geo_mom
168!f2py intent(in) :: io_viewer
169    write(io_viewer,100) "atom center of dipole origin:", idx_diporg
170    write(io_viewer,110) "coordinates of dipole origin:", dipole_origin
171    write(io_viewer,110) "scale constant:", scal_const
172    write(io_viewer,100) "order of multipole moments:", order_mom
173    write(io_viewer,100) "order of partial geometric derivatives on dipole origin:", &
174                         order_geo_mom
175    write(io_viewer,100) "order of electronic derivatives:", order_elec
176    return
177100 format("dump_carmom>> ",A,I8)
178110 format("dump_carmom>> ",A,3F16.8)
179  end subroutine dump_carmom
180
181  !> \brief dumps the information of Dirac delta function integrals
182  !> \author Bin Gao
183  !> \date 2011-03-24
184  !> \param order_elec is the order of electronic derivatives
185  !> \param idx_delta is the atomic center of delta function origin (<1 for non-atomic center)
186  !> \param delta_origin contains the coordinates of delta function origin
187  !> \param idx_diporg is the atomic center of dipole origin (<1 for non-atomic center)
188  !> \param dipole_origin contains the coordinates of dipole origin
189  !> \param scal_const is the scale constant for Dirac delta function
190  !> \param order_mom is the order of Cartesian multipole moments
191  !> \param order_geo_pot is the order of geometric derivatives on delta function origin
192  !> \param order_geo_mom is the order of geometric derivatives on dipole origin
193  !> \param io_viewer is the IO unit to dump the information
194  subroutine dump_delta(order_elec, idx_delta, delta_origin, idx_diporg,     &
195                        dipole_origin, scal_const, order_mom, order_geo_pot, &
196                        order_geo_mom, io_viewer)
197    use xkind
198    implicit none
199    integer, intent(in) :: order_elec
200    integer, intent(in) :: idx_delta
201    real(REALK), intent(in) :: delta_origin(3)
202    integer, intent(in) :: idx_diporg
203    real(REALK), intent(in) :: dipole_origin(3)
204    real(REALK), intent(in) :: scal_const
205    integer, intent(in) :: order_mom
206    integer, intent(in) :: order_geo_pot
207    integer, intent(in) :: order_geo_mom
208    integer, intent(in) :: io_viewer
209!f2py intent(in) :: order_elec
210!f2py intent(in) :: idx_delta
211!f2py intent(in) :: delta_origin
212!f2py intent(in) :: idx_diporg
213!f2py intent(in) :: dipole_origin
214!f2py intent(in) :: scal_const
215!f2py intent(in) :: order_mom
216!f2py intent(in) :: order_geo_pot
217!f2py intent(in) :: order_geo_mom
218!f2py intent(in) :: io_viewer
219    write(io_viewer,100) "atom center of Dirac delta function:", idx_delta
220    write(io_viewer,110) "coordinates of center of Dirac delta function:", &
221                         delta_origin
222    write(io_viewer,110) "scale constant:", scal_const
223    write(io_viewer,100)                                                 &
224      "order of partial geometric derivatives on Dirac delta function:", &
225      order_geo_pot
226    write(io_viewer,100) "atom center of dipole origin:", idx_diporg
227    write(io_viewer,110) "coordinates of dipole origin:", dipole_origin
228    write(io_viewer,100) "order of multipole moments:", order_mom
229    write(io_viewer,100) "order of partial geometric derivatives on dipole origin:", &
230                         order_geo_mom
231    write(io_viewer,100) "order of electronic derivatives:", order_elec
232    return
233100 format("dump_delta>> ",A,I8)
234110 format("dump_delta>> ",A,3F16.8)
235  end subroutine dump_delta
236
237  !> \brief dumps the information of nuclear attraction potential integrals
238  !> \author Bin Gao
239  !> \date 2011-03-24
240  !> \param order_elec is the order of electronic derivatives
241  !> \param idx_nucorg is the atomic center of nuclear potential origin (<1 for non-atomic center)
242  !> \param nucpot_origin contains the coordinates of nuclear potential origin
243  !> \param idx_diporg is the atomic center of dipole origin (<1 for non-atomic center)
244  !> \param dipole_origin contains the coordinates of dipole origin
245  !> \param scal_const is the scale constant for nuclear attraction potential
246  !> \param order_mom is the order of Cartesian multipole moments
247  !> \param order_geo_pot is the order of geometric derivatives on nuclear potential origin
248  !> \param order_geo_mom is the order of geometric derivatives on dipole origin
249  !> \param io_viewer is the IO unit to dump the information
250  subroutine dump_nucpot(order_elec, idx_nucorg, nucpot_origin, idx_diporg,   &
251                         dipole_origin, scal_const, order_mom, order_geo_pot, &
252                         order_geo_mom, io_viewer)
253    use xkind
254    implicit none
255    integer, intent(in) :: order_elec
256    integer, intent(in) :: idx_nucorg
257    real(REALK), intent(in) :: nucpot_origin(3)
258    integer, intent(in) :: idx_diporg
259    real(REALK), intent(in) :: dipole_origin(3)
260    real(REALK), intent(in) :: scal_const
261    integer, intent(in) :: order_mom
262    integer, intent(in) :: order_geo_pot
263    integer, intent(in) :: order_geo_mom
264    integer, intent(in) :: io_viewer
265!f2py intent(in) :: order_elec
266!f2py intent(in) :: idx_nucorg
267!f2py intent(in) :: nucpot_origin
268!f2py intent(in) :: idx_diporg
269!f2py intent(in) :: dipole_origin
270!f2py intent(in) :: scal_const
271!f2py intent(in) :: order_mom
272!f2py intent(in) :: order_geo_pot
273!f2py intent(in) :: order_geo_mom
274!f2py intent(in) :: io_viewer
275    write(io_viewer,100) "atom center of nuclear attraction potential:", idx_nucorg
276    write(io_viewer,110) "coordinates of center of nuclear attraction potential:", &
277                         nucpot_origin
278    write(io_viewer,110) "scale constant:", scal_const
279    write(io_viewer,100)                                                         &
280      "order of partial geometric derivatives on nuclear attraction potential:", &
281      order_geo_pot
282    write(io_viewer,100) "atom center of dipole origin:", idx_diporg
283    write(io_viewer,110) "coordinates of dipole origin:", dipole_origin
284    write(io_viewer,100) "order of multipole moments:", order_mom
285    write(io_viewer,100) "order of partial geometric derivatives on dipole origin:", &
286                         order_geo_mom
287    write(io_viewer,100) "order of electronic derivatives:", order_elec
288    return
289100 format("dump_nucpot>> ",A,I8)
290110 format("dump_nucpot>> ",A,3F16.8)
291  end subroutine dump_nucpot
292
293  !> \brief dumps the information of inverse square distance potential integrals
294  !> \author Bin Gao
295  !> \date 2011-03-24
296  !> \param order_elec is the order of electronic derivatives
297  !> \param idx_isdorg is the atomic center of inverse square distance potential
298  !>        origin (<1 for non-atomic center)
299  !> \param isdpot_origin contains the coordinates of inverse square distance potential origin
300  !> \param idx_diporg is the atomic center of dipole origin (<1 for non-atomic center)
301  !> \param dipole_origin contains the coordinates of dipole origin
302  !> \param scal_const is the scale constant for inverse square distance potential
303  !> \param order_mom is the order of Cartesian multipole moments
304  !> \param order_geo_pot is the order of geometric derivatives on inverse square distance potential origin
305  !> \param order_geo_mom is the order of geometric derivatives on dipole origin
306  !> \param io_viewer is the IO unit to dump the information
307  subroutine dump_isdpot(order_elec, idx_isdorg, isdpot_origin, idx_diporg,   &
308                         dipole_origin, scal_const, order_mom, order_geo_pot, &
309                         order_geo_mom, io_viewer)
310    use xkind
311    implicit none
312    integer, intent(in) :: order_elec
313    integer, intent(in) :: idx_isdorg
314    real(REALK), intent(in) :: isdpot_origin(3)
315    integer, intent(in) :: idx_diporg
316    real(REALK), intent(in) :: dipole_origin(3)
317    real(REALK), intent(in) :: scal_const
318    integer, intent(in) :: order_mom
319    integer, intent(in) :: order_geo_pot
320    integer, intent(in) :: order_geo_mom
321    integer, intent(in) :: io_viewer
322!f2py intent(in) :: order_elec
323!f2py intent(in) :: idx_isdorg
324!f2py intent(in) :: isdpot_origin
325!f2py intent(in) :: idx_diporg
326!f2py intent(in) :: dipole_origin
327!f2py intent(in) :: scal_const
328!f2py intent(in) :: order_mom
329!f2py intent(in) :: order_geo_pot
330!f2py intent(in) :: order_geo_mom
331!f2py intent(in) :: io_viewer
332    write(io_viewer,100) "atom center of inverse square distance potential:", idx_isdorg
333    write(io_viewer,110) "coordinates of center of inverse square distance potential:", &
334                         isdpot_origin
335    write(io_viewer,110) "scale constant:", scal_const
336    write(io_viewer,100)                                                              &
337      "order of partial geometric derivatives on inverse square distance potential:", &
338      order_geo_pot
339    write(io_viewer,100) "atom center of dipole origin:", idx_diporg
340    write(io_viewer,110) "coordinates of dipole origin:", dipole_origin
341    write(io_viewer,100) "order of multipole moments:", order_mom
342    write(io_viewer,100) "order of partial geometric derivatives on dipole origin:", &
343                         order_geo_mom
344    write(io_viewer,100) "order of electronic derivatives:", order_elec
345    return
346100 format("dump_isdpot>> ",A,I8)
347110 format("dump_isdpot>> ",A,3F16.8)
348  end subroutine dump_isdpot
349
350  !> \brief dumps the information of Gaussian charge potential integrals
351  !> \author Bin Gao
352  !> \date 2011-03-24
353  !> \param order_elec is the order of electronic derivatives
354  !> \param idx_gauorg is the atomic center of Gaussian charge potential origin (<1 for non-atomic center)
355  !> \param gaupot_origin contains the coordinates of Gaussian charge potential origin
356  !> \param gaupot_expt is the exponent used in the Gaussian broadening function of the charge
357  !> \param idx_diporg is the atomic center of dipole origin (<1 for non-atomic center)
358  !> \param dipole_origin contains the coordinates of dipole origin
359  !> \param scal_const is the scale constant for Gaussian charge potential
360  !> \param order_mom is the order of Cartesian multipole moments
361  !> \param order_geo_pot is the order of geometric derivatives on gaupot function origin
362  !> \param order_geo_mom is the order of geometric derivatives on dipole origin
363  !> \param io_viewer is the IO unit to dump the information
364  subroutine dump_gaupot(order_elec, idx_gauorg, gaupot_origin, gaupot_expt, &
365                         idx_diporg, dipole_origin, scal_const, order_mom,   &
366                         order_geo_pot, order_geo_mom, io_viewer)
367    use xkind
368    implicit none
369    integer, intent(in) :: order_elec
370    integer, intent(in) :: idx_gauorg
371    real(REALK), intent(in) :: gaupot_origin(3)
372    real(REALK), intent(in) :: gaupot_expt
373    integer, intent(in) :: idx_diporg
374    real(REALK), intent(in) :: dipole_origin(3)
375    real(REALK), intent(in) :: scal_const
376    integer, intent(in) :: order_mom
377    integer, intent(in) :: order_geo_pot
378    integer, intent(in) :: order_geo_mom
379    integer, intent(in) :: io_viewer
380!f2py intent(in) :: order_elec
381!f2py intent(in) :: idx_gauorg
382!f2py intent(in) :: gaupot_origin
383!f2py intent(in) :: gaupot_expt
384!f2py intent(in) :: idx_diporg
385!f2py intent(in) :: dipole_origin
386!f2py intent(in) :: scal_const
387!f2py intent(in) :: order_mom
388!f2py intent(in) :: order_geo_pot
389!f2py intent(in) :: order_geo_mom
390!f2py intent(in) :: io_viewer
391    write(io_viewer,100) "atom center of Gaussian charge potential:", idx_gauorg
392    write(io_viewer,110) "coordinates of center of Gaussian charge potential:", &
393                         gaupot_origin
394    write(io_viewer,110) "exponent in Gaussian broadening function of the charge", &
395                         gaupot_expt
396    write(io_viewer,110) "scale constant:", scal_const
397    write(io_viewer,100)                                                      &
398      "order of partial geometric derivatives on Gaussian charge potential:", &
399      order_geo_pot
400    write(io_viewer,100) "atom center of dipole origin:", idx_diporg
401    write(io_viewer,110) "coordinates of dipole origin:", dipole_origin
402    write(io_viewer,100) "order of multipole moments:", order_mom
403    write(io_viewer,100) "order of partial geometric derivatives on dipole origin:", &
404                         order_geo_mom
405    write(io_viewer,100) "order of electronic derivatives:", order_elec
406    return
407100 format("dump_gaupot>> ",A,I8)
408110 format("dump_gaupot>> ",A,3F16.8)
409  end subroutine dump_gaupot
410
411  !> \brief dumps the information of diamagnetic spin-orbit coupling integrals
412  !> \author Bin Gao
413  !> \date 2011-03-24
414  !> \param order_elec is the order of electronic derivatives
415  !> \param idx_dso contains the atomic centers of diamagnetic spin-orbit coupling
416  !>        origins (<1 for non-atomic center)
417  !> \param coord_dso contains the coordinates of diamagnetic spin-orbit coupling origins
418  !> \param scal_const is the scale constant for diamagnetic spin-orbit coupling
419  !> \param order_geo_dso contains the orders of geometric derivatives on diamagnetic
420  !>        spin-orbit coupling origins
421  !> \param io_viewer is the IO unit to dump the information
422  subroutine dump_dso(order_elec, idx_dso, coord_dso, scal_const, &
423                      order_geo_dso, io_viewer)
424    use xkind
425    implicit none
426    integer, intent(in) :: order_elec
427    integer, intent(in) :: idx_dso(2)
428    real(REALK), intent(in) :: coord_dso(3,2)
429    real(REALK), intent(in) :: scal_const
430    integer, intent(in) :: order_geo_dso(2)
431    integer, intent(in) :: io_viewer
432!f2py intent(in) :: order_elec
433!f2py intent(in) :: idx_dso
434!f2py intent(in) :: coord_dso
435!f2py intent(in) :: scal_const
436!f2py intent(in) :: order_geo_dso
437!f2py intent(in) :: io_viewer
438    write(io_viewer,100) "first atom center:", idx_dso(1)
439    write(io_viewer,110) "coordinates of the first atom center:", coord_dso(:,1)
440    write(io_viewer,100)                                                  &
441      "order of partial geometric derivatives on the first atom center:", &
442      order_geo_dso(1)
443    write(io_viewer,100) "second atom center:", idx_dso(2)
444    write(io_viewer,110) "coordinates of the second atom center:", coord_dso(:,2)
445    write(io_viewer,100)                                                   &
446      "order of partial geometric derivatives on the second atom center:", &
447      order_geo_dso(2)
448    write(io_viewer,110) "scale constant:", scal_const
449    write(io_viewer,100) "order of electronic derivatives:", order_elec
450    return
451100 format("dump_dso>> ",A,I8)
452110 format("dump_dso>> ",A,3F16.8)
453  end subroutine dump_dso
454
455  !> \brief dumps the information of effective core potential (ECP)
456  !> \author Bin Gao
457  !> \date 2011-03-24
458  !> \param ecp_local indicates if it is the local or non-local part
459  !> \param idx_ecp is the atomic index of ECP center
460  !> \param coord_ecp is the coordinates of ECP center
461  !> \param core_angl is the number of core electrons (local part) or angular number (non-local part)
462  !> \param num_gauss is the number of Gaussians combining the radical function
463  !> \param gau_rpower contains the powers of r in Gaussians
464  !> \param gau_exponent contains the exponents of Gaussians
465  !> \param gau_coeff contains the combination coefficients of Gaussians
466  !> \param io_viewer is the logical unit number of the viewer
467  subroutine dump_ecp(ecp_local, idx_ecp, coord_ecp, core_angl, num_gauss, &
468                      gau_rpower, gau_exponent, gau_coeff, io_viewer)
469    use xkind
470    implicit none
471    logical, intent(in) :: ecp_local
472    integer, intent(in) :: idx_ecp
473    real(REALK), intent(in) :: coord_ecp(3)
474    integer, intent(in) :: core_angl
475    integer, intent(in) :: num_gauss
476    integer, intent(in) :: gau_rpower(num_gauss)
477    real(REALK), intent(in) :: gau_exponent(num_gauss)
478    real(REALK), intent(in) :: gau_coeff(num_gauss)
479    integer, intent(in) :: io_viewer
480!f2py intent(in) :: ecp_local
481!f2py intent(in) :: idx_ecp
482!f2py intent(in) :: coord_ecp
483!f2py intent(in) :: core_angl
484!f2py intent(hide) :: num_gauss
485!f2py intent(in) :: gau_rpower
486!f2py intent(in) :: gau_exponent
487!f2py depend(num_gauss) :: gau_exponent
488!f2py intent(in) :: gau_coeff
489!f2py depend(num_gauss) :: gau_coeff
490!f2py intent(in) :: io_viewer
491    integer igau  !incremental recorder over Gaussians
492    if (ecp_local) then
493      write(io_viewer,100) "operator>> effective core potential (local part)"
494      write(io_viewer,100) "number of core electrons:", core_angl
495    else
496      write(io_viewer,100) "operator>> effective core potential (non-local part)"
497      write(io_viewer,100) "angular number:", core_angl
498    end if
499    write(io_viewer,100) "atom center of ECP:", idx_ecp
500    write(io_viewer,110) "coordinates of ECP center:", coord_ecp
501    write(io_viewer,100) "    powers      exponents      coefficients"
502    do igau = 1, num_gauss
503      write(io_viewer,120) -gau_rpower(igau), gau_exponent(igau), gau_coeff(igau)
504    end do
505    return
506100 format("dump_ecp>> ",A,I8)
507110 format("dump_ecp>> ",A,3F16.8)
508120 format("dump_ecp>> ",I4,2F16.8)
509  end subroutine dump_ecp
510
511  !> \brief dumps the information of model core potential (Version 1) -- potential
512  !> \author Bin Gao
513  !> \date 2011-03-24
514  !> \param idx_mcp is the atomic index of MCP1 center
515  !> \param coord_mcp is the coordinates of MCP1 center
516  !> \param charge_mcp is the effective core charge
517  !> \param num_gauss is the number of Gaussians in the model potential
518  !> \param gau_rpower contains the powers of r in Gaussians
519  !> \param gau_exponent contains the exponents of Gaussians
520  !> \param gau_coeff contains the combination coefficients of Gaussians
521  !> \param io_viewer is the logical unit number of the viewer
522  subroutine dump_mcp1_pot(idx_mcp, coord_mcp, charge_mcp, num_gauss, &
523                           gau_rpower, gau_exponent, gau_coeff, io_viewer)
524    use xkind
525    implicit none
526    integer, intent(in) :: idx_mcp
527    real(REALK), intent(in) :: coord_mcp(3)
528    real(REALK), intent(in) :: charge_mcp
529    integer, intent(in) :: num_gauss
530    integer, intent(in) :: gau_rpower(num_gauss)
531    real(REALK), intent(in) :: gau_exponent(num_gauss)
532    real(REALK), intent(in) :: gau_coeff(num_gauss)
533    integer, intent(in) :: io_viewer
534!f2py intent(in) :: idx_mcp
535!f2py intent(in) :: coord_mcp
536!f2py intent(in) :: charge_mcp
537!f2py intent(hide) :: num_gauss
538!f2py intent(in) :: gau_rpower
539!f2py intent(in) :: gau_exponent
540!f2py depend(num_gauss) :: gau_exponent
541!f2py intent(in) :: gau_coeff
542!f2py depend(num_gauss) :: gau_coeff
543!f2py intent(in) :: io_viewer
544    integer igau  !incremental recorder over Gaussians
545    write(io_viewer,100) "operator>> model core potential (Version 1) -- potential"
546    write(io_viewer,100) "atom center of MCP1:", idx_mcp
547    write(io_viewer,110) "coordinates of MCP1 center:", coord_mcp
548    write(io_viewer,110) "effective core charge:", charge_mcp
549    write(io_viewer,100) "    powers      exponents      coefficients"
550    do igau = 1, num_gauss
551      write(io_viewer,120) gau_rpower(igau), gau_exponent(igau), gau_coeff(igau)
552    end do
553    return
554100 format("dump_mcp1_pot>> ",A,I8)
555110 format("dump_mcp1_pot>> ",A,3F16.8)
556120 format("dump_mcp1_pot>> ",I4,2F16.8)
557  end subroutine dump_mcp1_pot
558
559  !> \brief dumps the information of model core potential (Version 1) -- core orbitals
560  !> \author Bin Gao
561  !> \date 2011-03-24
562  !> \param idx_mcp is the atomic index of MCP1 center
563  !> \param coord_mcp is the coordinates of MCP1 center
564  !> \param core_angnum is the angular number of core orbitals
565  !> \param param_eshfit is the parameter of energy shift operator
566  !> \param core_nexp is the number of exponents of core orbtials
567  !> \param core_exponent contains the orbital exponents
568  !> \param core_nshell is the number of shells of core orbitals
569  !> \param core_coeff contains the expansion coefficients
570  !> \param io_viewer is the logical unit number of the viewer
571  subroutine dump_mcp1_core(idx_mcp, coord_mcp, core_angnum, param_eshfit, &
572                            core_nexp, core_exponent, core_nshell,         &
573                            core_coeff, io_viewer)
574    use xkind
575    implicit none
576    integer, intent(in) :: idx_mcp
577    real(REALK), intent(in) :: coord_mcp(3)
578    integer, intent(in) :: core_angnum
579    real(REALK), intent(in) :: param_eshfit
580    integer, intent(in) :: core_nexp
581    real(REALK), intent(in) :: core_exponent(core_nexp)
582    integer, intent(in) :: core_nshell
583    real(REALK), intent(in) :: core_coeff(core_nshell,core_nexp)
584    integer, intent(in) :: io_viewer
585!f2py intent(in) :: idx_mcp
586!f2py intent(in) :: coord_mcp
587!f2py intent(in) :: core_angnum
588!f2py intent(in) :: param_eshfit
589!f2py intent(hide) :: core_nexp
590!f2py intent(in) :: core_exponent
591!f2py intent(hide) :: core_nshell
592!f2py intent(in) :: core_coeff
593!f2py depend(core_nexp) :: core_coeff
594!f2py intent(in) :: io_viewer
595    integer iexp  !incremental recorder over orbital exponents
596    write(io_viewer,100) "operator>> model core potential (Version 1) -- core orbitals"
597    write(io_viewer,100) "atom center of MCP1:", idx_mcp
598    write(io_viewer,110) "coordinates of MCP1 center:", coord_mcp
599    write(io_viewer,100) "angular number of core orbitals:", core_angnum
600    write(io_viewer,110) "parameter of energy shift operator:", param_eshfit
601    write(io_viewer,100) "number of exponents of core orbtials:", core_nexp
602    write(io_viewer,100) "number of shells of core orbitals:", core_nshell
603    write(io_viewer,100) "    exponents      coefficients"
604    do iexp = 1, core_nexp
605      write(io_viewer,120) core_exponent(iexp), core_coeff(1:min(4,core_nshell),iexp)
606      write(io_viewer,130) core_coeff(5:core_nshell,iexp)
607    end do
608    return
609100 format("dump_mcp1_core>> ",A,I8)
610110 format("dump_mcp1_core>> ",A,3F16.8)
611120 format("dump_dump_mcp1_core>> ",5F14.7)
612130 format("dump_dump_mcp1_core>> ",14X,4F14.7)
613  end subroutine dump_mcp1_core
614