Fixed paw ortho bug when using processor groups. D1dB_SumAll was being used instead of D1dB_Vector_SumAll..EJB

This commit is contained in:
Eric Bylaska 2023-07-20 16:47:38 -07:00
parent cbfe8f42cb
commit 3355e82837
8 changed files with 200 additions and 37 deletions

View file

@ -4,11 +4,10 @@
*
* ***********************************
* * *
* * D1dB_SumAll *
* * *
* * *
* * D1dB_SumAll *
* * *
* ***********************************
subroutine D1dB_SumAll(sum)
c implicit none
real*8 sum
@ -191,11 +190,10 @@ c implicit none
* ***********************************
* * *
* * D1dB_Vector_SumAll *
* * *
* * *
* * D1dB_Vector_SumAll *
* * *
* ***********************************
subroutine D1dB_Vector_SumAll(n,sum)
c implicit none
integer n

View file

@ -3935,10 +3935,9 @@ c ****************************************
c ****************************************
c * *
c * Dneall_m_scale_s21 *
c * Dneall_m_scale_s21 *
c * *
c ****************************************
subroutine Dneall_m_scale_s21(mb,dte,s21)
implicit none
integer mb
@ -4008,9 +4007,10 @@ c ****************************************
end
c ****************************************
c * *
c * Dneall_m_scale_s11 *
c * Dneall_m_scale_s11 *
c * *
c ****************************************
@ -4507,7 +4507,7 @@ c ****************************************
real*8 ddot
external ddot
if (parallelized) then
if (mb.eq.0) then
nn = ne(1) + ne(2)
@ -4545,8 +4545,8 @@ c ****************************************
end do
!$OMP END MASTER
!$OMP BARRIER
call D1dB_SumAll(nn*nprj,dbl_mb(sw1b(1)))
call D1dB_SumAll(nn*nprj,dbl_mb(sw2b(1)))
call D1dB_Vector_SumAll(nn*nprj,dbl_mb(sw1b(1)))
call D1dB_Vector_SumAll(nn*nprj,dbl_mb(sw2b(1)))
if (mparallelized) then
call Parallel_taskid(taskid)
@ -4609,6 +4609,7 @@ c ****************************************
if (.not.value)
> call errquit('Dneall_m_add_sw1sw2:pop stack',1,MA_ERR)
else
if (mb.eq.0) then
ms1 = 1
ms2 = ispin

View file

@ -3,9 +3,9 @@
*
* ***********************************
* * *
* * control_read *
* * *
* * *
* * control_read *
* * *
* ***********************************
logical function control_read(code_in,rtdb)
implicit none

View file

@ -1460,9 +1460,9 @@ c end if
* ***************************
* * *
* * ion_chargeexist *
* * *
* * *
* * ion_chargeexist *
* * *
* ***************************
logical function ion_chargeexist()
implicit none
@ -1475,9 +1475,9 @@ c end if
end
* ***************************
* * *
* * ion_nion *
* * *
* * *
* * ion_nion *
* * *
* ***************************
integer function ion_nion()
implicit none

View file

@ -4,13 +4,13 @@
OBJ_OPTIMIZE = \
cpsd_input.o \
cpmd_input.o \
cpsdv5.o \
cpmdv5.o \
cpsdv5.o \
cpmdv5.o \
mmsdv1.o \
inner_loop_md.o \
mm_loop.o \
psi_lmbda.o psi_lmbda_sic.o psi_lmbda_paw.o psi_lmbda2.o psi_lmbda3.o \
psi_lmbda_omp.o \
inner_loop_md.o \
mm_loop.o \
psi_lmbda.o psi_lmbda_sic.o psi_lmbda_paw.o psi_lmbda2.o psi_lmbda3.o \
psi_lmbda_omp.o psi_lmbda_paw1.o \
fei_output.o \
dipole_motion.o \
pspw_et.o \
@ -23,9 +23,9 @@
LIBRARY = libpspw.a
USES_BLAS = psi_lmbda.f psi_lmbda2.f psi_lmbda3.f psi_lmbda_sic.f psi_lmbda_paw.f \
USES_BLAS = psi_lmbda.f psi_lmbda2.f psi_lmbda3.f psi_lmbda_sic.f psi_lmbda_paw.f psi_lmbda_paw1.f \
inner_loop.F inner_loop_md.F cpsdv5.F cpmdv5.F psi_lmbda2.f psi_lmbda_omp.F \
mmsdv1.F inner_loop_qmmm_step.F cpmd_qmmm_start.F cpmd_qmmm_stop.F mm_loop.F pspw_et.F
mmsdv1.F inner_loop_qmmm_step.F cpmd_qmmm_start.F cpmd_qmmm_stop.F mm_loop.F pspw_et.F
LIB_DEFINES =

View file

@ -670,13 +670,13 @@ c call Pack_cc_Sum(1,psi2(1,n),psi1(1,n),psi2(1,n))
call Dneall_f_Sortho(ms,psi2,psi_r,npack1)
end if
else
call psi_lmbda_paw(ispin,ne,nemax,npack1,psi_r,psi2,
> dte,
> lmd,dbl_mb(tmp_L(1)),ierr)
call psi_lmbda_paw(ispin,ne,nemax,npack1,psi_r,psi2,dte,
> lmd,dbl_mb(tmp_L(1)),ierr)
c call psi_lmbda_paw1(ispin,ne,nemax,npack1,psi_r,psi2,dte,
c > lmd,dbl_mb(tmp_L(1)),ierr)
end if
else if (fractional) then
call psi_lmbda2(ispin,ne,nemax,npack1,psi1,psi2,
> dte,occ1,
call psi_lmbda2(ispin,ne,nemax,npack1,psi1,psi2,dte,occ1,
> lmd,dbl_mb(tmp_L(1)),ierr)
!else if (sic.or.hfx) then
else if (sic) then

View file

@ -2,7 +2,7 @@ c $Id: psi_lmbda_paw.f 21177 2011-10-10 17:09:43Z bylaska $
* ***********************************
* * *
* * psi_lmbda_paw *
* * psi_lmbda_paw *
* * *
* ***********************************

View file

@ -0,0 +1,164 @@
*
* $Id$
*
* ********************************************
* * *
* * psi_lmbda_paw1 *
* * *
* ********************************************
subroutine psi_lmbda_paw1(ispin,ne,nemaxq,npack1,
> psi1,psi2,
> dte,
> lmbda,tmp,ierr)
implicit none
integer ispin,ne(2),nemaxq,npack1
complex*16 psi1(npack1,nemaxq)
complex*16 psi2(npack1,nemaxq)
real*8 dte
real*8 lmbda(*)
real*8 tmp(*)
integer ierr
integer MASTER
parameter (MASTER=0)
* **** local variables ****
logical notgram
integer taskid
integer ms,it
integer nn
integer s11,s12,s21,s22,st1,st2,sa1,sa0
real*8 adiff
* :::: iteration limit and tolerence for non-liner equations ::::
integer itrlmd
real*8 convg,convg2
parameter (itrlmd=120, convg=1.0d-15,convg2=1.0e-12)
* **** external functions ****
real*8 Dneall_m_dmax
external Dneall_m_dmax
integer Parallel_threadid
external Parallel_threadid
call nwpw_timing_start(3)
c call Dneall_m_size(1,nn)
c s11 = 0*nn + 1
c s21 = 1*nn + 1
c s22 = 2*nn + 1
c s12 = 3*nn + 1
c sa0 = 4*nn + 1
c sa1 = 5*nn + 1
c st1 = 6*nn + 1
c st2 = 7*nn + 1
c call dcopy(8*nn,0.0d0,0,tmp,1)
*:::::::::::::::::::::: Lagrangian multipliers ::::::::::::::::::::::
DO 640 ms=1,ispin
notgram=.true.
IF(ne(ms).le.0) GO TO 640
call Dneall_m_size(ms,nn)
s11 = 0*nn + 1
s21 = 1*nn + 1
s22 = 2*nn + 1
s12 = 3*nn + 1
sa0 = 4*nn + 1
sa1 = 5*nn + 1
st1 = 6*nn + 1
st2 = 7*nn + 1
ccccccccc Uncomment in this block and then run ch3cl-pspw-sd.nw cccccccccccc
call Dneall_ffm_sym_Multiply(ms,psi2,psi2,npack1,tmp(s22))
call Dneall_ffm_Multiply(ms,psi2,psi1,npack1,tmp(s21))
call Dneall_ffm_Multiply(ms,psi1,psi2,npack1,tmp(s12))
call Dneall_ffm_sym_Multiply(ms,psi1,psi1,npack1,tmp(s11))
!call Dne_ffm_combo_sym_Multiply(ms,psi1,psi2,npack1,tmp(s11))
* *** add the extra paw overlap part ***
call psp_add_paw_extra_overlap1(ms,psi2,tmp(s22))
call psp_add_paw_extra_overlap2(ms,psi2,psi1,tmp(s21))
call psp_add_paw_extra_overlap2(ms,psi1,psi2,tmp(s12))
call psp_add_paw_extra_overlap1(ms,psi1,tmp(s11))
ccccccccc Uncomment in this block and then run ch3cl-pspw-sd.nw cccccccccccc
* ***** scale the overlap matrices ****
call Dneall_m_scale_s22(ms,dte,tmp(s22))
call Dneall_m_scale_s21(ms,dte,tmp(s21))
call Dneall_m_scale_s21(ms,dte,tmp(s12))
call Dneall_m_scale_s11(ms,dte,tmp(s11))
!call dcopy(nn,tmp(s22),1,tmp(sa0),1)
call Parallel_shared_vector_copy(.true.,nn,tmp(s22),tmp(sa0))
!call Parallel_shared_vector_copy(.true.,nn,tmp(s21),tmp(s12))
do it=1,itrlmd
!call dcopy(nn,tmp(s22),1,tmp(sa1),1)
call Parallel_shared_vector_copy(.true.,nn,tmp(s22),tmp(sa1))
call Dneall_mmm_Multiply0(ms,
> tmp(s21),tmp(sa0),1.0d0,
> tmp(sa1),1.0d0)
call Dneall_mmm_Multiply0(ms,
> tmp(sa0),tmp(s12),1.0d0,
> tmp(sa1),1.0d0)
call Dneall_mmm_Multiply0(ms,
> tmp(s11),tmp(sa0),1.0d0,
> tmp(st1),0.0d0)
call Dneall_mmm_Multiply(ms,
> tmp(sa0),tmp(st1),1.0d0,
> tmp(sa1),1.0d0)
!call dcopy(nn,tmp(sa1),1,tmp(st1),1)
call Parallel_shared_vector_copy(.true.,nn,tmp(sa1),tmp(st1))
call daxpy_omp(nn,(-1.0d0),tmp(sa0),1,tmp(st1),1)
adiff = Dneall_m_dmax(ms,tmp(st1))
if(adiff.lt.convg) GO TO 630
if (adiff.gt.1.0d10) go to 620
!call dcopy(nn,tmp(sa1),1,tmp(sa0),1)
call Parallel_shared_vector_copy(.true.,nn,tmp(sa1),tmp(sa0))
end do
if (adiff.lt.convg2) GO TO 630
620 continue
ierr=10
call Parallel_taskid(taskid)
if (taskid.eq.MASTER) then
WRITE(6,*)
> 'Warning: Lagrange Multiplier tolerance too high:',adiff
WRITE(6,*) ' +Try using a smaller time step'
WRITE(6,*) ' +Gram-Schmidt being performed, spin:',ms
end if
c call Dneall_f_ortho(ms,psi2,npack1)
call Dneall_f_GramSchmidt(ms,psi2,npack1)
notgram = .false.
630 continue
* ::::::::::::::::: correction due to the constraint :::::::::::::::::
if (notgram) then
call Dneall_fmf_Multiply(ms,
> psi1,npack1,
> tmp(sa1), dte,
> psi2,1.0d0)
end if
call Dneall_mm_Expand(ms,tmp(sa1),lmbda)
640 continue
c*::::::::::::::::: correction due to the constraint :::::::::::::::::
c call Dneall_fmf_Multiply(0,
c > psi1,npack1,
c > lmbda, dte,
c > psi2,1.0d0)
call nwpw_timing_end(3)
return
end