[Q-e-developers] Re: About Charge Density

Paolo Giannozzi p.giannozzi at gmail.com
Mon Jun 15 14:57:27 CEST 2015


Symmetrization transforms a sum over the Irreducible Brillouin Zone (IBZ)
only:
F_u(r) = \sum_{k \in IBZ} w_k f_k(r) (where w_k are symmetry weights) into
the desired quantity
F(r) = \sum_{k \in BZ} f_k(r), where the sum is performed over all
k-points. It can be done in
real or in reciprocal space. If your additional term contains a sum over
k-points, you should
make the sum over the IBZ (that is, over the available k-points and
symmetry weight), add
your term before the call to "symrho"; otherwise, after it.

Paolo

On Sun, Jun 14, 2015 at 5:40 AM, Jiqiang Li <jqli14 at fudan.edu.cn> wrote:

> Dear Paolo,
>        I am afraid the potinit.f90 initializes the potential by calling
> subroutine "v_of_rho" at very beginning of a calculation. However, I want
> to calculate the potential by new charge density in every SCF loop, |i>
> changes with the iterations, but the array psi_m(r)*psi_n(r) is fixed. so I
> focus on the electron.f90 and modified the subroutine "sum_band".  I am
> not sure the correct place to add the additional term of charge density.
> I have some puzzle as mentioned in  the previous email:
>
>       (1)The subroutine *symrho* is marked as "symmetrize rho(G)", what
> does the "symmetrize charge density in G space' mean? and what's its
> consequence? Should I add the term before this symmetrization process or
> after it?
>
>       (2)I add the code after line 443 in subroutine "sum_band" in
> QE(v5.1) as follow, ( the band index of |m> and |n> are (nbnd -level_ccmg1)
> and (nbnd-level_ccmg2), respectively. )
>
> IF (loc_ccmg) THEN
>     !
>     !
>     IF ( nks > 1 ) REWIND( iunigk )
>     ik_ccmg = 1
>     !
>     IF (lsda) current_spin = isk(ik_ccmg)
>     npw = ngk(ik_ccmg)
>     w1_ccmg = 2.0 / omega
>     !
>     IF (nks > 1) THEN
>         READ( iunigk) igk
>
>         CALL get_buffer(evc, nwordwfc, iunwfc, ik_ccmg)  ! copy evc(1:nwordwfc) from the ik-th record in file with a unit of iunwfc, i.e. wavefucntion
>     END IF
>     !
>     IF (nkb > 0) CALL init_us_2(npw, igk, xk(1, ik_ccmg), vkb)
>     !
>     IF (noncolin) THEN
>     !
>     ELSE
>          IF ( io_wave_ccmg == 'write' ) THEN
>         psic(:) = (0.D0, 0.D0)
>         psic(nls(igk(1:npw))) =  evc(1:npw, nbnd-level_ccmg1) + &
>
>                                 (0.D0, 1.D0) * evc(1:npw, nbnd -level_ccmg2)
>         psic(nlsm(igk(1:npw))) = CONJG(evc(1:npw, nbnd -level_ccmg1) - &
>
>                                 (0.D0, 1.D0) * evc(1:npw, nbnd -level_ccmg2))
>              CALL invfft('Wave', psic, dffts)
>
>              CALL get_rho_gamma_ccmg(rho_r_ccmg(:,current_spin), dffts % nnr, psic, w1_ccmg)
>          END IF
>      ENDIF ! noncolinear or colinear
>  ENDIF  ! loc_ccmg or not
>  !
>  IF (loc_ccmg) THEN
>      IF (noncolin) THEN
>
>      ! the case of noncolinear has not been supported so far, remaining further work.
>      ELSE
>           IF ( io_wave_ccmg == 'write' ) THEN
>               CALL MPI_COMM_RANK(MPI_COMM_WORLD, myid_ccmg, ierr_ccmg)
>               WRITE( file_ccmg,'(I4)') myid_ccmg
>               file_ccmg = TRIM(ADJUSTL(file_ccmg))
>               file_ccmg = 'wave_ccmg/'//'rho'//TRIM(file_ccmg)//'.dat'
>               OPEN(unit = 12306, file = file_ccmg)
>               DO ig = 1, size(rho_r_ccmg,1)
>                   WRITE(12306,*) rho_r_ccmg(ig,current_spin)
>               END DO
>               CLOSE(12306)
>           END IF
>           ! if reading wavefunction.
>           IF ( io_wave_ccmg == 'read' ) THEN
>               CALL MPI_COMM_RANK(MPI_COMM_WORLD, myid_ccmg, ierr_ccmg)
>               WRITE( file_ccmg,'(I4)') myid_ccmg
>               file_ccmg = TRIM(ADJUSTL(file_ccmg))
>               file_ccmg = 'wave_ccmg/'//'rho'//TRIM(file_ccmg)//'.dat'
>               INQUIRE(file = file_ccmg, exist = alive_ccmg )
>
>               IF ( alive_ccmg == .false. ) CALL errore( 'sum_band', 'rho.dat does not exsit!', 1 )
>               OPEN(unit = 12306, file = file_ccmg)
>               DO ig = 1, size(rho_r_ccmg,1)
>                   READ(12306,*) rho_r_ccmg(ig,current_spin)
>               END DO
>               CLOSE(12306)
>
>                rho%of_r(:, current_spin)=rho%of_r(:, current_spin)+rho_r_ccmg(:,current_spin)
>           END IF
>       ENDIF ! noncolinear or colinear
>   ENDIF  ! loc_ccmg or not
>   !
>
>  Is the place correct?
>
> Jiqiang
>
>
>
> *From:* Paolo Giannozzi <p.giannozzi at gmail.com>
> *Date:* 2015-06-12 04:21
> *To:* Jiqiang Li <jqli14 at fudan.edu.cn>
> *Subject:* Re: About charge density
> On Thu, Jun 11, 2015 at 5:49 PM, Jiqiang Li <jqli14 at fudan.edu.cn> wrote:
>
>
>>       (3) Actually I want to realize these process:  First I run a whole
>> ordinary SCF calculation and at the end output the  wavefunction on the
>> disk. Then I run another new SCF calculation, but in this new calculation  a
>> new charge density, rho_new = rho + <m|n> = \sum_i(occupation) { <i|i>} +
>> <m|n>, is treated as input charge density in SCF loop . The term <m|n> is
>> fixed as a constant array  read from disk in every SCF loop. So Where
>> should I add the code to realize rho_new=rho + <m|n> ?
>>
>>
> In subroutine PW/src/potinit.f90, which either reads from file or
> generates the starting charge
> density (not the potential, actually: it is computed at the end by calling
> routine "v_of_rho").
> You might also use routine "plugin_scf_potential" to add your preferred
> term to the charge density.
>
> Paolo
>
>
> --
> Paolo Giannozzi, Dept. Chemistry&Physics&Environment,
> Univ. Udine, via delle Scienze 208, 33100 Udine, Italy
> hone +39-0432-558216, fax +39-0432-558222
>
> ------------------------------
>
> _______________________________________________
> Q-e-developers mailing list
> Q-e-developers at qe-forge.org
> http://qe-forge.org/mailman/listinfo/q-e-developers
>
>


-- 
Paolo Giannozzi, Dept. Chemistry&Physics&Environment,
Univ. Udine, via delle Scienze 208, 33100 Udine, Italy
hone +39-0432-558216, fax +39-0432-558222
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.quantum-espresso.org/pipermail/developers/attachments/20150615/6999cbed/attachment.html>


More information about the developers mailing list