Update solids example to use kpoint symmetry.

examples/ewf/solids/01-simple-sym.py

@@ -29,7 +29,8 @@
 emb.kernel()
 
 # Some versions of pyscf may encounter an error when generating IAOs with space_group_symmetry=True.
-# To avoid this, turn off space group symmetry before generating IAOs, as below:
+# To avoid this, either update to a version more recent that 2.3, or turn off space group symmetry before generating
+# IAOs, as in the lines below.
 # emb.mf.mol.space_group_symmetry = False
 # emb.kernel()
 

examples/ewf/solids/20-mRPA-interactions.py

@@ -10,28 +10,35 @@
 
 cell = pyscf.pbc.gto.Cell()
 cell.a, cell.atom = solids.graphene()
-cell.basis = 'cc-pVDZ'
+cell.basis = 'sto-3g'
 cell.output = 'pyscf.out'
 cell.dimension = 2
+cell.space_group_symmetry = True
+cell.symmorphic = True
 cell.build()
 
 # HF
-hf = pyscf.pbc.scf.RHF(cell)
-hf = hf.density_fit()#auxbasis='cc-pVDZ-ri')
+kmesh = [2,2,1]
+kpts = cell.make_kpts(kmesh, space_group_symmetry=True, time_reversal_symmetry=True, symmorphic=True)
+hf = pyscf.pbc.scf.KRHF(cell, cell.make_kpts([2,2,1]))
+hf = hf.density_fit()
 hf.kernel()
-
+# This may be required to avoid issues discussed in 01-simple-sym.py.
+hf.mol.space_group_symmetry = False
 # Run embedded
-emb_bare = vayesta.ewf.EWF(hf, bath_options=dict(bathtype="rpa", threshold=1e-1), screening=None)
+emb_bare = vayesta.ewf.EWF(hf, bath_options=dict(bathtype="rpa", threshold=1e-2), screening=None)
 emb_bare.kernel()
-
-emb_mrpa = vayesta.ewf.EWF(hf, bath_options=dict(bathtype="rpa", threshold=1e-1), screening="mrpa")
+# Run calculation using screened interactions and cumulant correction for nonlocal energy.
+emb_mrpa = vayesta.ewf.EWF(hf, bath_options=dict(bathtype="rpa", threshold=1e-2), screening="mrpa",
+                           ext_rpa_correction="cumulant")
 emb_mrpa.kernel()
 
 # Reference full system CCSD:
-cc = pyscf.pbc.cc.CCSD(hf)
+cc = pyscf.pbc.cc.KCCSD(hf)
 cc.kernel()
 
-print("E(HF)=              %+16.8f Ha" % hf.e_tot)
-print("E(Emb. bare CCSD)=  %+16.8f Ha" % emb_bare.e_tot)
-print("E(Emb. mRPA CCSD)=  %+16.8f Ha" % emb_mrpa.e_tot)
-print("E(CCSD)=            %+16.8f Ha" % cc.e_tot)
+print("Error(HF)=                     %+16.8f Ha" % (hf.e_tot - cc.e_tot))
+print("Error(Emb. bare CCSD)=         %+16.8f Ha" % (emb_bare.e_tot - cc.e_tot))
+print("Error(Emb. mRPA CCSD)=         %+16.8f Ha" % (emb_mrpa.e_tot - cc.e_tot))
+print("Error(Emb. mRPA CCSD + ΔE_k)=  %+16.8f Ha" % (emb_mrpa.e_tot+emb_mrpa.e_nonlocal-cc.e_tot))
+print("E(CCSD)=                       %+16.8f Ha" % cc.e_tot)