1! --- 2! Copyright (C) 1996-2016 The SIESTA group 3! This file is distributed under the terms of the 4! GNU General Public License: see COPYING in the top directory 5! or http://www.gnu.org/copyleft/gpl.txt . 6! See Docs/Contributors.txt for a list of contributors. 7! --- 8 module m_born_charge 9 10 private 11 public :: born_charge 12 13 CONTAINS 14 15 subroutine born_charge() 16 use precision 17 use siesta_options, only:dx 18 use sparse_matrices 19 use parallel, only: Node, Nodes, IOnode 20 USE m_ksvinit, only: nkpol, kpol, wgthpol 21 use siesta_geom, only: na_u, na_s, xa, scell, ucell, isa 22 use atomlist, only: no_u, indxuo, iphorb, iaorb, lasto 23 use atomlist, only: rmaxo, no_s, no_l 24 use m_spin, only: nspin, SPpol, NonCol, SpOrb ! CC RC Added 25 use m_ksv, only: KSV_pol 26 use siesta_geom, only: shape 27#ifdef MPI 28 use m_mpi_utils, only: globalize_sum 29#endif 30 31 implicit none 32 33 real(dp) :: polxyz(3,nspin) ! Automatic, small 34 real(dp) :: polR(3,nspin) ! Automatic, small 35 36 real(dp) :: qspin(4) ! Local variable 37 38 integer :: j, ind, io, ispin 39#ifdef MPI 40 real(dp):: qtmp(4) ! Temporary for globalized spin charge 41#endif 42 43! Computes Born-effective charges by finite differences, 44! using the "force constant calculation" mode of operation. 45! It is called only for those steps (even) in which the dx 46! is positive. 47! 48 if (NonCol .or. SpOrb) then 49 if ( IONode ) then 50 write(*,'(/,a)') 51 . 'born_charge: implemented only for spin-polarized ' 52 write(*,'(a)') 53 . 'born_charge: or paramagnetic calculations' 54 end if 55 56 return 57 58 endif 59 60 if (nkpol.lt.1) then 61 if (IOnode) write(6,'(/,a,f12.6)') 62 . 'siesta: specify polarization grid for BC calculation' 63 if (IOnode) write(6,'(a,f12.6)') 64 . 'siesta: The Born charge matrix will not be calculated' 65 RETURN 66 endif 67 if (IOnode) write(6,'(/,a,f12.6)') 68 . 'siesta: Calculating polarization. ' 69 70! Find total population of spin up and down 71 if (nspin .ge. 2) then 72 do ispin = 1,nspin 73 qspin(ispin) = 0.0_dp 74 do io = 1,no_l 75 do j = 1,numh(io) 76 ind = listhptr(io) + j 77 qspin(ispin) = qspin(ispin) 78 . + Dscf(ind,ispin)*S(ind) 79 enddo 80 enddo 81 enddo 82#ifdef MPI 83! Global reduction of spin components 84 call globalize_sum(qspin(1:nspin),qtmp(1:nspin)) 85 qspin(1:nspin) = qtmp(1:nspin) 86#endif 87 endif 88 if (nkpol.gt.0) then 89 ! Note use of xa instead of xa_last, since this routine 90 ! is called at every geometry step, before moving the atoms. 91 ! 92 call KSV_pol(na_u, na_s, xa, rmaxo, scell, ucell, 93 . no_u, no_l, no_s, nspin, qspin, 94 . maxnh, nkpol, numh, listhptr, listh, 95 . H, S, xijo, indxuo, isa, iphorb, 96 . iaorb, lasto, shape, 97 . nkpol,kpol,wgthpol, polR, polxyz) 98 endif 99 if (nkpol.gt.0.and.IOnode) then 100 call obc( polR, ucell, dx, nspin) 101 endif 102 103 end subroutine born_charge 104 105 end module m_born_charge 106