phase0_2014.03.01.patch
phase0_2014.03.01.patch
—
differences between files,
34 KB (34884 bytes)
ファイルコンテンツ
diff -uprN phase0_2014.03/src_phase/Initialization.F90 phase0_2014.03.01/src_phase/Initialization.F90
--- phase0_2014.03/src_phase/Initialization.F90 2014-11-30 23:11:10.000000000 +0900
+++ phase0_2014.03.01/src_phase/Initialization.F90 2015-06-03 09:36:31.813693486 +0900
@@ -160,7 +160,9 @@ contains
subroutine aavers
include 'version.h' ! svn_revision
character(len=72) :: vers, system, codename
- write(vers,'("Revision:",i5," -- ORG_Parallel --")') svn_revision
+! write(vers,'("Revision:",i5," -- ORG_Parallel --")') svn_revision
+ write(vers,'("phase/0 2014.03.01 Revision:",i5," -- ORG_Parallel --")') svn_revision
+
codename = 'phaseUnif'
system = ''
diff -uprN phase0_2014.03/src_phase/NEB.F90 phase0_2014.03.01/src_phase/NEB.F90
--- phase0_2014.03/src_phase/NEB.F90 2014-11-29 14:33:06.000000000 +0900
+++ phase0_2014.03.01/src_phase/NEB.F90 2015-06-03 09:35:09.548581471 +0900
@@ -1250,6 +1250,8 @@ subroutine neb_force
real(8), allocatable :: fx0(:), fy0(:), fz0(:)
real(8), allocatable :: fx1(:), fy1(:), fz1(:)
+ integer :: ci_target=-1
+
if (mype==0) write(nfneb,*) 'neb force...'
allocate( k_spring(neb%number_of_images-1) )
@@ -1261,6 +1263,17 @@ subroutine neb_force
end do
max_energy = maxval( neb%image(:)%energy )
+ if( ( neb%step >= neb%cond%ci_neb_start) .and. ( neb%cond%climbing_image ) ) then
+ do k=2,neb%number_of_images-1
+ if (neb%image(k)%energy==max_energy) then
+ ci_target=k
+ write(nfneb,*) 'climing image : ',ci_target
+ exit
+ endif
+ enddo
+ endif
+
+
do k = 2, neb%number_of_images-1
@@ -1335,7 +1348,7 @@ subroutine neb_force
fz_spring(:) = k_spring0(image_id) * ( rz2(:) - rz0(:) ) &
+ k_spring0(image_id - 1 ) * ( rz1(:) - rz0(:) )
- if( ( neb%step < neb%cond%ci_neb_start ) .or. ( .not. neb%cond%climbing_image) ) then
+ if(k/=ci_target)then
! get true force perpendicular to local tangent
@@ -1459,7 +1472,7 @@ subroutine neb_force
write(nfneb,'(a31,3e11.3)') ' max force : ', &
maxval( abs( fx1 ) ), maxval( abs( fy1 ) ), maxval( abs( fz1 ) )
endif
- elseif( ( neb%step >= neb%cond%ci_neb_start) .and. ( neb%cond%climbing_image ) ) then
+ else
! get true force along local tangent
@@ -1965,6 +1978,7 @@ subroutine write_result(inc_energy)
open(nfdynm,file=F_DYNM)
open(nfenf,file=F_ENF)
call m_IS_wd_speciesname_etc(nfdynm)
+ call m_IS_wd_speciesname_etc(nfnebdynm)
write(nfenf,*) '#step image image_distance energy force_org force_neb force_normal'
sum_d = 0.0d0
do i=1, neb%number_of_images
@@ -1988,14 +2002,21 @@ subroutine write_result(inc_energy)
maxval(abs(neb%image(i)%force(:,:))), &
maxval(abs(neb%image(i)%true_force(:,:)))
endif
+ write(nfnebdynm,'(a,i5,i8,a)') ' cps and forc at (nebstep, iter_total=',neb%step,iteration,')'
write(nfdynm,'(a,i5,i8,a)') ' cps and forc at (nebstep, iter_total=',neb%step,iteration,')'
do j = 1, neb%image(i)%num_atom
- write(nfnebdynm,'(3i5,3f20.10)') neb%step, i, j, neb%image(i)%cps(j,1:3)
+! write(nfnebdynm,'(3i5,3f20.10)') neb%step, i, j, neb%image(i)%cps(j,1:3)
+ write(nfnebdynm,'(" ",i4,3f15.9,3f12.6)') j &
+ &, neb%image(i)%cps(j,1),neb%image(i)%cps(j,2),neb%image(i)%cps(j,3) &
+ &, neb%image(i)%force(j,1), neb%image(i)%force(j,2), neb%image(i)%force(j,3)
write(nfdynm,'(" ",i4,3f15.9,3f12.6)') j &
&, neb%image(i)%cps(j,1),neb%image(i)%cps(j,2),neb%image(i)%cps(j,3) &
&, neb%image(i)%force(j,1), neb%image(i)%force(j,2), neb%image(i)%force(j,3)
end do
end do
+
+ if(inc_energy) write(nfnebenf,*)
+
close(nfdynm)
close(nfenf)
! temporary
diff -uprN phase0_2014.03/src_phase/Preparation.F90 phase0_2014.03.01/src_phase/Preparation.F90
--- phase0_2014.03/src_phase/Preparation.F90 2014-12-01 14:04:42.000000000 +0900
+++ phase0_2014.03.01/src_phase/Preparation.F90 2015-06-03 09:35:09.548581471 +0900
@@ -151,7 +151,7 @@ subroutine Preparation()
use m_Dipole, only:m_Dipole_calc
use m_Wannier, only:m_Wan_gen_weight
use m_FiniteElectricField, only:m_FEF_init
- use m_IterationNumbers, only : m_Iter_rd_iteration_numbers,iteration_unit_cell
+ use m_IterationNumbers, only : m_Iter_rd_iteration_numbers,iteration_unit_cell,m_Iter_cmix_reset
#ifdef _INCLUDE_EXX_
use m_ES_ExactExchange, only: m_ES_EXX_init
#endif
@@ -561,6 +561,8 @@ subroutine Preparation()
if ( sw_phonon == OFF ) call m_IS_gnrt_supercell_symmetry(paramset,nfout)
+ call m_Iter_cmix_reset()
+
contains
subroutine fft_box_finding_way(outer_or_inner)
integer,intent(out):: outer_or_inner
diff -uprN phase0_2014.03/src_phase/m_Charge_Density.F90 phase0_2014.03.01/src_phase/m_Charge_Density.F90
--- phase0_2014.03/src_phase/m_Charge_Density.F90 2014-12-01 20:34:03.000000000 +0900
+++ phase0_2014.03.01/src_phase/m_Charge_Density.F90 2015-06-03 09:35:09.549581969 +0900
@@ -5694,8 +5694,8 @@ contains
end subroutine m_CD_rspace_charge
- subroutine m_CD_get_rspace_charge(nfout,na,nb,nc,rho)
- integer,intent(in) :: na,nb,nc,nfout
+ subroutine m_CD_get_rspace_charge(nfout,na,nb,nc,rho,is)
+ integer,intent(in) :: na,nb,nc,nfout,is
real(kind=8), dimension(na,nb,nc), intent(out) :: rho
integer,dimension(3) :: boxsize
real(kind=DP),dimension(3,3) :: cellsize
@@ -5703,12 +5703,12 @@ contains
integer :: iloop
call m_CD_alloc_rspace_charge()
- do iloop=1,nspin
- call map_valence_charge_to_fft_box(iloop)
+ !do iloop=1,nspin
+ call map_valence_charge_to_fft_box(is)
call m_FFT_CD_inverse0(nfout,afft)
!!$ call m_FFT_CD0(nfout,afft,INVERSE)
call chgr()
- enddo
+ !enddo
call m_CD_dealloc_rspace_charge()
contains
diff -uprN phase0_2014.03/src_phase/m_ES_ortho.F90 phase0_2014.03.01/src_phase/m_ES_ortho.F90
--- phase0_2014.03/src_phase/m_ES_ortho.F90 2014-11-29 14:33:03.000000000 +0900
+++ phase0_2014.03.01/src_phase/m_ES_ortho.F90 2015-06-03 09:36:08.291938506 +0900
@@ -1554,11 +1554,19 @@ contains
!diagonal
do i1 = i, NB_END
if(mode == ORTHONORMALIZATION .or. mode == NORMALIZATION) then
- call WSW_t_g(ik,i1,mod_pot,fr,psi_t,np_g1k_x,neg,kimg,bpr_t,bpi_t)
+ if(mod_pot == VANDERBILT_TYPE) then
+ call WSW_t_g(ik,i1,mod_pot,fr,psi_t,np_g1k_x,neg,kimg,bpr_t,bpi_t)
+ else
+ call WSW_t_g(ik,i1,mod_pot,fr,psi_t,np_g1k_x,neg,kimg)
+ endif
!!$ ! fr = 1/dsqrt(<Psi(i)|S|Psi(i)>)
- if(dabs(fr-1.d0) > DELTA)&
- & call normalize_bp_and_psi_t_g(ik,i1,fr,mod_pot &
- & , psi_t,np_g1k_x,neg,kimg,bpr_t,bpi_t)
+ if(dabs(fr-1.d0) > DELTA .and. mod_pot == VANDERBILT_TYPE )then
+ call normalize_bp_and_psi_t_g(ik,i1,fr,mod_pot &
+ , psi_t,np_g1k_x,neg,kimg,bpr_t,bpi_t)
+ else if (dabs(fr-1.d0) > DELTA ) then
+ call normalize_bp_and_psi_t_g(ik,i1,fr,mod_pot &
+ , psi_t,np_g1k_x,neg,kimg)
+ endif
! |Psi(i)> = |Psi(i)> * fr, <beta|Psi(i)> = <beta|Psi(i)> * fr
if(mod_pot == VANDERBILT_TYPE) call cp_bpr2bprtw(i1,neg) ! bpr_t -> bpr_tw1, bpr_tw2
end if
@@ -1566,11 +1574,15 @@ contains
if(i1 == neg) cycle
call cp_psi2psii_g(ik,i1) ! psi_t(:,i,:) -> psi_ir,psi_ii
call W1SW2_t_r_g(ik,i1,NB_END,mod_pot,psi_t,np_g1k_x,neg,kimg) ! -> p1Sp2
- if(mod_pot == VANDERBILT_TYPE) &
- & call cp_bpr2bpi_g(kimg_t_wk,i1,bpr_t,bpi_t) ! bpr_t(:,i1)->bp_ir,bp_ii
- call modify_bp_and_psi_t_r_g(ik,i1,NB_END,mod_pot &
+ if(mod_pot == VANDERBILT_TYPE) then
+ call cp_bpr2bpi_g(kimg_t_wk,i1,bpr_t,bpi_t) ! bpr_t(:,i1)->bp_ir,bp_ii
+ call modify_bp_and_psi_t_r_g(ik,i1,NB_END,mod_pot &
& ,psi_t,np_g1k_x,neg,kimg,bpr_t,bpi_t)
! psi_t, bpr_t, pbi_t, p1Sp2 -> psi_t, bpr_t, bpi_t
+ else
+ call modify_bp_and_psi_t_r_g(ik,i1,NB_END,mod_pot &
+ & ,psi_t,np_g1k_x,neg,kimg)
+ endif
end if
end do ! i1-loop
if(mode /= NORMALIZATION) then
@@ -1598,9 +1610,14 @@ contains
end if
if(nrank_e > 1) call mpi_allreduce(MPI_IN_PLACE, p1Sp2_NB, NB*NB*kimg_t_wk &
& ,mpi_double_precision,mpi_sum,mpi_k_world(myrank_k),ierr)
- call modify_bp_and_psi_t_r_blk_g(ik,i,i2,p1Sp2_NB,NB_END,NB_END2,kimg_t_wk &
- & , mod_pot,psi_t,psi_t_dia,np_g1k_x,neg,kimg,bpr_t,bpi_t,bpr_t_dia,bpi_t_dia)
+ if(mod_pot == VANDERBILT_TYPE) then
+ call modify_bp_and_psi_t_r_blk_g(ik,i,i2,p1Sp2_NB,NB_END,NB_END2,kimg_t_wk &
+ & , mod_pot,psi_t,psi_t_dia,np_g1k_x,neg,kimg,bpr_t,bpi_t,bpr_t_dia,bpi_t_dia)
! psi_t, bpr_t, pbi_t, p1Sp2 -> psi_t, bpr_t, bpi_t
+ else
+ call modify_bp_and_psi_t_r_blk_g(ik,i,i2,p1Sp2_NB,NB_END,NB_END2,kimg_t_wk &
+ & , mod_pot,psi_t,psi_t_dia,np_g1k_x,neg,kimg)
+ endif
end if
end do ! i2-loop
end if ! if(mode /= NORMALIZATION)
@@ -6621,8 +6638,13 @@ contains
if(mod_pot == VANDERBILT_TYPE) &
& call cp_bpr2bpi_g(kimg_t_wk,i1,bpr_t,bpi_t) ! -> bp_ir, bp_ii
call W1SW2_t_r_g(ik,i1,NB_END,mod_pot,phi_t,np_g1k_x,neg,kimg) ! -> p1Sp2
- call modify_bp_and_psi_t_r_g(ik,i1,NB_END,mod_pot &
- & ,phi_t,np_g1k_x,neg,kimg,phifr_t,phifi_t) ! psi_t, bpr_t, pbi_t, p1Sp2 -> psi_t, bpr_t, bpi_t
+ if(mod_pot == VANDERBILT_TYPE) then
+ call modify_bp_and_psi_t_r_g(ik,i1,NB_END,mod_pot &
+ & ,phi_t,np_g1k_x,neg,kimg,phifr_t,phifi_t) ! psi_t, bpr_t, pbi_t, p1Sp2 -> psi_t, bpr_t, bpi_t
+ else
+ call modify_bp_and_psi_t_r_g(ik,i1,NB_END,mod_pot &
+ & ,phi_t,np_g1k_x,neg,kimg) ! psi_t, bpr_t, pbi_t, p1Sp2 -> psi_t, bpr_t, bpi_t
+ endif
end do ! i1-loop
if (sw_timing_2ndlevel == ON) call tstatc0_end(id_sname2)
!lower
@@ -6654,9 +6676,14 @@ contains
! ==============================================================================
if(nrank_e > 1) call mpi_allreduce(MPI_IN_PLACE, p1Sp2_NB, NB*NB*kimg_t_wk &
& ,mpi_double_precision,mpi_sum,mpi_k_world(myrank_k),ierr)
- call modify_bp_and_psi_t_r_blk_g(ik,i,i2,p1Sp2_NB,NB_END,NB_END2,kimg_t_wk &
- & , mod_pot,phi_t,psi_t_dia,np_g1k_x,neg,kimg,phifr_t,phifi_t,bpr_t_dia,bpi_t_dia)
+ if(mod_pot==VANDERBILT_TYPE)then
+ call modify_bp_and_psi_t_r_blk_g(ik,i,i2,p1Sp2_NB,NB_END,NB_END2,kimg_t_wk &
+ & , mod_pot,phi_t,psi_t_dia,np_g1k_x,neg,kimg,phifr_t,phifi_t,bpr_t_dia,bpi_t_dia)
! phi_t, phifr_t, phifi_t, p1Sp2 -> phi_t, phifr_t, phifi_t
+ else
+ call modify_bp_and_psi_t_r_blk_g(ik,i,i2,p1Sp2_NB,NB_END,NB_END2,kimg_t_wk &
+ & , mod_pot,phi_t,psi_t_dia,np_g1k_x,neg,kimg)
+ endif
end do ! i2-loop
if (sw_timing_2ndlevel == ON) call tstatc0_end(id_sname3)
end do ! i-loop
diff -uprN phase0_2014.03/src_phase/m_Kpoints.F90 phase0_2014.03.01/src_phase/m_Kpoints.F90
--- phase0_2014.03/src_phase/m_Kpoints.F90 2014-11-30 22:32:16.000000000 +0900
+++ phase0_2014.03.01/src_phase/m_Kpoints.F90 2015-06-03 09:35:09.563588967 +0900
@@ -882,8 +882,10 @@ contains
write(nfout,'(" << trmat >>")')
write(nfout,'(3f8.4)') ((trmat(i,j),j=1,3),i=1,3)
end if
-
- call gen_kpoint_in_full_BZ( preallocation, kv3_fbz, vkxyz_fbz )
+! === modified by T. Y. 2015/01/16 ===
+!!$ call gen_kpoint_in_full_BZ( preallocation, kv3_fbz, vkxyz_fbz )
+ call gen_kpoint_in_full_BZ( preallocation)
+! ====================================
! ========================================================================= 12.0A
else if(way_ksample == MESH) then
@@ -2853,10 +2855,13 @@ contains
end subroutine m_Kp_dealloc
! === KT_add ==== 2014/09/30
- subroutine gen_kpoint_in_full_BZ( paramset, kv3_fbz, vkxyz_fbz )
+! === modified by T. Y. 2015/01/16 ===
+!!$ subroutine gen_kpoint_in_full_BZ( paramset, kv3_fbz, vkxyz_fbz )
+ subroutine gen_kpoint_in_full_BZ( paramset)
logical, intent(in) :: paramset
- integer, intent(in) :: kv3_fbz
- real(kind=DP), intent(inout) :: vkxyz_fbz(kv3_fbz,3,CARTS)
+!!$ integer, intent(in) :: kv3_fbz
+!!$ real(kind=DP), intent(inout) :: vkxyz_fbz(kv3_fbz,3,CARTS)
+! ====================================
integer :: i, ik
real(kind=DP) :: b1(3), b2(3), b3(3)
@@ -2898,6 +2903,10 @@ contains
real(kind=DP) :: mtr(3,3)
real(kind=DP), parameter :: lambda = 1.d0+1.d-8
real(kind=DP) :: gg, gg_ref, ko(3)
+! === modified by T. Y. 2015/01/16 ===
+ integer,parameter :: max_icount = 10000
+ integer :: icount
+! ====================================
do i=1,3
do j=1,3
@@ -2909,6 +2918,7 @@ contains
ko(1:3) = vkxyz_fbz(ik,1:3,BUCS)
do i=1,nface
gg_ref = dot_product(iface_bz(:,i),matmul(mtr,iface_bz(:,i))) * 0.5d0
+ icount = 0
do
gg = dot_product(iface_bz(:,i),matmul(mtr,ko))
if(gg > gg_ref*lambda) then
@@ -2916,6 +2926,14 @@ contains
else
exit
end if
+! === modified by T. Y. 2015/01/16 ===
+ icount = icount+1
+ if(icount > max_icount) then
+ write(nfout,'(" icount > max_icount in <<move_kpt_to_inside_fbz>>")')
+ write(nfout,'(" icount = ",i8)') icount
+ stop ' stop at <<move_kpt_to_inside_fbz>>. Icount is too large'
+ end if
+! ====================================
end do
end do
vkxyz_fbz(ik,1:3,BUCS) = ko(1:3)
diff -uprN phase0_2014.03/src_phase/m_Total_Energy.F90 phase0_2014.03.01/src_phase/m_Total_Energy.F90
--- phase0_2014.03/src_phase/m_Total_Energy.F90 2014-11-29 14:33:05.000000000 +0900
+++ phase0_2014.03.01/src_phase/m_Total_Energy.F90 2015-06-03 09:35:09.563588967 +0900
@@ -1925,6 +1925,13 @@ contains
etotal0 = ekinet+ehartr+exc+elocal+enonlc+eewald-epc
endif
! ================= 13.0S
+#ifdef __FUJITSU
+ call mpi_bcast(etotal0,1,mpi_double_precision,0,mpi_comm_group,ierr)
+#else
+#ifdef __ETOTAL_UNIFY_
+ call mpi_bcast(etotal0,1,mpi_double_precision,0,mpi_comm_group,ierr)
+#endif
+#endif
if(sw_dipole_correction == ON) etotal0 = etotal0 + edip
if(sw_hubbard == ON) etotal0 = etotal0 + ehub0
diff -uprN phase0_2014.03/src_phase/m_XC_Potential.F90 phase0_2014.03.01/src_phase/m_XC_Potential.F90
--- phase0_2014.03/src_phase/m_XC_Potential.F90 2014-11-30 22:32:17.000000000 +0900
+++ phase0_2014.03.01/src_phase/m_XC_Potential.F90 2015-06-03 09:35:09.564589467 +0900
@@ -1917,11 +1917,18 @@ contains
allocate(dfddrho_vdw(nfftx,nffty,nfftz));dfddrho_vdw=0.d0
allocate(chgrhr_red(nfftcd));chgrhr_red=0.0d0
allocate(grad_rho_red(nfftcd));grad_rho_red=0.d0
- do i=ista_fftph,iend_fftph
- chgrhr_red(i) = chgrhr_l(i,ispin)
- grad_rho_red(i) = grad_rho(i,ispin)
- enddo
- call mpi_allreduce(MPI_IN_PLACE,chgrhr_red,nfftcd*ispin,mpi_double_precision,mpi_sum,&
+ if(ispin==1)then
+ do i=ista_fftph,iend_fftph
+ chgrhr_red(i) = chgrhr_l(i,1)
+ grad_rho_red(i) = grad_rho(i,1)
+ enddo
+ else
+ do i=ista_fftph,iend_fftph
+ chgrhr_red(i) = chgrhr_l(i,1)+chgrhr_l(i,2)
+ grad_rho_red(i) = grad_rho(i,1)+grad_rho(i,2)
+ enddo
+ endif
+ call mpi_allreduce(MPI_IN_PLACE,chgrhr_red,nfftcd,mpi_double_precision,mpi_sum,&
& mpi_cdfft_world(myrank_ggacmp),ierr)
call mpi_allreduce(MPI_IN_PLACE,grad_rho_red,nfftcd,mpi_double_precision,mpi_sum,&
& mpi_cdfft_world(myrank_ggacmp),ierr)
@@ -1931,8 +1938,8 @@ contains
do iz=1,nfftz
ixyz = (iz-1)*nffty*nfftx+(iy-1)*nfftx+ix
if(ixyz<ista_fftph.or.ixyz>iend_fftph) cycle
- dF_drho(ixyz,ispin) = dF_drho(ixyz,ispin)+dfdrho_vdw(ix,iy,iz)
- dF_dgradrho(ixyz,ispin) = dF_dgradrho(ixyz,ispin)+dfddrho_vdw(ix,iy,iz)
+ dF_drho(ixyz,1:ispin) = dF_drho(ixyz,1:ispin)+dfdrho_vdw(ix,iy,iz)
+ dF_dgradrho(ixyz,1:ispin) = dF_dgradrho(ixyz,1:ispin)+dfddrho_vdw(ix,iy,iz)
enddo
enddo
enddo
diff -uprN phase0_2014.03/src_phase/m_vdWDF.F90 phase0_2014.03.01/src_phase/m_vdWDF.F90
--- phase0_2014.03/src_phase/m_vdWDF.F90 2014-11-29 14:33:06.000000000 +0900
+++ phase0_2014.03.01/src_phase/m_vdWDF.F90 2015-06-03 09:35:09.580597466 +0900
@@ -295,14 +295,14 @@ module m_vdWDF
firstcall = .false.
end subroutine initialize_vdwdf_scf
- subroutine initialize_vdwdf_oneshot()
+ subroutine initialize_vdwdf_oneshot(is)
use progress_bar, only : set_printable
+ integer, intent(in) :: is
integer :: i,cix,ciy,ciz,nrxyz
real(kind=DP) :: q
call set_printable(printable)
if(printable) write(nfout,'(a)') 'initialization ...'
-
call do_cell_params()
rinplw = 1.d0/(dble(na*nb*nc))
@@ -312,7 +312,7 @@ module m_vdWDF
maxk = dble(nr12)/r12max
call alloc_vdw()
- call m_CD_get_rspace_charge(nfout,na,nb,nc,rho)
+ call m_CD_get_rspace_charge(nfout,na,nb,nc,rho,is)
call print_vdw_parameters()
call derivation(na,nb,nc,aa,rho,dv,grad)
diff -uprN phase0_2014.03/src_phase/vdW.F90 phase0_2014.03.01/src_phase/vdW.F90
--- phase0_2014.03/src_phase/vdW.F90 2014-11-29 14:33:05.000000000 +0900
+++ phase0_2014.03.01/src_phase/vdW.F90 2015-06-03 09:35:09.592603467 +0900
@@ -56,30 +56,36 @@ end subroutine vdW_scf
subroutine vdW_oneshot()
use m_Const_Parameters, only : DP
- use m_Control_Parameters, only : printable
+ use m_Control_Parameters, only : printable,nspin
use m_Files, only : nfout
use m_Parallelization, only : mype
use m_Timing, only : tstatc0_begin,tstatc0_end
use m_Total_Energy, only : etotal
use m_vdWDF, only: initialize_vdwdf_oneshot,build_theta,vdWdf_core,corrections,Ecnl_12,Ecnl_3,Ecnl_3s,finalize_vdwdf
implicit none
+ integer :: is
integer :: id_sname = -1
real(kind=DP) :: Ecnl
call tstatc0_begin('vdW_oneshot ',id_sname,1)
if(mype==0) then
write(nfout,'(a)') "** 'oneshot' calculation of the vdW-interaction **"
endif
- call initialize_vdwdf_oneshot()
+ Ecnl=0.d0
+ is=1
+ if(nspin>1) is=-1
+
+ call initialize_vdwdf_oneshot(is)
call build_theta()
call vdWdf_core()
+ Ecnl = Ecnl+Ecnl_12
call corrections()
- Ecnl = (Ecnl_12 + Ecnl_3 - Ecnl_3s)
+ Ecnl = Ecnl+Ecnl_3 - Ecnl_3s
+ call finalize_vdwdf()
etotal = etotal + Ecnl
if(printable) then
write(nfout,'(a,f20.10,a)') 'vdW energy : ',Ecnl, ' hartree'
write(nfout,'(a,f20.10,a)') '--> total energy : ',etotal,' hartree'
endif
- call finalize_vdwdf()
call tstatc0_end(id_sname)
end subroutine vdW_oneshot
#endif
diff -uprN phase0_2014.03/src_phase_3d/Initialization.F90 phase0_2014.03.01/src_phase_3d/Initialization.F90
--- phase0_2014.03/src_phase_3d/Initialization.F90 2014-11-30 23:16:06.000000000 +0900
+++ phase0_2014.03.01/src_phase_3d/Initialization.F90 2015-06-03 09:37:52.228877999 +0900
@@ -160,7 +160,8 @@ contains
subroutine aavers
include 'version.h' ! svn_revision
character(len=72) :: vers, system, codename
- write(vers,'("Revision:",i5, " --- 3D_Parallel --")') svn_revision
+! write(vers,'("Revision:",i5, " --- 3D_Parallel --")') svn_revision
+ write(vers,'("phase/0 2014.03.01 Revision:",i5," --- 3D_Parallel --")') svn_revision
codename = 'phaseUnif'
system = ''
diff -uprN phase0_2014.03/src_phase_3d/NEB.F90 phase0_2014.03.01/src_phase_3d/NEB.F90
--- phase0_2014.03/src_phase_3d/NEB.F90 2014-11-30 12:39:44.000000000 +0900
+++ phase0_2014.03.01/src_phase_3d/NEB.F90 2015-06-03 09:35:09.594604467 +0900
@@ -1250,6 +1250,8 @@ subroutine neb_force
real(8), allocatable :: fx0(:), fy0(:), fz0(:)
real(8), allocatable :: fx1(:), fy1(:), fz1(:)
+ integer :: ci_target=-1
+
if (mype==0) write(nfneb,*) 'neb force...'
allocate( k_spring(neb%number_of_images-1) )
@@ -1261,6 +1263,17 @@ subroutine neb_force
end do
max_energy = maxval( neb%image(:)%energy )
+ if( ( neb%step >= neb%cond%ci_neb_start) .and. ( neb%cond%climbing_image ) ) then
+ do k=2,neb%number_of_images-1
+ if (neb%image(k)%energy==max_energy) then
+ ci_target=k
+ write(nfneb,*) 'climing image : ',ci_target
+ exit
+ endif
+ enddo
+ endif
+
+
do k = 2, neb%number_of_images-1
@@ -1335,7 +1348,7 @@ subroutine neb_force
fz_spring(:) = k_spring0(image_id) * ( rz2(:) - rz0(:) ) &
+ k_spring0(image_id - 1 ) * ( rz1(:) - rz0(:) )
- if( ( neb%step < neb%cond%ci_neb_start ) .or. ( .not. neb%cond%climbing_image) ) then
+ if(k/=ci_target)then
! get true force perpendicular to local tangent
@@ -1459,7 +1472,7 @@ subroutine neb_force
write(nfneb,'(a31,3e11.3)') ' max force : ', &
maxval( abs( fx1 ) ), maxval( abs( fy1 ) ), maxval( abs( fz1 ) )
endif
- elseif( ( neb%step >= neb%cond%ci_neb_start) .and. ( neb%cond%climbing_image ) ) then
+ else
! get true force along local tangent
@@ -1965,6 +1978,7 @@ subroutine write_result(inc_energy)
open(nfdynm,file=F_DYNM)
open(nfenf,file=F_ENF)
call m_IS_wd_speciesname_etc(nfdynm)
+ call m_IS_wd_speciesname_etc(nfnebdynm)
write(nfenf,*) '#step image image_distance energy force_org force_neb force_normal'
sum_d = 0.0d0
do i=1, neb%number_of_images
@@ -1988,14 +2002,21 @@ subroutine write_result(inc_energy)
maxval(abs(neb%image(i)%force(:,:))), &
maxval(abs(neb%image(i)%true_force(:,:)))
endif
+ write(nfnebdynm,'(a,i5,i8,a)') ' cps and forc at (nebstep, iter_total=',neb%step,iteration,')'
write(nfdynm,'(a,i5,i8,a)') ' cps and forc at (nebstep, iter_total=',neb%step,iteration,')'
do j = 1, neb%image(i)%num_atom
- write(nfnebdynm,'(3i5,3f20.10)') neb%step, i, j, neb%image(i)%cps(j,1:3)
+! write(nfnebdynm,'(3i5,3f20.10)') neb%step, i, j, neb%image(i)%cps(j,1:3)
+ write(nfnebdynm,'(" ",i4,3f15.9,3f12.6)') j &
+ &, neb%image(i)%cps(j,1),neb%image(i)%cps(j,2),neb%image(i)%cps(j,3) &
+ &, neb%image(i)%force(j,1), neb%image(i)%force(j,2), neb%image(i)%force(j,3)
write(nfdynm,'(" ",i4,3f15.9,3f12.6)') j &
&, neb%image(i)%cps(j,1),neb%image(i)%cps(j,2),neb%image(i)%cps(j,3) &
&, neb%image(i)%force(j,1), neb%image(i)%force(j,2), neb%image(i)%force(j,3)
end do
end do
+
+ if(inc_energy) write(nfnebenf,*)
+
close(nfdynm)
close(nfenf)
! temporary
diff -uprN phase0_2014.03/src_phase_3d/Preparation.F90 phase0_2014.03.01/src_phase_3d/Preparation.F90
--- phase0_2014.03/src_phase_3d/Preparation.F90 2014-12-01 14:05:27.000000000 +0900
+++ phase0_2014.03.01/src_phase_3d/Preparation.F90 2015-06-03 09:35:09.594604467 +0900
@@ -150,7 +150,7 @@ subroutine Preparation()
& ,m_ES_cp_iconv, m_ES_add_neg
use m_Wannier, only:m_Wan_gen_weight
use m_FiniteElectricField, only:m_FEF_init
- use m_IterationNumbers, only : m_Iter_rd_iteration_numbers,iteration_unit_cell
+ use m_IterationNumbers, only : m_Iter_rd_iteration_numbers,iteration_unit_cell,m_Iter_cmix_reset
#ifdef _INCLUDE_EXX_
use m_ES_ExactExchange, only: m_ES_EXX_init
#endif
@@ -485,6 +485,8 @@ subroutine Preparation()
if ( sw_phonon == OFF ) call m_IS_gnrt_supercell_symmetry(paramset,nfout)
+ call m_Iter_cmix_reset()
+
contains
subroutine fft_box_finding_way(outer_or_inner)
integer,intent(out):: outer_or_inner
diff -uprN phase0_2014.03/src_phase_3d/m_ES_ortho.F90 phase0_2014.03.01/src_phase_3d/m_ES_ortho.F90
--- phase0_2014.03/src_phase_3d/m_ES_ortho.F90 2014-11-30 10:11:34.000000000 +0900
+++ phase0_2014.03.01/src_phase_3d/m_ES_ortho.F90 2015-06-03 09:37:30.810177002 +0900
@@ -1084,10 +1084,18 @@ contains
#endif
do i1 = L_NB_STA, L_NB_END
if(mode == ORTHONORMALIZATION .or. mode == NORMALIZATION) then
- call WSW_t_g(ik,i1,mod_pot,fr,psi_t,np_g1k_x,np_e,kimg,bpr_t,bpi_t) ! fr = 1/dsqrt(<Psi(i)|S|Psi(i)>)
- if(dabs(fr-1.d0) > DELTA) &
- & call normalize_bp_and_psi_t_g(ik,i1,fr,mod_pot &
- & ,psi_t,np_g1k_x,np_e,kimg,bpr_t,bpi_t)
+ if(mod_pot == VANDERBILT_TYPE) then
+ call WSW_t_g(ik,i1,mod_pot,fr,psi_t,np_g1k_x,np_e,kimg,bpr_t,bpi_t) ! fr = 1/dsqrt(<Psi(i)|S|Psi(i)>)
+ else
+ call WSW_t_g(ik,i1,mod_pot,fr,psi_t,np_g1k_x,np_e,kimg) ! fr = 1/dsqrt(<Psi(i)|S|Psi(i)>)
+ endif
+ if(dabs(fr-1.d0) > DELTA .and. mod_pot == VANDERBILT_TYPE) then
+ call normalize_bp_and_psi_t_g(ik,i1,fr,mod_pot &
+ & ,psi_t,np_g1k_x,np_e,kimg,bpr_t,bpi_t)
+ else if(dabs(fr-1.d0) > DELTA) then
+ call normalize_bp_and_psi_t_g(ik,i1,fr,mod_pot &
+ & ,psi_t,np_g1k_x,np_e,kimg)
+ endif
! |Psi(i)> = |Psi(i)> * fr, <beta|Psi(i)> = <beta|Psi(i)> * fr
if(mod_pot == VANDERBILT_TYPE) call cp_bpr2bprtw(i1,L_NB_END) ! bpr_t -> bpr_tw1, bpr_tw2
end if
@@ -1095,10 +1103,14 @@ contains
if(i1 == neg) cycle
call cp_psi2psii_g(ik,i1) ! psi_t(:,i1,:) -> psi_ir,psi_ii
call W1SW2_t_r_g(ik,i1,L_NB_END,mod_pot,psi_t,np_g1k_x,np_e,kimg) ! -> p1Sp2
- if(mod_pot == VANDERBILT_TYPE) &
- & call cp_bpr2bpi_g(kimg_t_wk, i1,bpr_t,bpi_t) ! -> bp_ir, bp_ii
- call modify_bp_and_psi_t_r_g(ik,i1,L_NB_END,mod_pot &
+ if(mod_pot == VANDERBILT_TYPE) then
+ call cp_bpr2bpi_g(kimg_t_wk, i1,bpr_t,bpi_t) ! -> bp_ir, bp_ii
+ call modify_bp_and_psi_t_r_g(ik,i1,L_NB_END,mod_pot &
& ,psi_t,np_g1k_x,np_e,kimg,bpr_t,bpi_t) ! psi_t, bpr_t, pbi_t, p1Sp2 -> psi_t, bpr_t, bpi_t
+ else
+ call modify_bp_and_psi_t_r_g(ik,i1,L_NB_END,mod_pot &
+ & ,psi_t,np_g1k_x,np_e,kimg) ! psi_t, bpr_t, pbi_t, p1Sp2 -> psi_t, bpr_t, bpi_t
+ endif
end if
end do ! i1-loop
#ifdef __TIMER_DO__
@@ -1233,10 +1245,14 @@ contains
local_block = nbsn(nbs)
L_NB_STA = nbsn_sta(local_block)
L_NB_END = nbsn_end(local_block)
- if(mod_pot == VANDERBILT_TYPE) call cp_bpr2bprtw(L_NB_STA,L_NB_END) ! bpr_t -> bpr_tw1, bpr_tw2
-
- call W1SW2_t_r_block_g(ik,i,L_NB_STA,p1Sp2_NB,NB_END,L_NB_END,kimg_t_wk &
- & , mod_pot,psi_t,psi_t_dia,np_g1k_x,np_e,kimg,bpr_tw1_dia,bpi_tw1_dia) ! -> p1Sp2
+ if(mod_pot == VANDERBILT_TYPE) then
+ call cp_bpr2bprtw(L_NB_STA,L_NB_END) ! bpr_t -> bpr_tw1, bpr_tw2
+ call W1SW2_t_r_block_g(ik,i,L_NB_STA,p1Sp2_NB,NB_END,L_NB_END,kimg_t_wk &
+ & , mod_pot,psi_t,psi_t_dia,np_g1k_x,np_e,kimg,bpr_tw1_dia,bpi_tw1_dia) ! -> p1Sp2
+ else
+ call W1SW2_t_r_block_g(ik,i,L_NB_STA,p1Sp2_NB,NB_END,L_NB_END,kimg_t_wk &
+ & , mod_pot,psi_t,psi_t_dia,np_g1k_x,np_e,kimg) ! -> p1Sp2
+ endif
icount = icount + 1
!-$ tune FUJITSU for block ----->>
@@ -1338,9 +1354,14 @@ contains
! end do ; enddo ; enddo
! endif
- call modify_bp_and_psi_t_r_blk_g(ik,i,L_NB_STA,p1Sp2_NB,NB_END,L_NB_END,kimg_t_wk &
- & , mod_pot,psi_t,psi_t_dia,np_g1k_x,np_e,kimg,bpr_t,bpi_t,bpr_t_dia,bpi_t_dia)
+ if(mod_pot == VANDERBILT_TYPE) then
+ call modify_bp_and_psi_t_r_blk_g(ik,i,L_NB_STA,p1Sp2_NB,NB_END,L_NB_END,kimg_t_wk &
+ & , mod_pot,psi_t,psi_t_dia,np_g1k_x,np_e,kimg,bpr_t,bpi_t,bpr_t_dia,bpi_t_dia)
! psi_t, bpr_t, pbi_t, p1Sp2 -> psi_t, bpr_t, bpi_t
+ else
+ call modify_bp_and_psi_t_r_blk_g(ik,i,L_NB_STA,p1Sp2_NB,NB_END,L_NB_END,kimg_t_wk &
+ & , mod_pot,psi_t,psi_t_dia,np_g1k_x,np_e,kimg)
+ endif
end do ! i2-loop2
#ifdef __TIMER_DO__
call timer_end(545)
@@ -3700,9 +3721,13 @@ contains
do i1 = L_NB_STA, L_NB_END
if(mode == ORTHONORMALIZATION .or. mode == NORMALIZATION) then
call WSW_t_g(ik,i1,mod_pot,fr,phi_t,np_g1k_x,np_e,kimg,phifr_t,phifi_t)
- if(dabs(fr-1.d0) > DELTA) &
- & call normalize_bp_and_psi_t_g(ik,i1,fr,mod_pot &
- & , phi_t,np_g1k_x,np_e,kimg,phifr_t,phifi_t)
+ if(dabs(fr-1.d0) > DELTA .and. mod_pot==VANDERBILT_TYPE) then
+ call normalize_bp_and_psi_t_g(ik,i1,fr,mod_pot &
+ & , phi_t,np_g1k_x,np_e,kimg,phifr_t,phifi_t)
+ else if (dabs(fr-1.d0) > DELTA) then
+ call normalize_bp_and_psi_t_g(ik,i1,fr,mod_pot &
+ & , phi_t,np_g1k_x,np_e,kimg)
+ endif
! |Phi(i)> = |Phi(i)> * fr, <beta|Phi(i)> = <beta|Phi(i)> * fr
!x!!$ if(mod_pot == VANDERBILT_TYPE) call cp_bpr2bprtw(i1,L_NB_END)
end if
@@ -3713,8 +3738,13 @@ contains
if(mod_pot == VANDERBILT_TYPE) &
& call cp_bpr2bpi_g(kimg_t_wk,i1,bpr_t,bpi_t) ! -> bp_ir, bp_ii
call W1SW2_t_r_g(ik,i1,L_NB_END,mod_pot,phi_t,np_g1k_x,np_e,kimg) ! -> p1Sp2
- call modify_bp_and_psi_t_r_g(ik,i1,L_NB_END,mod_pot &
- & ,phi_t,np_g1k_x,np_e,kimg,phifr_t,phifi_t)
+ if(mod_pot == VANDERBILT_TYPE) then
+ call modify_bp_and_psi_t_r_g(ik,i1,L_NB_END,mod_pot &
+ & ,phi_t,np_g1k_x,np_e,kimg,phifr_t,phifi_t)
+ else
+ call modify_bp_and_psi_t_r_g(ik,i1,L_NB_END,mod_pot &
+ & ,phi_t,np_g1k_x,np_e,kimg)
+ endif
! psi_t, bpr_t, pbi_t, p1Sp2 -> psi_t, bpr_t, bpi_t
end do ! i1-loop
if (sw_timing_2ndlevel == ON) call tstatc0_end(id_sname2)
@@ -3783,8 +3813,13 @@ contains
p1Sp2_NB(j_NB,i_NB,iq) = p1Sp2_t1_NB(j_NB,i_NB,iq,icount)
end do; enddo ; enddo
- call modify_bp_and_psi_t_r_blk_g(ik,i,L_NB_STA,p1Sp2_NB,NB_END,L_NB_END,kimg_t_wk &
- & , mod_pot,phi_t,psi_t_dia,np_g1k_x,np_e,kimg,phifr_t,phifi_t,bpr_t_dia,bpi_t_dia)
+ if(mod_pot==VANDERBILT_TYPE)then
+ call modify_bp_and_psi_t_r_blk_g(ik,i,L_NB_STA,p1Sp2_NB,NB_END,L_NB_END,kimg_t_wk &
+ & , mod_pot,phi_t,psi_t_dia,np_g1k_x,np_e,kimg,phifr_t,phifi_t,bpr_t_dia,bpi_t_dia)
+ else
+ call modify_bp_and_psi_t_r_blk_g(ik,i,L_NB_STA,p1Sp2_NB,NB_END,L_NB_END,kimg_t_wk &
+ & , mod_pot,phi_t,psi_t_dia,np_g1k_x,np_e,kimg)
+ endif
! phi_t, phifr_t, phifi_t, p1Sp2 -> phi_t, phifr_t, phifi_t
end do ! i2-loop
if (sw_timing_2ndlevel == ON) call tstatc0_end(id_sname3)
diff -uprN phase0_2014.03/src_phase_3d/m_Kpoints.F90 phase0_2014.03.01/src_phase_3d/m_Kpoints.F90
--- phase0_2014.03/src_phase_3d/m_Kpoints.F90 2014-11-30 23:16:06.000000000 +0900
+++ phase0_2014.03.01/src_phase_3d/m_Kpoints.F90 2015-06-03 09:35:09.594604467 +0900
@@ -878,12 +878,10 @@ contains
& gen_name_in_carts, &
& knv3_t_fbz, kv3_fbz, vkxyz_fbz, to_ibz_from_fbz_for_kpoint )
!-(b_Kpoints)
- if(ipri_kp_t >= 1) then
- write(nfout,'(" << trmat >>")')
- write(nfout,'(3f8.4)') ((trmat(i,j),j=1,3),i=1,3)
- end if
-
- call gen_kpoint_in_full_BZ( preallocation, kv3_fbz, vkxyz_fbz )
+! === modified by T. Y. 2015/01/16 ===
+!!$ call gen_kpoint_in_full_BZ( preallocation, kv3_fbz, vkxyz_fbz )
+ call gen_kpoint_in_full_BZ( preallocation)
+! ====================================
! ========================================================================= 12.0A
else if(way_ksample == MESH) then
@@ -2853,10 +2851,13 @@ contains
end subroutine m_Kp_dealloc
! === KT_add ==== 2014/09/30
- subroutine gen_kpoint_in_full_BZ( paramset, kv3_fbz, vkxyz_fbz )
+! === modified by T. Y. 2015/01/16 ===
+!!$ subroutine gen_kpoint_in_full_BZ( paramset, kv3_fbz, vkxyz_fbz )
+ subroutine gen_kpoint_in_full_BZ( paramset)
logical, intent(in) :: paramset
- integer, intent(in) :: kv3_fbz
- real(kind=DP), intent(inout) :: vkxyz_fbz(kv3_fbz,3,CARTS)
+!!$ integer, intent(in) :: kv3_fbz
+!!$ real(kind=DP), intent(inout) :: vkxyz_fbz(kv3_fbz,3,CARTS)
+! ====================================
integer :: i, ik
real(kind=DP) :: b1(3), b2(3), b3(3)
diff -uprN phase0_2014.03/src_phase_3d/m_Total_Energy.F90 phase0_2014.03.01/src_phase_3d/m_Total_Energy.F90
--- phase0_2014.03/src_phase_3d/m_Total_Energy.F90 2014-11-30 10:11:34.000000000 +0900
+++ phase0_2014.03.01/src_phase_3d/m_Total_Energy.F90 2015-06-03 09:35:09.595604967 +0900
@@ -498,6 +498,13 @@ contains
etotal0 = ekinet+ehartr+exc+elocal+enonlc+eewald-epc
endif
! ================= 13.0S
+#ifdef __FUJITSU
+ call mpi_bcast(etotal0,1,mpi_double_precision,0,mpi_comm_group,ierr)
+#else
+#ifdef __ETOTAL_UNIFY_
+ call mpi_bcast(etotal0,1,mpi_double_precision,0,mpi_comm_group,ierr)
+#endif
+#endif
if(sw_dipole_correction == ON) etotal0 = etotal0 + edip
if(sw_hubbard == ON) etotal0 = etotal0 + ehub0
diff -uprN phase0_2014.03/src_phase_3d/m_XC_Potential.F90 phase0_2014.03.01/src_phase_3d/m_XC_Potential.F90
--- phase0_2014.03/src_phase_3d/m_XC_Potential.F90 2014-11-30 23:16:06.000000000 +0900
+++ phase0_2014.03.01/src_phase_3d/m_XC_Potential.F90 2015-06-03 09:35:09.595604967 +0900
@@ -1686,6 +1686,10 @@ contains
& .and. xctype /= 'ggapbe' .and. xctype /= 'ldapbe' &
& .and. xctype /= 'katopbe' .and. xctype /= 'ggapbek' &
& .and. xctype /= 'ggapbe'.and. xctype /= 'ggapbex' &
+!================= rpbe, revpbe
+ & .and. xctype /= 'rpbe' &
+ & .and. xctype /= 'revpbe' &
+!================= rpbe, revpbe
#ifdef DISABLE_VDWDF
& .and. xctype /= 'ggapbey') then
#else