[gpaw-users] LCAO-TDDFT propagation on a charged system
Varadharajan Srinivasan
varadharajan.srinivasan at gmail.com
Wed Aug 25 08:27:42 CEST 2021
Dear Tuomas,
Thank you for your reply and my apologies for the delay in writing back. I
realised that KS decomposition may be more complicated so for now I used
the fix on the density calculation.
I will work on correcting these for KS decomposition as well.
I applied the patch you had supplied for the time propagation and it worked
fine. I will try and compile some examples including the systems we were
interested in to support this.
Once the corrections are incorporated I will also, as you suggested, post a
merge request.
Best,
Vardha.
On Tue, Aug 17, 2021 at 2:09 PM Tuomas Rossi <tuomas.rossi at aalto.fi> wrote:
> Dear Vardha,
>
>
> On 13.8.2021 15.01, Varadharajan Srinivasan wrote:
> > Dear Tuomas,
> >
> > I think their some related issue downstream as well when the density
> > matrix is processed for KS decomposition and density calculations. These
> > codes currently raise the error
> >
> > 'K-points are not fully supported'
> >
> > which I think is arising because thecheck that is done is
> >
> > if len(paw.wfs.kpt_u) > 1
> >
> > However, if I understand correctly in the spin-polarized (nspins=2)
> > case, even for gamma point only calculations, one kpoint object is
> > assigned for each spin.
> >
> > Is this correct?
>
> Yes, kpt_u in the code lists both k-points and spins. But, the KS
> decomposition code supports at the moment only spin-paired gamma point
> calculations.
>
>
> > I tested a temporary bypass for this which seemed to work
> >
> > gamma_only = len(self.kpt_u) == 1 or (len(self.kpt_u) == 2 and
> > paw.wfs.nspins == 2)
> >
> > assert gamma_only, 'K-points not implemented'
> >
> > and also modified the get_density method in lcaotddft/densitymatrix.py
> > to loop over the two spins to compute the total density at gamma point
> >
> > rho_G = density.gd.zeros()## Total density
> > rho_Gs = density.gd.zeros() ## Density from one spin at a time
> >
> > for iu in range(2):
> > kpt = wfs.kpt_u[iu]
> > assert kpt.q == 0
> > wfs.basis_functions.construct_density(rho_MM.astype(wfs.dtype),
> > rho_Gs, kpt.q)
>
> Here the same rho_MM is used for all spins/k-points, but the correct one
> for the index u from rho_uMM should be used.
>
>
> > rho_G += rho_Gs
> >
> > I haven't done extensive tests yet but it didn't throw any more errors.
> >
> > Do let me know if I am completely off here.
>
> This looks like a good start. Similar generalizations would need to be
> done throughout the KS decomposition code (see the issue with rho_MM
> pointed above) and then the code tested that it works correctly. You are
> welcome to submit a merge request with the updates if you like to work
> them out.
>
> Regarding the propagation scheme, the fix in
> https://gitlab.com/gpaw/gpaw/-/merge_requests/910 works now fine in
> simple tests. It would be great if you could test that it works
> expectedly in your systems too.
>
>
> Best regards,
>
> Tuomas
>
>
>
> > Thanks,
> >
> > Vardha.
> >
> >
> > On Tue, Aug 10, 2021 at 3:09 PM Varadharajan Srinivasan
> > <varadharajan.srinivasan at gmail.com
> > <mailto:varadharajan.srinivasan at gmail.com>> wrote:
> >
> > Dear Tuomas,
> >
> > Thanks for confirming and for fixing the issue. I am happy to test
> > it on some small systems if that could help.
> >
> > Best,
> > Vardha
> >
> > On Tuesday, August 10, 2021, Tuomas Rossi <tuomas.rossi at aalto.fi
> > <mailto:tuomas.rossi at aalto.fi>> wrote:
> >
> > Dear Vardha,
> >
> > It seems that the spin-polarized calculation was broken for LCAO
> > RT-TDDFT. Thank you for reporting the issue!
> >
> > A fix for the spin-polarized calculation is work in progress in
> > https://gitlab.com/gpaw/gpaw/-/merge_requests/910
> > <https://gitlab.com/gpaw/gpaw/-/merge_requests/910>>
> Before the fix, odd-electron systems would work only in
> > spin-paired mode (resulting in partial occupations).
> >
> > Best regards,
> >
> > Tuomas
> >
> >
> > On 5.8.2021 10.11, Varadharajan Srinivasan via gpaw-users wrote:
> >
> > Dear all,
> >
> > I was wondering whether there are any limitations on using
> > the LCAO-based RT-TDDFT on a singly charged (odd-electron)
> > system? I ran an example for Na2- and, while I was able to
> > perform the ground state calculation, the time-dependent run
> > crashed with the following error:
> > vt_G = hamiltonian.vt_sG[kpt.s]
> > IndexError: list index out of range
> >
> > I am able to run RT-TDDFT on the same system in FD mode just
> > fine. My ground state and tddft input scripts are given
> > below. I am using GPAW 21.6.0. Any help would be
> > appreciated. I apologize if this issue has been answered
> > before but I couldn't find any posts. Also, not sure if
> > there are any obvious mistakes I am making here.
> >
> > Thanks,
> > Vardha.
> >
> > GS:
> >
> > from ase import Atoms
> > from gpaw import GPAW
> > from gpaw.poisson import PoissonSolver
> > from gpaw.poisson_extravacuum import ExtraVacuumPoissonSolver
> >
> > # Sodium atom chain
> > atoms = Atoms('Na2',
> > positions=[[i * 3.0, 0, 0] for i in
> range(2)])
> > atoms.center(vacuum=6.0)
> >
> > # Use an advanced Poisson solver
> > ps = ExtraVacuumPoissonSolver(gpts=(512, 256, 256),
> >
> > poissonsolver_large=PoissonSolver(),
> > coarses=2,
> >
> > poissonsolver_small=PoissonSolver())
> >
> > # Ground-state calculation
> > calc = GPAW(mode='lcao', h=0.3, basis='pvalence.dz
> > <http://pvalence.dz/>> <http://pvalence.dz/
> > <http://pvalence.dz/>> xc='LDA', charge = -1.0,
> > setups={'Na': '1'}, spinpol = True,
> > poissonsolver=ps,
> > convergence={'density': 1e-12},
> > txt='gs.out')
> > atoms.calc = calc
> > energy = atoms.get_potential_energy()
> > calc.write('gs.gpw', mode='all')
> >
> > TD:from gpaw.lcaotddft import LCAOTDDFT
> > from gpaw.lcaotddft.dipolemomentwriter import
> DipoleMomentWriter
> > from gpaw.lcaotddft.wfwriter import WaveFunctionWriter
> >
> > # Read the ground-state file
> > td_calc = LCAOTDDFT('gs.gpw', txt='td.out')
> >
> > # Attach any data recording or analysis tools
> > DipoleMomentWriter(td_calc, 'dm.dat')
> > WaveFunctionWriter(td_calc, 'wf.ulm')
> >
> > # Kick and propagate
> > td_calc.absorption_kick([1e-5, 0., 0.])
> > td_calc.propagate(20, 1500)
> >
> > # Save the state for restarting later
> > td_calc.write('td.gpw', mode='all')
> >
> >
> >
> > _______________________________________________
> > gpaw-users mailing list
> > gpaw-users at listserv.fysik.dtu.dk
> > <mailto:gpaw-users at listserv.fysik.dtu.dk>
> > https://listserv.fysik.dtu.dk/mailman/listinfo/gpaw-users
> > <https://listserv.fysik.dtu.dk/mailman/listinfo/gpaw-users>>
>
> >
> >
>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://listserv.fysik.dtu.dk/pipermail/gpaw-users/attachments/20210825/03a64fa7/attachment-0001.htm>
More information about the gpaw-users
mailing list