From 38cf4205177fc8362b400f064dbc0c614910a6c1 Mon Sep 17 00:00:00 2001
From: Daniele Rapetti <daniele.rapetti@sissa.it>
Date: Thu, 6 Jul 2023 16:21:50 +0200
Subject: [PATCH] pw now should give the correct masses to plumed

adding a comment to explain the new piece of code

corrected F90 error and patched the wanted version of qe
---
 .../qespresso-6.2.diff/PW/src/plugin_ext_forces.f90 | 13 +++++++++----
 .../qespresso-7.0.diff/PW/src/plugin_ext_forces.f90 | 13 ++++++++++---
 2 files changed, 19 insertions(+), 7 deletions(-)

diff --git a/patches/qespresso-6.2.diff/PW/src/plugin_ext_forces.f90 b/patches/qespresso-6.2.diff/PW/src/plugin_ext_forces.f90
index ac9925473c..4df0bbaab7 100644
--- a/patches/qespresso-6.2.diff/PW/src/plugin_ext_forces.f90
+++ b/patches/qespresso-6.2.diff/PW/src/plugin_ext_forces.f90
@@ -18,19 +18,21 @@ SUBROUTINE plugin_ext_forces()
   USE plugin_flags
   !
   USE cell_base,        ONLY : alat, at
-  USE ions_base,        ONLY : tau, nat,amass
+  USE ions_base,        ONLY : tau, nat, amass, ityp
   USE force_mod,        ONLY : force,sigma
   USE control_flags,    ONLY : istep
   USE ener,             ONLY : etot 
   !
   IMPLICIT NONE
   !
-  INTEGER:: i,j
+  INTEGER:: i,j,ia
   REAL(DP) :: at_plumed(3,3)
   REAL(DP) :: virial(3,3)
   REAL(DP) :: volume
   REAL(DP), ALLOCATABLE :: tau_plumed(:,:)
+  REAL(DP) :: masses_plumed(nat)
   !
+  masses_plumed = 0.0_DP
   IF(use_plumed) then
     IF(ionode)THEN
       at_plumed=alat*at;  ! the cell, rescaled properly
@@ -43,9 +45,12 @@ SUBROUTINE plugin_ext_forces()
              -at_plumed(1,2)*at_plumed(3,3)*at_plumed(2,1) &
              -at_plumed(1,3)*at_plumed(3,1)*at_plumed(2,2) 
       virial=-sigma*volume
-
+      ! the masses in QE are stored per type, see q-e//Modules/ions_base.f90
+      do ia=1,nat
+        masses_plumed(ia)=amass(ityp(ia))
+      end do
       CALL plumed_f_gcmd("setStep"//char(0),istep)
-      CALL plumed_f_gcmd("setMasses"//char(0),amass)
+      CALL plumed_f_gcmd("setMasses"//char(0),masses_plumed)
       CALL plumed_f_gcmd("setForces"//char(0),force)
       CALL plumed_f_gcmd("setPositions"//char(0),tau_plumed)
       CALL plumed_f_gcmd("setBox"//char(0),at_plumed)
diff --git a/patches/qespresso-7.0.diff/PW/src/plugin_ext_forces.f90 b/patches/qespresso-7.0.diff/PW/src/plugin_ext_forces.f90
index 52d35bc5aa..5d4762fb42 100644
--- a/patches/qespresso-7.0.diff/PW/src/plugin_ext_forces.f90
+++ b/patches/qespresso-7.0.diff/PW/src/plugin_ext_forces.f90
@@ -18,7 +18,7 @@ SUBROUTINE plugin_ext_forces()
   USE plugin_flags,     ONLY : use_plumed
   !
   USE cell_base,        ONLY : alat, at
-  USE ions_base,        ONLY : tau, nat,amass
+  USE ions_base,        ONLY : tau, nat, amass, ityp
   USE force_mod,        ONLY : force,sigma
   USE control_flags,    ONLY : istep
   USE ener,             ONLY : etot 
@@ -27,12 +27,14 @@ SUBROUTINE plugin_ext_forces()
   !
   IMPLICIT NONE
   !
-  INTEGER:: i,j
+  INTEGER:: i,j,ia
   REAL(DP) :: at_plumed(3,3)
   REAL(DP) :: virial(3,3)
   REAL(DP) :: volume
   REAL(DP), ALLOCATABLE :: tau_plumed(:,:)
+  REAL(DP) :: masses_plumed(nat)
   !
+  masses_plumed = 0.0_DP
   IF(use_plumed) then
     IF(ionode)THEN
       at_plumed=alat*at;  ! the cell, rescaled properly
@@ -46,8 +48,13 @@ SUBROUTINE plugin_ext_forces()
              -at_plumed(1,3)*at_plumed(3,1)*at_plumed(2,2) 
       virial=-sigma*volume
 
+      ! the masses in QE are stored per type, see q-e//Modules/ions_base.f90
+      do ia=1,nat
+        masses_plumed(ia)=amass(ityp(ia))
+      end do
+
       CALL plumed_f_gcmd("setStep"//char(0),istep)
-      CALL plumed_f_gcmd("setMasses"//char(0),amass)
+      CALL plumed_f_gcmd("setMasses"//char(0),masses_plumed)
       CALL plumed_f_gcmd("setForces"//char(0),force)
       CALL plumed_f_gcmd("setPositions"//char(0),tau_plumed)
       CALL plumed_f_gcmd("setBox"//char(0),at_plumed)