<div dir="ltr"><div><div>Looks fine to me. It would be preferrable to have the same set of environ routines and of "hooks" for both CP and PW, but this is not going to be easy. <br><br></div>Dead code and redundant stuff can (should) be removed<br><br></div>Paolo<br></div><div class="gmail_extra"><br><div class="gmail_quote">On Tue, Oct 6, 2015 at 12:12 PM, Oliviero Andreussi <span dir="ltr"><<a href="mailto:oliviero.andreussi@usi.ch" target="_blank">oliviero.andreussi@usi.ch</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Dear QE-developers and, in particular, CP-developers,<br>
<br>
I am in the process of interfacing with CP the module we are developing,<br>
Environ (<a href="http://www.quantum-environment.org" rel="noreferrer" target="_blank">www.quantum-environment.org</a>), to model environment effects in<br>
first-principles simulations. I have been following the same strategy<br>
that was followed with PW, i.e. via the so-called plugins structure,<br>
which is meant to provide a more general platform for developers to<br>
introduce inside the code additional terms to the<br>
potential/energy/forces that may depend on the electronic density of the<br>
system (and thus need to be updated at each SCF/CP step). The idea is<br>
that a few calls are inserted in appropriate places inside the main<br>
programs (PW, and now CP) to plugin_#something# subroutines, which are<br>
empty subroutines that can then be patched by the external modules.<br>
<br>
Following what was done for PW (mostly by Layla and myself), I have now<br>
modified CP to allow the use of plugins. In particular, I have inserted<br>
in a few subroutines of CP a few calls to new plugin_#something#<br>
subroutines, I have created this new empy files and I have added them to<br>
the CPV/src/Makefile. Clearly, these modifications alone will/should not<br>
change the main program in any sensible way, apart from making it more<br>
ugly. Nonetheless, before commiting these changes, I would like to be<br>
sure that no harm is done and that the other developers of CP are ok<br>
with my modifications. At the bottom of this email I copy the diff of my<br>
development versions with respect to the svn files, the main changes are<br>
last, in the subroutine vofrho_x in file vofrho.f90. Please let me know<br>
if you have any comment/suggestion/request, if no remarks are done, I<br>
will commit what I have by the end of the week.<br>
<br>
A few side notes:<br>
<br>
1) I also created a plugin_utilities.f90 file which contains a few<br>
subroutines that do some reciprocal-space calculations or are interfaces<br>
for subroutines that do so (compute an electrostatic potential, a<br>
gradient, forces, etc.). Although there are already subroutines that do<br>
something similar, unfortunately their structures are not general enough<br>
to be used by my modules and I had to create more general versions or<br>
interfaces that could be called with less arguments. In PW this<br>
subroutines were included in the appropriate files, (e.g. v_h_of_rho_r<br>
inside PW/src/v_of_rho.f90 to compute a Hartree potential from a density<br>
in real-space, or the external_gradient interface in<br>
PW/src/grad_corr.f90 to compute a gradient without the need of<br>
specifying reciprocal-space information). In CP this could also be done,<br>
but I am not sure it is really useful, so I decided to put everything<br>
inside a specific file. This way, it should be easier to clean it up or<br>
remove it, in case no one will need it anymore.<br>
<br>
2) While looking for subroutines to do reciprocal-space calculations, I<br>
found the vofps_x and vofloc_x routines in potentials.f90 that are never<br>
called by CP and that are potentially buggy (the loops over g vectors<br>
run on the dense grid up to ngm, while rhops and vps are only defined on<br>
the smooth grid up to ngms). Maybe they should be removed from the code<br>
or fixed, if used by anyone.<br>
<br>
3) Environ contains also the quantum-enthalpy and quantum-surface terms<br>
developed by Matteo Cococcioni and co-workers and that are already<br>
inside the CP code (vol_clu.f90 and related calls). In principles, it<br>
would be possible to remove these parts from the code and only leave<br>
them for Environ users.<br>
<br>
Best,<br>
<br>
Oliviero Andreussi<br>
Senior Postdoctoral Associate<br>
École Polytechnique Fédérale de Lausanne (EPFL) and<br>
Università della Svizzera Italiana (USI) of Lugano, Switzerland<br>
<a href="mailto:oliviero.andreussi@usi.ch">oliviero.andreussi@usi.ch</a><br>
<a href="mailto:oliviero.andreussi@epfl.ch">oliviero.andreussi@epfl.ch</a><br>
<br>
<br>
Here is the list of empty, patchable, subroutines I plan to add to CPV/src/<br>
<br>
?       plugin_clean.f90<br>
?       plugin_init_ions.f90<br>
?       plugin_get_potential.f90<br>
?       plugin_print_info.f90<br>
?       plugin_add_potential.f90<br>
?       plugin_print_energies.f90<br>
?       plugin_energy.f90<br>
?       plugin_int_forces.f90<br>
?       plugin_init_cell.f90<br>
?       plugin_read_input.f90<br>
?       plugin_utilities.f90<br>
?       plugin_clock.f90<br>
?       plugin_init_base.f90<br>
<br>
and here is the diff of the files I plan to commit with respect to SVN<br>
updated now<br>
<br>
===================================================================<br>
--- cpr.f90    (revision 11775)<br>
+++ cpr.f90    (working copy)<br>
@@ -295,6 +295,10 @@<br>
       IF ( ( tfor .OR. tfirst ) .AND. tefield ) CALL efield_update( eigr )<br>
       IF ( ( tfor .OR. tfirst ) .AND. tefield2 ) CALL efield_update2(<br>
eigr )<br>
       !<br>
+     ! ... pass ions information to plugins<br>
+     !<br>
+     CALL plugin_init_ions( tau0 )<br>
+     !<br>
       IF ( lda_plus_u ) then<br>
          ! forceh    ! Forces on ions due to Hubbard U<br>
          forceh=0.0d0<br>
@@ -801,6 +805,8 @@<br>
             IF ( tefield )  CALL efield_update( eigr )<br>
             IF ( tefield2 ) CALL efield_update2( eigr )<br>
             !<br>
+           CALL plugin_init_ions( tau0 )<br>
+           !<br>
             lambdam = lambda<br>
             !<br>
             CALL move_electrons( nfi, tfirst, tlast, bg(:,1), bg(:,2),<br>
bg(:,3),&<br>
@@ -1128,6 +1134,8 @@<br>
    !<br>
    CALL print_clock( 'ALLTOALL' )<br>
    !<br>
+  CALL plugin_clock()<br>
+  !<br>
    CALL mp_report()<br>
    !<br>
  END SUBROUTINE terminate_run<br>
<br>
===================================================================<br>
--- fromscra.f90    (revision 11775)<br>
+++ fromscra.f90    (working copy)<br>
@@ -116,6 +116,10 @@<br>
         CALL phbox( taub, iverbosity, eigrb )<br>
      END IF<br>
      !<br>
+    !     pass ions informations to plugins<br>
+    !<br>
+    CALL plugin_init_ions( tau0 )<br>
+    !<br>
      !     wfc initialization with random numbers<br>
      !<br>
      CALL wave_rand_init( cm_bgrp )<br>
<br>
===================================================================<br>
--- restart_sub.f90    (revision 11775)<br>
+++ restart_sub.f90    (working copy)<br>
@@ -155,6 +155,8 @@<br>
     IF ( tefield  ) CALL efield_berry_setup( eigr, tau0 )<br>
     IF ( tefield2 ) CALL efield_berry_setup2( eigr, tau0 )<br>
     !<br>
+   CALL plugin_init_ions( tau0 )<br>
+   !<br>
     edft%eself = eself<br>
     !<br>
     IF( tzerop .or. tzeroe .or. tzeroc ) THEN<br>
<br>
===================================================================<br>
--- input.f90    (revision 11775)<br>
+++ input.f90    (working copy)<br>
@@ -708,6 +708,11 @@<br>
        ! force pairing<br>
<br>
        force_pairing_ = force_pairing<br>
+<br>
+      ! ... having set all input keywords, read plugins' input file(s)<br>
+<br>
+      CALL plugin_read_input()<br>
+<br>
        !<br>
        ! ... the 'ATOMIC_SPECIES' card must be present, check it<br>
<br>
@@ -1133,6 +1138,8 @@<br>
        !<br>
        !   CALL sic_info()  ! maybe useful<br>
        !<br>
+      CALL plugin_print_info( )<br>
+      !<br>
        IF(tefield) call efield_info( )<br>
        IF(tefield2) call efield_info2( )<br>
<br>
===================================================================<br>
--- init_run.f90    (revision 11775)<br>
+++ init_run.f90    (working copy)<br>
@@ -117,6 +117,10 @@<br>
    !<br>
    CALL init_dimensions()<br>
    !<br>
+  ! ... initialization of plugin variables and arrays<br>
+  !<br>
+  CALL plugin_init_base()<br>
+  !<br>
    ! ... initialize atomic positions and cell<br>
    !<br>
    CALL init_geometry()<br>
<br>
===================================================================<br>
--- init.f90    (revision 11775)<br>
+++ init.f90    (working copy)<br>
@@ -419,5 +419,9 @@<br>
        !<br>
        call gcalb ( )<br>
        !<br>
+      !   pass new cell parameters to plugins<br>
+      !<br>
+      CALL plugin_init_cell( )<br>
+      !<br>
        return<br>
      end subroutine newinit_x<br>
<br>
===================================================================<br>
--- energies.f90    (revision 11775)<br>
+++ energies.f90    (working copy)<br>
@@ -197,6 +197,8 @@<br>
               IF( textfor ) WRITE( stdout, 16 ) eextfor<br>
            END IF<br>
            !<br>
+          CALL plugin_print_energies()<br>
+          !<br>
  1         FORMAT(6X,'                total energy = ',F18.10,' Hartree<br>
a.u.')<br>
  2         FORMAT(6X,'              kinetic energy = ',F18.10,' Hartree<br>
a.u.')<br>
  3         FORMAT(6X,'        electrostatic energy = ',F18.10,' Hartree<br>
a.u.')<br>
<br>
===================================================================<br>
--- stop_run.f90    (revision 11775)<br>
+++ stop_run.f90    (working copy)<br>
@@ -26,6 +26,8 @@<br>
    !<br>
    IF ( lconstrain ) CALL deallocate_constraint()<br>
    !<br>
+  CALL plugin_clean()<br>
+  !<br>
    CALL mp_global_end()<br>
    !<br>
  END SUBROUTINE stop_run<br>
<br>
===================================================================<br>
--- vofrho.f90    (revision 11775)<br>
+++ vofrho.f90    (working copy)<br>
@@ -62,6 +62,8 @@<br>
        USE tsvdw_module,     ONLY: EtsvdW,UtsvdW,FtsvdW,HtsvdW<br>
        USE mp_global,        ONLY: me_image<br>
        USE exx_module,       ONLY: dexx_dh, exxalfa  ! exx_wf related<br>
+<br>
+      USE plugin_variables, ONLY: plugin_etot<br>
<br>
        IMPLICIT NONE<br>
  !<br>
@@ -143,6 +145,16 @@<br>
  !<br>
        if (abivol.or.abisur) call vol_clu(rhor,rhog,sfac,nfi)<br>
        !<br>
+      !     compute plugin contributions to the potential, add it later<br>
+      !<br>
+      CALL plugin_get_potential(rhor,nfi)<br>
+      !<br>
+      !     compute plugin contributions to the energy<br>
+      !<br>
+      plugin_etot = 0.0_dp<br>
+      !<br>
+      CALL plugin_energy(rhor,plugin_etot)<br>
+      !<br>
        ttsic = ( ABS( self_interaction ) /= 0 )<br>
        !<br>
        IF( ttsic ) ALLOCATE( self_vloc( ngm ) )<br>
@@ -413,6 +425,11 @@<br>
          END IF<br>
          !<br>
        END IF<br>
+      !<br>
+      !     add plugin contributions to potential here...<br>
+      !<br>
+      CALL plugin_add_potential( rhor )<br>
+      !<br>
  !<br>
  !     rhor contains the xc potential in r-space<br>
  !<br>
@@ -499,6 +516,11 @@<br>
              fion = fion + fion1<br>
              !fion=fion+FtsvdW<br>
           END IF<br>
+         !<br>
+         !     plugin patches on internal forces<br>
+         !<br>
+         CALL plugin_int_forces(fion)<br>
+<br>
        END IF<br>
<br>
        DEALLOCATE( fion1 )<br>
@@ -625,6 +647,10 @@<br>
        if (abivol) etot = etot + P_ext*volclu<br>
        if (abisur) etot = etot + Surf_t*surfclu<br>
        !<br>
+      !     Add possible external contribution from plugins to the energy<br>
+      !<br>
+      etot = etot + plugin_etot<br>
+      !<br>
        IF( tpre ) THEN<br>
           !<br>
           detot6 = dekin6 + dh6 + dps6 + dsr6<br>
<br>
<br>
_______________________________________________<br>
Q-e-developers mailing list<br>
<a href="mailto:Q-e-developers@qe-forge.org">Q-e-developers@qe-forge.org</a><br>
<a href="http://qe-forge.org/mailman/listinfo/q-e-developers" rel="noreferrer" target="_blank">http://qe-forge.org/mailman/listinfo/q-e-developers</a><br>
</blockquote></div><br><br clear="all"><br>-- <br><div class="gmail_signature"><div dir="ltr"><div><div dir="ltr"><span><span><font color="#888888">Paolo Giannozzi, Dept. Chemistry&Physics&Environment,<br>
Univ. Udine, via delle Scienze 208, 33100 Udine, Italy<br>
Phone <a href="tel:%2B39-0432-558216" value="+390432558216" target="_blank">+39-0432-558216</a>, fax <a href="tel:%2B39-0432-558222" value="+390432558222" target="_blank">+39-0432-558222</a></font></span></span></div></div></div></div>
</div>