Added the first tutorial for plugins

docs/source/appendix/contributors.rst

@@ -46,6 +46,7 @@ The following people have contributed to this project (by alphabetical order):
 * Nicolas Renon
 * Lorenzo Tenti
 * Julien Toulouse
+* Diata Traoré
 * Mikaël Véril
 
 

docs/source/index.rst

@@ -39,9 +39,10 @@
    programmers_guide/programming
    programmers_guide/ezfio
    programmers_guide/plugins
+   programmers_guide/plugins_tuto_intro
+   programmers_guide/plugins_tuto_I
    programmers_guide/new_ks
    programmers_guide/index
-   programmers_guide/plugins
 
 
 .. toctree::
@@ -52,5 +53,6 @@
    appendix/benchmarks
    appendix/license
    appendix/contributors
+   appendix/references
 
 

docs/source/programmers_guide/plugins_tuto_I.rst

@@ -0,0 +1 @@
+.. include:: ../../../plugins/local/tuto_plugins/tuto_I/tuto_I.rst

docs/source/programmers_guide/plugins_tuto_intro.rst

@@ -0,0 +1 @@
+.. include:: ../../../plugins/README.rst

plugins/README.rst

@@ -0,0 +1,131 @@
+==============================
+Tutorial for creating a plugin
+==============================
+
+Introduction: what is a plugin, and what tutorial will be about ?
+=================================================================
+
+The |QP| is split into two kinds of routines/global variables (i.e. *providers*): 
+ 1)  the **core modules** locatedin qp2/src/, which contains all the bulk of a quantum chemistry software (integrals, matrix elements between Slater determinants, linear algebra routines, DFT stuffs etc..)
+ 2) the **plugins** which are external routines/*providers* connected to the qp2/src/ routines/*providers*.
+
+More precisely, a **plugin** of the |QP| is a directory where you can create routines, 
+providers and executables that use all the global variables/functions/routines already created 
+in the modules of qp2/src or in other plugins. 
+
+Instead of giving a theoretical lecture on what is a plugin, 
+we will go through a series of examples that allow you to do the following thing: 
+
+1)   print out **one- and two-electron integrals** on the AO/MO basis, creates two providers which manipulate these objects, print out these providers, 
+
+2)  browse the **Slater determinants stored** in the |EZFIO| wave function and compute their matrix elements, 
+
+3) build the **Hamiltonian matrix** and **diagonalize** it either with **Lapack or Davidson**,
+
+4)  print out the **one- and two-electron rdms**, 
+
+5)   obtain the **AOs** and **MOs** on the **DFT grid**, together with the **density**,
+
+How the tutorial will be done
+-----------------------------
+
+This tuto is as follows: 
+
+ 1)  you **READ THIS FILE UNTIL THE END** in order to get the big picture and vocabulary, 
+
+ 2) you go to the directory :file:`qp2/plugins/tuto_plugins/` and you will find detailed tutorials for each of the 5 examples. 
+
+Creating a plugin: the basic
+----------------------------
+
+The first thing to do is to be in the QPSH mode: you execute the qp2/bin/qpsh script that essentially loads all 
+the environement variables and allows for the completion of command lines in bash (that is an AMAZING feature :) 
+
+Then, you need to known **where** you want to create your plugin, and what is the **name** of the plugin. 
+
+.. important::
+
+  The plugins are **NECESSARILY** located in qp2/plugins/, and from there you can create any structures of directories.
+
+
+Ex: If you want to create a plugin named "my_fancy_plugin" in the directory plugins/plugins_test/, 
+this goes with the command 
+
+.. code:: bash
+
+    qp plugins create -n my_fancy_plugin -r plugins_test/
+
+Then, to create the plugin of your dreams, the two questions you need to answer are the following: 
+
+1) What do I **need** to compute what I want, which means what are the **objects** that I need ?
+
+   There are two kind of objects:
+
+   + the *routines/functions*: 
+
+     Ex: Linear algebra routines, integration routines etc ...
+
+   + the global variables which are called the *providers*: 
+
+     Ex: one-electron integrals, Slater determinants, density matrices etc ...
+
+2) **Where do I find** these objects ? 
+
+   The objects (routines/functions/providers) are necessarily created in other *modules/plugins*.  
+
+.. seealso::
+
+   The routine :c:func:`lapack_diagd` (which diagonalises a real hermitian matrix) is located in the file 
+   :file:`qp2/src/utils/linear_algebra.irp.f` 
+   therefore it "belongs" to the module :ref:`module_utils`
+
+   The routine :c:func:`ao_to_mo` (which converts a given matrix A from the AO basis to the MO basis) is located in the file 
+   :file:`qp2/src/mo_one_e_ints/ao_to_mo.irp.f`
+   therefore it "belongs" to the module :ref:`module_mo_one_e_ints`
+
+   The provider :c:data:`ao_one_e_integrals` (which is the integrals of one-body part of H on the AO basis) is located in the file 
+   :file:`qp2/src/ao_one_e_ints/ao_one_e_ints.irp.f` 
+   therefore it belongs to the module :ref:`module_ao_one_e_ints`
+
+   The provider :c:data:`one_e_dm_mo_beta_average` (which is the state average beta density matrix on the MO basis) is located in the file 
+   :file:`qp2/src/determinants/density_matrix.irp.f`
+   therefore it belongs to the module :ref:`module_determinants`
+
+To import all the variables that you need, you just need to write the name of the plugins in the :file:`NEED` file .
+
+To import all the variables/routines of the module :ref:`module_utils`, :ref:`module_determinants` and :ref:`module_mo_one_e_ints`, the  :file:`NEED` file you will need is simply the following:
+
+.. code:: bash
+
+ cat NEED
+
+ utils
+ determinants
+ mo_one_e_ints
+
+
+.. important::
+
+  There are **many** routines/providers in the core modules of QP. 
+
+  Nevertheless, as everything is coded with the |IRPF90|, you can use the following amazing tools: :command:`irpman`
+
+  :command:`irpman` can be used in command line in bash to obtain all the info on a routine or variable ! 
+
+
+Example: execute the following command line : 
+
+.. code:: bash
+
+  irpman ao_one_e_integrals
+
+Then all the information you need on :c:data:`ao_one_e_integrals` will appear on the screen. 
+This includes
+
+ - **where** the provider is created, (*i.e.* the actual file where the provider is designed)
+ - the **type** of the provider (*i.e.* a logical, integer etc ...)
+ - the **dimension** if it is an array, 
+ - what other *providers* are **needed** to build this provider, 
+ - what other *providers* **need** this provider. 
+
+

plugins/local/cipsi_tc_bi_ortho/selection.irp.f

@@ -960,7 +960,7 @@ subroutine fill_buffer_double(i_generator, sp, h1, h2, bannedOrb, banned, fock_d
 !          endif
           e_pert(istate) = 0.25 * val / delta_E
 !          e_pert(istate) = 0.5d0 * (tmp - delta_E)
-          if(dsqrt(dabs(tmp)).gt.1.d-4.and.dabs(alpha_h_psi).gt.1.d-4)then
+          if(dsqrt(tmp).gt.1.d-4.and.dabs(psi_h_alpha).gt.1.d-4)then
            coef(istate)   = e_pert(istate) / psi_h_alpha
           else
            coef(istate)   = alpha_h_psi / delta_E

plugins/local/tuto_plugins/H2.xyz

@@ -0,0 +1,6 @@
+2
+H2, equilibrium geometry
+H   0.0   0.0    0.
+H   0.0   0.0    0.74
+
+

plugins/local/tuto_plugins/n2.xyz

@@ -0,0 +1,4 @@
+2
+N2 Geo: Experiment Mult: 1 symmetry: 14
+N    0.0    0.0    0.5488
+N    0.0    0.0    -0.5488

plugins/local/tuto_plugins/tuto_I/print_one_e_h.irp.f

@@ -0,0 +1,20 @@
+program my_program_to_print_stuffs
+  implicit none
+  BEGIN_DOC
+! TODO : Put the documentation of the program here
+  END_DOC
+ integer :: i,j
+ print*,'AO integrals '
+ do i = 1, ao_num
+  do j = 1, ao_num
+   print*,j,i,ao_one_e_integrals(j,i)
+  enddo
+ enddo
+
+ print*,'MO integrals '
+ do i = 1, mo_num
+  do j = 1, mo_num
+   print*,j,i,mo_one_e_integrals(j,i)
+  enddo
+ enddo
+end

plugins/local/tuto_plugins/tuto_I/print_traces_on_e.irp.f

@@ -0,0 +1,24 @@
+program my_program
+  implicit none
+  BEGIN_DOC
+! This program is there essentially to show how one can use providers in programs 
+  END_DOC
+ integer :: i,j
+ double precision :: accu
+ print*,'Trace on the AO basis '
+ print*,trace_ao_one_e_ints
+ print*,'Trace on the AO basis after projection on the MO basis'
+ print*,trace_ao_one_e_ints_from_mo
+ print*,'Trace of MO integrals '
+ print*,trace_mo_one_e_ints
+ print*,'ao_num = ',ao_num
+ print*,'mo_num = ',mo_num
+ if(ao_num .ne. mo_num)then
+  print*,'The AO basis and MO basis are different ...'
+  print*,'Trace on the AO basis should not be the same as Trace of MO integrals'
+  print*,'Only the second one must be equal to the trace on the MO integrals'
+ else 
+  print*,'The AO basis and MO basis are the same !'
+  print*,'All traces should coincide '
+ endif
+end

plugins/local/tuto_plugins/tuto_I/print_two_e_h.irp.f

@@ -0,0 +1,32 @@
+program my_program_to_print_stuffs
+  implicit none
+  BEGIN_DOC
+! TODO : Put the documentation of the program here
+  END_DOC
+ integer :: i,j,k,l
+ double precision :: integral 
+ double precision :: get_ao_two_e_integral, get_two_e_integral ! declaration of the functions 
+ print*,'AO integrals, physicist notations : <i j|k l>'
+ do i = 1, ao_num
+  do j = 1, ao_num
+   do k = 1, ao_num
+    do l = 1, ao_num
+     integral = get_ao_two_e_integral(i, j, k, l, ao_integrals_map)
+     print*,i,j,k,l,integral 
+    enddo
+   enddo
+  enddo
+ enddo
+
+ print*,'MO integrals, physicist notations : <i j|k l>'
+ do i = 1, mo_num
+  do j = 1, mo_num
+   do k = 1, mo_num
+    do l = 1, mo_num
+     integral = get_two_e_integral(i, j, k, l, mo_integrals_map)
+     print*,i,j,k,l,integral 
+    enddo
+   enddo
+  enddo
+ enddo
+end

plugins/local/tuto_plugins/tuto_I/traces_one_e.irp.f

@@ -0,0 +1,111 @@
+
+! This file is an example of the kind of manipulations that you can do with providers 
+!
+
+!!!!!!!!!!!!!!!!!!!!!!!!!! Main providers useful for the program !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+!!!               type             name 
+BEGIN_PROVIDER [ double precision, trace_mo_one_e_ints]
+ implicit none
+ BEGIN_DOC
+! trace_mo_one_e_ints = Trace of the one-electron integrals on the MO basis
+!
+!                     = sum_i mo_one_e_integrals(i,i)
+ END_DOC
+ integer :: i
+ trace_mo_one_e_ints = 0.d0
+ do i = 1, mo_num
+  trace_mo_one_e_ints += mo_one_e_integrals(i,i)
+ enddo
+END_PROVIDER  
+
+BEGIN_PROVIDER [ double precision, trace_ao_one_e_ints]
+ implicit none
+ BEGIN_DOC
+! trace_ao_one_e_ints = Trace of the one-electron integrals on the AO basis taking into account the non orthogonality
+!
+! Be aware that the trace of an operator in a non orthonormal basis is Tr(A S^{-1}) = \sum_{m,n}(A_mn S^{-1}_mn) 
+! 
+! WARNING: it is equal to the trace on the MO basis if and only if the AO basis and MO basis 
+!          have the same number of functions
+ END_DOC
+ integer :: i,j
+ double precision, allocatable :: inv_overlap_times_integrals(:,:) ! = h S^{-1}
+ allocate(inv_overlap_times_integrals(ao_num,ao_num))
+ ! routine that computes the product of two matrices, you can check it with 
+ ! irpman get_AB_prod
+ call get_AB_prod(ao_one_e_integrals,ao_num,ao_num,s_inv,ao_num,inv_overlap_times_integrals)
+ ! Tr(inv_overlap_times_integrals) = Tr(h S^{-1})
+ trace_ao_one_e_ints = 0.d0
+ do i = 1, ao_num
+  trace_ao_one_e_ints += inv_overlap_times_integrals(i,i)
+ enddo
+ !
+ ! testing the formula Tr(A S^{-1}) = \sum_{m,n}(A_mn S^{-1}_mn)
+ double precision :: test
+ test = 0.d0
+ do i = 1, ao_num
+  do j = 1, ao_num
+   test += ao_one_e_integrals(j,i) * s_inv(i,j)
+  enddo
+ enddo
+ if(dabs(accu - trace_ao_one_e_ints).gt.1.d-12)then
+  print*,'Warning ! '
+  print*,'Something is wrong because Tr(AB) \ne sum_{mn}A_mn B_nm'
+ endif
+END_PROVIDER
+
+BEGIN_PROVIDER [ double precision, trace_ao_one_e_ints_from_mo]
+ implicit none
+ BEGIN_DOC
+! trace_ao_one_e_ints_from_mo = Trace of the one-electron integrals on the AO basis after projection on the MO basis 
+!
+!                             = Tr([SC h {SC}^+] S^{-1})
+!
+!                             = Be aware that the trace of an operator in a non orthonormal basis is = Tr(A S^{-1}) where S is the metric
+! Must be equal to the trace_mo_one_e_ints 
+ END_DOC
+ integer :: i
+ double precision, allocatable :: inv_overlap_times_integrals(:,:)
+ allocate(inv_overlap_times_integrals(ao_num,ao_num))
+ ! Using the provider ao_one_e_integrals_from_mo = [SC h {SC}^+] 
+ call get_AB_prod(ao_one_e_integrals_from_mo,ao_num,ao_num,s_inv,ao_num,inv_overlap_times_integrals)
+ ! inv_overlap_times_integrals = [SC h {SC}^+] S^{-1}
+ trace_ao_one_e_ints_from_mo = 0.d0
+ ! Computing the trace
+ do i = 1, ao_num
+  trace_ao_one_e_ints_from_mo += inv_overlap_times_integrals(i,i)
+ enddo
+END_PROVIDER  
+
+!!!!!!!!!!!!!!!!!!!!!!!!!!! Additional providers to check some stuffs !!!!!!!!!!!!!!!!!!!!!!!!!
+
+BEGIN_PROVIDER [ double precision, ao_one_e_int_no_ov_from_mo, (ao_num, ao_num) ]
+ BEGIN_DOC
+ ! ao_one_e_int_no_ov_from_mo = C mo_one_e_integrals C^T 
+ !
+ ! WARNING : NON EQUAL TO ao_one_e_integrals due to the non orthogonality 
+ END_DOC
+ call mo_to_ao_no_overlap(mo_one_e_integrals,mo_num,ao_one_e_int_no_ov_from_mo,ao_num) 
+END_PROVIDER 
+
+BEGIN_PROVIDER [ double precision, ao_one_e_int_no_ov_from_mo_ov_ov, (ao_num, ao_num)]
+ BEGIN_DOC
+ ! ao_one_e_int_no_ov_from_mo_ov_ov = S ao_one_e_int_no_ov_from_mo S = SC mo_one_e_integrals (SC)^T
+ !
+ ! EQUAL TO ao_one_e_integrals ONLY IF ao_num = mo_num
+ END_DOC
+ double precision, allocatable :: tmp(:,:)
+ allocate(tmp(ao_num, ao_num))
+ call get_AB_prod(ao_overlap,ao_num,ao_num,ao_one_e_int_no_ov_from_mo,ao_num,tmp)
+ call get_AB_prod(tmp,ao_num,ao_num,ao_overlap,ao_num,ao_one_e_int_no_ov_from_mo_ov_ov)
+END_PROVIDER 
+
+BEGIN_PROVIDER [ double precision, c_t_s_c, (mo_num, mo_num)]
+ implicit none
+ BEGIN_DOC
+! C^T S C = should be the identity 
+ END_DOC 
+ call get_AB_prod(mo_coef_transp,mo_num,ao_num,S_mo_coef,mo_num,c_t_s_c)
+END_PROVIDER 
+