[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