Li2 ENS PDFT InversionΒΆ
[1]:
import numpy as np
import matplotlib.pyplot as plt
from CADMium import Pssolver, Psgrid, Partition, Inverter
import CADMium
a = 5.122/2
Za, Zb = 3,3
pol = 2
#Set up grid
NP = 7
NM = [4,4]
L = np.arccosh(15/a)
loc = np.array(range(-4,5)) #Stencil outline
grid = Psgrid(NP, NM, a, L, loc)
grid.initialize()
# ALPHA FRAGMENT
Nmo_a = [[1,2]]; Nmo_A = [[2,1]]
N_a = [[1,2]]; N_A = [[2,1]]
nu_a = 0.5
#Fragment b electrons
Nmo_b = [[1,2]]; Nmo_B = [[2,1]]
N_b = [[1,2]]; N_B = [[2,1]]
nu_b = 0.5
#Molecular elctron configuration
Nmo_m = [[3,3]]
N_m = [[3,3]]
part = Partition(grid, Za, Zb, pol, [Nmo_a, Nmo_A], [N_a, N_A], nu_a, [Nmo_b, Nmo_B], [N_b, N_B], nu_b, { "AB_SYM" : True,
# "ENS_SPIN_SYM" : False,
"interaction_type" : "dft",
"kinetic_part_type" : "inversion",
"hxc_part_type" : "exact",
# "k_family" : "gga",
# "ke_func_id" : 500,
})
print("be careful! Treating ensembles with nu_x=1.0 will break things")
#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" : "orbitalinvert",
"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})
be careful! Treating ensembles with nu_x=1.0 will break things
----> Begin SCF calculation for *Isolated* Fragments
Total Energy (a.u.) Inversion
__________________ ____________________________________
Iteration A B iters optimality res
___________________________________________________________________________________________
1 -8.62301 -8.62301 1.000e+00
2 -7.59515 -7.59515 1.530e-01
3 -7.38724 -7.38724 5.758e-02
4 -7.34928 -7.34928 2.753e-02
5 -7.34331 -7.34331 1.291e-02
6 -7.34300 -7.34300 6.027e-03
7 -7.34348 -7.34348 2.683e-03
8 -7.34364 -7.34364 1.377e-03
9 -7.34379 -7.34379 6.411e-04
10 -7.34386 -7.34386 2.991e-04
11 -7.34388 -7.34388 1.403e-04
12 -7.34390 -7.34390 6.598e-05
13 -7.34390 -7.34390 3.116e-05
14 -7.34390 -7.34390 1.481e-05
15 -7.34390 -7.34390 7.069e-06
16 -7.34390 -7.34390 3.397e-06
17 -7.34390 -7.34390 1.641e-06
18 -7.34390 -7.34390 7.943e-07
19 -7.34390 -7.34390 3.854e-07
20 -7.34390 -7.34390 1.874e-07
21 -7.34390 -7.34390 9.139e-08
22 -7.34390 -7.34390 4.474e-08
23 -7.34390 -7.34390 2.198e-08
24 -7.34390 -7.34390 1.082e-08
25 -7.34390 -7.34390 5.334e-09
----> Begin SCF calculation for *Interacting* Fragments
Total Energy (a.u.) Inversion
__________________ ____________________________________
Iteration A B iters optimality res
___________________________________________________________________________________________
1 -7.32470 -7.32470 10 +7.065e-15 +1.000e+00
2 -7.32892 -7.32892 11 +9.499e-15 +4.577e-02
3 -7.34091 -7.34091 9 +5.380e-15 +2.897e-02
4 -7.34086 -7.34086 10 +9.518e-15 +2.348e-02
5 -7.33601 -7.33601 8 +1.155e-13 +9.184e-03
6 -7.33587 -7.33587 7 +9.878e-13 +1.038e-02
7 -7.33853 -7.33853 7 +7.301e-15 +1.744e-03
8 -7.33926 -7.33926 6 +7.274e-15 +4.218e-03
9 -7.33811 -7.33811 6 +3.540e-14 +1.495e-03
10 -7.33739 -7.33739 6 +8.000e-15 +1.895e-03
11 -7.33772 -7.33772 6 +5.755e-15 +7.635e-04
12 -7.33819 -7.33819 5 +5.893e-15 +6.080e-04
13 -7.33821 -7.33821 5 +7.807e-15 +5.148e-04
14 -7.33801 -7.33801 5 +9.397e-15 +1.349e-04
15 -7.33791 -7.33791 4 +8.040e-14 +2.092e-04
16 -7.33795 -7.33795 4 +9.107e-14 +4.116e-05
17 -7.33801 -7.33801 4 +7.421e-15 +7.296e-05
18 -7.33802 -7.33802 4 +5.245e-15 +2.668e-05
19 -7.33801 -7.33801 4 +5.714e-15 +2.499e-05
20 -7.33800 -7.33800 4 +5.232e-15 +1.361e-05
21 -7.33800 -7.33800 4 +8.038e-15 +8.586e-06
22 -7.33800 -7.33800 4 +5.876e-15 +7.432e-06
23 -7.33800 -7.33800 3 +8.268e-15 +1.865e-06
24 -7.33800 -7.33800 3 +7.925e-15 +3.569e-06
25 -7.33800 -7.33800 3 +7.593e-15 +7.636e-07
26 -7.33800 -7.33800 3 +6.049e-15 +1.358e-06
27 -7.33800 -7.33800 3 +1.059e-14 +4.485e-07
28 -7.33800 -7.33800 2 +5.002e-13 +4.670e-07
29 -7.33800 -7.33800 2 +2.260e-13 +2.700e-07
30 -7.33800 -7.33800 2 +1.560e-13 +1.276e-07
31 -7.33800 -7.33800 2 +3.124e-14 +1.397e-07
32 -7.33800 -7.33800 2 +1.249e-14 +2.971e-08
33 -7.33800 -7.33800 2 +5.239e-15 +6.000e-08
34 -7.33800 -7.33800 2 +1.284e-14 +1.438e-08
[2]:
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: 5.122
Fragment Energy: -14.675995635145341
Partition Energy: -1.8071121141333901
Vnn Energy 1.757126122608356
Total Energy: -14.725981626670373