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