H2 FOO PDFT InversionΒΆ
[3]:
import numpy as np
import matplotlib.pyplot as plt
from CADMium import Pssolver, Psgrid, Partition, Inverter
import CADMium
a = 1.446/2
Za, Zb = 1,1
pol = 2
#Set up grid
NP = 7
NM = [4,4]
L = np.arccosh(10/a)
loc = np.array(range(-4,5)) #Stencil outline
grid = Psgrid(NP, NM, a, L, loc)
grid.initialize()
# ALPHA FRAGMENT
Nmo_a = [[1 , 1]];
N_a = [[0.5 , 0.5]];
nu_a = 1.0
#Fragment b electrons
Nmo_b = [[1 , 1]];
N_b = [[0.5 , 0.5]];
nu_b = 1.0
#Molecular elctron configuration
Nmo_m = [[1,1]]
N_m = [[1,1]]
part = Partition(grid, Za, Zb, pol, Nmo_a, N_a, nu_a, Nmo_b, N_b, nu_b, { "AB_SYM" : True,
"fractional" : True,
# "ENS_SPIN_SYM" : False,
"interaction_type" : "dft",
"kinetic_part_type" : "libxcke",
"hxc_part_type" : "exact",
"k_family" : "gga",
"ke_func_id" : 500,
})
#Setup inverter object
mol_solver = Pssolver(grid, Nmo_m, N_m)
part.inverter = Inverter(grid, mol_solver, { "AB_SYM" : True,
# "ENS_SPIN_SYM" : False,
"use_iterative" : False,
"invert_type" : "wuyang",
"DISP" : False,
})
part.optPartition.isolated = True
part.scf({"disp" : True,
"alpha" : [0.6],
"e_tol" : 1e-8})
part.optPartition.isolated = False
part.scf({"disp" : True,
"alpha" : [0.6],
"max_iter" : 200,
"e_tol" : 2e-8,
"iterative" : False,
"continuing" : True})
----> Begin SCF calculation for *Isolated* Fragments
Total Energy (a.u.)
__________________
Iteration A B res
_______________________________________________________
1 -0.48020 -0.48020 1.000e+00
2 -0.45519 -0.45519 2.286e-01
3 -0.44824 -0.44824 6.783e-02
4 -0.44635 -0.44635 1.880e-02
5 -0.44584 -0.44584 5.089e-03
6 -0.44571 -0.44571 1.351e-03
7 -0.44567 -0.44567 3.488e-04
8 -0.44566 -0.44566 8.599e-05
9 -0.44566 -0.44566 1.935e-05
10 -0.44566 -0.44566 3.488e-06
11 -0.44566 -0.44566 1.179e-06
12 -0.44566 -0.44566 4.803e-07
13 -0.44566 -0.44566 2.149e-07
14 -0.44566 -0.44566 1.172e-07
15 -0.44566 -0.44566 5.666e-08
16 -0.44566 -0.44566 2.583e-08
17 -0.44566 -0.44566 1.141e-08
18 -0.44566 -0.44566 4.941e-09
----> Begin SCF calculation for *Interacting* Fragments
Total Energy (a.u.)
__________________
Iteration A B res
_______________________________________________________
1 -0.30808 -0.30808 1.000e+00
2 -0.39247 -0.39247 2.612e-01
3 -0.41055 -0.41055 5.422e-02
4 -0.41386 -0.41386 9.816e-03
5 -0.41447 -0.41447 1.741e-03
6 -0.41460 -0.41460 3.968e-04
7 -0.41463 -0.41463 1.683e-04
8 -0.41463 -0.41463 7.136e-05
9 -0.41464 -0.41464 3.026e-05
10 -0.41464 -0.41464 1.280e-05
11 -0.41464 -0.41464 5.404e-06
12 -0.41464 -0.41464 2.288e-06
13 -0.41464 -0.41464 9.689e-07
14 -0.41464 -0.41464 4.105e-07
15 -0.41464 -0.41464 1.740e-07
16 -0.41464 -0.41464 7.379e-08
17 -0.41464 -0.41464 3.137e-08
18 -0.41464 -0.41464 1.345e-08
[4]:
print("Separation Distance:", 2*a)
print("Fragment Energy:", part.E.Ef)
print("Partition Energy:", part.E.Ep)
print("Vnn Energy", part.E.Vnn)
print("Total Energy:", part.E.E)
Separation Distance: 1.446
Fragment Energy: -0.8292816151184538
Partition Energy: -0.9999712791390211
Vnn Energy 0.6915629322268326
Total Energy: -1.1376899620306422
[ ]: