1import numpy as np
2from ase.build import molecule
3from ase.parallel import parprint, world
4
5from gpaw import GPAW
6from gpaw.test import equal
7from gpaw.cluster import Cluster
8from gpaw.analyse.overlap import Overlap
9from gpaw.lrtddft import LrTDDFT
10
11"""Evaluate the overlap between two independent calculations
12
13Differences are forced by different eigensolvers and differing number
14of Kohn-Sham states.
15"""
16
17
18def test_overlap(in_tmp_dir):
19    h = 0.4
20    box = 2
21    nbands = 4
22    txt = '-'
23    txt = None
24    np.set_printoptions(precision=3, suppress=True)
25
26    H2 = Cluster(molecule('H2'))
27    H2.minimal_box(box, h)
28
29    c1 = GPAW(h=h, txt=txt, eigensolver='dav', nbands=nbands,
30              convergence={'eigenstates': nbands})
31    c1.calculate(H2)
32    lr1 = LrTDDFT(c1)
33
34    parprint('sanity --------')
35    ov = Overlap(c1).pseudo(c1)
36    parprint('pseudo(normalized):\n', ov)
37    ov = Overlap(c1).pseudo(c1, False)
38    parprint('pseudo(not normalized):\n', ov)
39    ov = Overlap(c1).full(c1)
40    parprint('full:\n', ov)
41    equal(ov[0], np.eye(ov[0].shape[0], dtype=ov.dtype), 1e-10)
42
43    def show(c2):
44        c2.calculate(H2)
45        ov = Overlap(c1).pseudo(c2)
46        parprint('wave function overlap (pseudo):\n', ov)
47        ov = Overlap(c1).full(c2)
48        parprint('wave function overlap (full):\n', ov)
49        lr2 = LrTDDFT(c2)
50        ovkss = lr1.kss.overlap(ov[0], lr2.kss)
51        parprint('KSSingles overlap:\n', ovkss)
52        ovlr = lr1.overlap(ov[0], lr2)
53        parprint('LrTDDFT overlap:\n', ovlr)
54
55    parprint('cg --------')
56    c2 = GPAW(h=h, txt=txt, eigensolver='cg', nbands=nbands + 1,
57              convergence={'eigenstates': nbands + 1})
58    show(c2)
59
60    parprint('spin --------')
61    H2.set_initial_magnetic_moments([1, -1])
62    c2 = GPAW(h=h, txt=txt, spinpol=True, nbands=nbands + 1,
63              parallel={'domain': world.size},
64              convergence={'eigenstates': nbands + 1})
65    H2.set_initial_magnetic_moments([0, 0])
66    try:
67        show(c2)
68    except AssertionError:
69        parprint('Not ready')
70    if 1:
71        ov = Overlap(c1).pseudo(c2, otherspin=1)
72        parprint('wave function overlap (pseudo other spin):\n', ov)
73
74    parprint('k-points --------')
75
76    H2.set_pbc([1, 1, 1])
77    c1 = GPAW(h=h, txt=txt, nbands=nbands,
78              kpts=(1, 1, 3),
79              convergence={'eigenstates': nbands})
80    c1.calculate(H2)
81    c2 = GPAW(h=h, txt=txt, nbands=nbands + 1,
82              kpts=(1, 1, 3),
83              convergence={'eigenstates': nbands + 1})
84    try:
85        show(c2)
86    except (AssertionError, IndexError) as e:
87        parprint('Not ready', e)
88