Цель задания: рассчитать электронную и полную энергии молекулы водорода, построить 1D и 2D плот орбиталей.
from hfscf import * #испанский скрипт для расчета энергии и формы орбиталей молекулы водорода с помощью метода Хартри-Фока
import numpy as np
#Скрипт для метода самосогласованного поля
def SCF (r = 1.4632, Z=[1, 1], b1 = GTO["H"], b2 = GTO["H"], b = 3 , vbs=False):
#r - электронные степени свободы, R - ядерные степени свободы
R = [0, r]
#матрица перекрывания
if vbs: print ("*) Generating overlap matrix S.")
s_scf = S(R, b1, b2, b)
#гамильтониан
if vbs: print ("\n*) Generating Hamiltonian H.")
h_scf = H(R, Z, b1, b2, b)
#диагонализируем матрицу S и находим матрицу X
if vbs: print ("\n*) Diagonalizing matrix S and finding diagonal matrix X.")
X = diagon(s_scf)
Xa = X.getH()
#создаём матрицу плотности P
if vbs: print("\n*) Creating density matrix P.")
p_scf = np.matrix([[0,0],[0,0]], dtype=np.float64)
#итеративно оптимизируем матрицу P
if vbs: print("\n*) Starting with the SCF.")
for iteration in range(50):
#конструируем матрицу Фока
#F = H + G
if vbs: print("\n**) Generating the Fock matrix: calculating integral of two electrons.")
g_scf = G(r, p_scf, b1, b2, b)
f_scf = h_scf + g_scf
#смена базиса - генерируем матрицу F'
#F' = X_adj * F * X
if vbs: print("**) Changing the base of F.")
f_tra = Xa * f_scf * X
#диагонализируем F' и строим C'
if vbs: print("**) Diagonalizing F' and generating C'.")
c_tra = diagon2(f_tra)
#выводим матрицу коэффициента C
#C = X * C'
if vbs: print("**) Constructing coefficient matrix C.")
c_scf = X * c_tra
#из C получаем матрицу P
if vbs: print("**) Recalculating density matrix P.")
p_temp = P(c_scf)
print("\nFinished the " + str(iteration + 1) + ". iteration.\n")
#проверяем сходимость
if np.linalg.norm(p_temp - p_scf) < 1E-4:
#рассчитайтываем электронную энергию молекулы и полную энергию молекулы
energy_electron = ener_elec(p_temp, h_scf, f_scf)
print(energy_electron)
energy_sum = ener_tot(r, Z, 0)
print(energy_sum)
#строим графики 1D и 2D
orbital(c_tra, r, b1, b2, b)
orbital2D(c_scf, X, f_scf, r, b1, b2, b)
print("\n\n-->The self-consistent field IF has converged!")
return {"S" : s_scf, "H" : h_scf, "X" : X, "F" : f_scf, "C" : c_scf, "P" : p_temp}
else:
p_scf = p_temp
print("\n\n-->The self-consistent field has NOT converged!\nReview assumptions.")
return {"S" : s_scf, "H" : h_scf , "X" : X , "F" : f_scf ,"C" : c_scf, "P" : p_temp}
#Слайды из презентации
from IPython.display import Image, display
Image(filename='МСП.png')
SCF(vbs = True)
*) Generating overlap matrix S. *) Generating Hamiltonian H. *) Diagonalizing matrix S and finding diagonal matrix X. *) Creating density matrix P. *) Starting with the SCF. **) Generating the Fock matrix: calculating integral of two electrons. **) Changing the base of F. **) Diagonalizing F' and generating C'. **) Constructing coefficient matrix C. **) Recalculating density matrix P. Finished the 1. iteration. **) Generating the Fock matrix: calculating integral of two electrons. **) Changing the base of F. **) Diagonalizing F' and generating C'. **) Constructing coefficient matrix C. **) Recalculating density matrix P. Finished the 2. iteration. **) Generating the Fock matrix: calculating integral of two electrons. **) Changing the base of F. **) Diagonalizing F' and generating C'. **) Constructing coefficient matrix C. **) Recalculating density matrix P. Finished the 3. iteration. **) Generating the Fock matrix: calculating integral of two electrons. **) Changing the base of F. **) Diagonalizing F' and generating C'. **) Constructing coefficient matrix C. **) Recalculating density matrix P. Finished the 4. iteration. **) Generating the Fock matrix: calculating integral of two electrons. **) Changing the base of F. **) Diagonalizing F' and generating C'. **) Constructing coefficient matrix C. **) Recalculating density matrix P. Finished the 5. iteration. **) Generating the Fock matrix: calculating integral of two electrons. **) Changing the base of F. **) Diagonalizing F' and generating C'. **) Constructing coefficient matrix C. **) Recalculating density matrix P. Finished the 6. iteration. **) Generating the Fock matrix: calculating integral of two electrons. **) Changing the base of F. **) Diagonalizing F' and generating C'. **) Constructing coefficient matrix C. **) Recalculating density matrix P. Finished the 7. iteration. -1.7974485548087507 0.683433570256971
-->The self-consistent field IF has converged!
{'S': matrix([[0.99999999, 0.63749012], [0.63749012, 0.99999999]]), 'H': matrix([[-1.09920375, -0.91954841], [-0.91954841, -1.09920375]]), 'X': matrix([[ 0.55258063, 1.17442445], [ 0.55258063, -1.17442445]]), 'F': matrix([[-0.34684107, -0.57771139], [-0.57771139, -0.34684028]]), 'C': matrix([[ 0.55258114, -1.17442421], [ 0.55258013, 1.17442468]]), 'P': matrix([[0.61069182, 0.61069071], [0.61069071, 0.6106896 ]])}
Итого скрипт выдаёт графики и значения электронной (-1.7974485548087507) и полной (0.683433570256971) энергий, а также 1D и 2D графики.