Matriz de superposición (molécula H2, STO-3G)

Question

Me gustaría implementar el método SCF utilizando la Química Cuántica Moderna (Szabo). Pero cuando quiero calcular la matriz de superposición para H2 (base STO-3G) mis resultados no son iguales a la verdadera matriz de superposición que hay en el libro. Sin embargo, cuando comparo mis resultados en otra base (STO-1G) obtengo la matriz de superposición exacta. A continuación, mi código para el cálculo de los valores de superposición.

class gauss():

  def __init__(self, alpha, Ra, k = 1):

    self.alpha = alpha
    self.Ra = np.array(Ra)
    self.k = k

  def __mul__(self, other):
    'calcular integral'
    longitud = (self.Ra - other.Ra)
    escalar = np.dot(longitud, longitud)

    return np.exp(-self.alpha*other.alpha*escalar/(self.alpha + other.alpha))*np.pi/(self.alpha + other.alpha)**(3/2)

def superposicion(primero, segundo):

  s11 = 0
  for i in primero:
    for j in primero:

        s11 += i * j * i.k * j.k

  s12 = 0
  for i in primero:
    for j in segundo:

        s12 += i * j * i.k * j.k

  return s11, s12

## primitivas gaussianas:
s1 = gauss(0.168855404, (0.0, 0.0, 0.0), 0.444634542)
s2 = gauss(0.623913730, (0.0, 0.0, 0.0), 0.535328142) 
s3 = gauss(3.42525091, (0.0, 0.0, 0.0), 0.154328967)   
s1b = gauss(0.168855404, (1.4, 0.0, 0.0), 0.444634542)
s2b = gauss(0.623913730, (1.4, 0.0, 0.0), 0.535328142)
s3b = gauss(3.42525091, (1.4, 0.0, 0.0), 0.154328967)

primero = [s1, s2, s3]
segundo = [s1b, s2b, s3b]

s11, s12 = superposicion(primero, segundo)

print(s12 / s11) # imprimir s12 normalizado

0.781180936263708

pero el valor exacto del libro = 0.6593 (1.4 Bohr)

y si elijo la base como:

s = gauss(0.4166127, (0.0, 0.0, 0.0))
sb = gauss(0.4166127, (1.4, 0.0, 0.0))

Obtengo el valor exacto:

0.6647924143035879

del libro: 0.6648

Answer 1

Tienes dos errores en tu código, uno menor y uno mayor.

El error menor está en tu función de multiplicación. El valor que devuelves debería tener $\pi^{3/2}$, pero de la forma en que está dispuesto ahora, el exponente solo se está aplicando al denominador de esa expresión.

El error mayor es cómo estás normalizando tus integrales de superposición. No puedes simplemente dividir $\langle1|2\rangle$ por $\langle1|1\rangle$ para normalizarlo. Debes normalizar explícitamente cada una de las funciones de base gaussianas. Para hacer esto, puedes incluir un factor de ((4*self.alpha*other.alpha)/np.pi**2)**(3/4) en el valor que estás devolviendo de tu multiplicación (en términos de legibilidad de tu código, puede ser más sabio definir variables dentro de la función como constante_de_proporcionalidad, integral_gaussiana y normalizacion y luego devolver el producto de estas). O podrías simplificar y simplemente multiplicar el exponencial por ((4*self.alpha*other.alpha)/(self.alpha+other.alpha)**2)**(3/4). Esto también significa que solo necesitas la segunda parte de tu función de superposición, ya que la multiplicación se encargará de la normalización.

La razón por la que tu código funciona para STO-1G es que dividir por $\langle1|1\rangle$ accidentalmente da la normalización correcta ya que solo hay una constante de normalización asociada con la única función gaussiana de STO-1G.