INTEGRACION MONTE CARLO EN EL METODO DE ROOTHAAN

Manuel FERNANDEZ NUÑEZ Departamento de Química Física Universidad de Cádiz y Rosa RODRIGUEZ HUERTAS Departamento de Matemáticas Universidad de Cádiz

FACULTAD DE CIENCIAS, Apartado 40 11510 PUERTO REAL (CADIZ), España

Resumen En el presente trabajo se presenta un procedimiento para la evaluación de las integrales mono y bielectrónicas del método de Roothaan por técnicas Monte Carlo. Tras describir el método de muestreo y la optimización de la función densidad empleada, se analizan los resultados de operar en el caso de algunos sistemas bielectrónicos. La comparación de los resultados obtenidos con los procedentes de integraciones exactas prueba que el procedimiento propuesto es lo bastante eficaz como para ser de utilidad práctica, por ejemplo en el tratamiento de problemas con funciones de base no habituales.

Abstract In the present paper we present a procedure to evaluate the mono- and bielectronic integrals arising from the Roothaan method by means of Monte Carlo integration . The sampling method used, the optimization of the density function, and some applications to several bielectronic systems are discussed. The comparison of obtained results to these coming from exact integration prove the utility of the procedure, that could be useful in treating calculations with non standard basis functions.

1. INTRODUCCION La integración Monte Carlo ha sido poco usada hasta el momento en problemas de Química Cuántica debido, probablemente, a la dificultad de obtener con ella resultados suficientemente precisos con tiempos de cálculo aceptables. En efecto, es bien sabido [1] que el tiempo de cálculo en la integración Monte Carlo crece -en principio- con el cuadrado de la precisión deseada. Sin embargo, esta clase de integración tiene a su favor, por lo menos, tres ventajas importantes: a) La evolución de los medios de cálculo a disposición de los químicos cuánticos. Cálculos prácticamente impensables hace solo una década, pueden realizarse hoy incluso en el laboratorio mas modesto. b) La dificultad de la integración Monte Carlo resulta relativamente independiente del número de dimensiones de la integral considerada. Esto hace a este tipo de integración especialmente atractiva para tratar las integrales bielectrónicas o, mejor aún, las integrales multielectrónicas como las que aparecen en los métodos de minimización de la varianza [2]. c) Salvo en algunos casos muy especiales, la dificultad de la integración Monte Carlo depende muy poco del tipo de integrando. Por ejemplo, desde este punto de vista las integrales que aparecen empleando una base gaussiana y las que aparecen con una base de Slater, tienen la mísma dificultad. De estos tres factores positivos, nos parece especialmente importante el tercero, ya que implica la posibilidad de estandarizar los cálculos quimicocuánticos con funciones de base arbitrarias, o con potenciales distintos del electrostático, tan costosos en la actualidad incluso a niveles exploratorios- por la dificultad de obtener las fórmulas de integración necesarias para llevarlos a cabo. Por este motivo, nos pareció interesante realizar una investigación acerca de las posibilidades de aplicar este método de integración. El trabajo incluye un estudio del cálculo de las integrales monoelectrónicas y bielectrónicas, seguido de un análisis de la aplicación de los resultados en algunos cálculos Roothaan sencillos, que ponen de manifiesto que los resultados son suficientemente precisos incluso con tiempos de cálculo bastante moderados. Conviene resaltar que estos ejemplos se han llevado a cabo con ordenadores personales, por lo que en caso de disponerse de medios de cálculo de mayor entidad podrían tratarse, de la mísma forma, problemas muchísimo mas complejos.

2 RESUMEN DE LOS FUNDAMENTOS DEL METODO DE INTEGRACION MONTE CARLO La integración MC puede basarse en la siguiente propiedad [3]: Si {ri} son variables aleatorias idénticamente distribuidas con una función densidad w(r), normalizada en un dominio D, y g es una función de ri entonces la media (esperanza matemática) cumple: 1 N r  r E  ∑ g(r ) =Eg(r ) i   i  N i

(1)

siendo la varianza: 1 N r  1 r var  ∑ g(r ) = var g(r ) i i   N N  i 

Si

r g(ri )

es contínua:

(2)

E [g(rri)]= ∫D g(rri)w(rri)dτ

(3)

1 N  E  ∑ g(rri) = ∫D g(rri)w(rri)dτ N i   

(4)

y por tanto:

En consecuencia, una forma de calcular la integral: F(D)= ∫D f(rr) dτ

(5)

consiste en descomponer f(r) en el producto de una función auxiliar g(r) por una función densidad w(r), muestrear aleatoriamente ri con la densidad w(r) y estimar F(D) mediante la media de los valores g(ri), ya que esta media constituye una estimación centrada de F(D). El error estandar de esta estimación es [4]: ε=

varg σ = N N

(6)

lo que indica que la disminución del error en el cálculo de la integral al aumentar el número de puntos empleados para calcularla, es lenta. El camino para aumentar la eficiencia del método de Monte Carlo en el cálculo de las integrales pasa entonces por la elección de una descomposición del integrando f(r) = g(r) · w(r) que sea adecuada para que la constante σ de la ecuación (6) resulte lo menor posible.

Es bien sabido [5], que si w(r) es proporcional a f(r) resulta var [ g ] = 0 con lo que, idealmente, convendría emplear como densidad el resultado de normalizar el integrando. Así sería σ = 0 y se obtendría el resultado deseado de forma óptima. Desgraciadamente este camino queda muy limitado en la práctica por la dificultad de encontrar el método de muestreo adecuado para la densidad w(r), de manera que lo más que puede hacerse es utilizar una densidad que siendo fácil de muestrear, se parezca en lo posible al resultado de normalizar f(r). Cuando la función densidad es constante, el muestreo es un problema bien establecido [6]. Existen programas eficientes para generar números aleatorios 0 < ηi < 1 , por ejemplo los suministrados en [7]. Sin embargo es importante llamar la atención sobre un problema que presentan muchos de estos programas: En realidad no suministran números realmente aleatorios, sino pseudoaleatorios. Ello se traduce en que a partir de la generación de cierta cantidad de datos, típicamente 105-106, cumplen su periodo y empiezan a repetir datos, con lo que no permiten aumentar tanto como se quiera el número de puntos empleados para calcular la integral tratada. Este defecto lo tienen los programas RAN1 y RAN2 de la referencia [7], que no permiten emplear mas allá de unos 7·105 puntos por integral. El programa RAN3, en cambio, hemos comprobado que permite emplear no menos de 109 puntos, lo cual es ya más que suficiente para nuestros propósitos. Los programas generadores de números aleatorios uniformemente distribuidos ηi sirven también, de forma indirecta, para generar números aleatorios con otras distribuciones, xi . Mediante la relación [8] : xi ∫ w(x)dx =η i a

(7)

que genera, a partir de un número aleatorio uniforme ηi , otro xi con una función densidad w(x) definida en a < x < b, La obtención del valor xi a partir de la integral (7) no siempre es factible. Precisamente por eso, al plantear las integraciones via Monte Carlo no se puede utilizar una densidad proporcional al integrando, como sería de desear, ya que ello implicaría saber resolver analíticamente la propia integral propuesta (y sobraría el Monte Carlo). Hay que utilizar una densidad "lo mejor posible", entre las que puedan ser muestreadas, vía fórmula (7) o procedimientos alternativos como el propuesto en [9].

3. INTEGRALES MONOELECTRONICAS Las integrales de tipo: I=∫∫∫ fµ fν rn + 2 senθdr dθ dφ

(8)

han sido calculadas empleando una función densidad: w(r, θ, φ)=

p −pr senθ 4π e

(9)

a la cual corresponden, aplicando (6), las fórmulas de muestreo: r =− i

ln(η ) i p

(10)

θ= arc cos(1−2η ) i

(11)

φ = 2πη i i

(12)

Una ventaja importante de la función densidad empleada es que resuelve de manera sencilla el problema de manipular el intervalo infinito en la variable radial. Otra ventaja, no menos importante, es que permite optimizar el parámetro p de manera que la eficiencia de la integración sea máxima. En efecto, si la integral a evaluar es: I=∫ F(x)dx

(13)

y se emplea una densidad de tipo p·exp(-px), la varianza de la estimación MC de la integral es: 2  pexp(−px)dx−  ∫

  F(x) varI= ∫   p • exp(−px) 

2  F(x) pexp(−px)  p • exp(−px) 

2 F2(x) =∫ dx − [∫ F(x)dx] p • exp(−px)

(14)

Minimizando esta varianza respecto al parámetro p, se obtiene: 1 2 F (x)exp(px)dx= ∫ xF2(x)exp(px)dx p∫

(15)

Aplicando esta condición en el caso en que la integral sea del tipo (8) con orbitales de base de tipo Slater, con lo que: nµ +nν +n +2 0 exp− (α + α )r  F(r) = r µ ν  

se deduce que, para este caso: p opt =

αµ + αν nµ + nν + nο + 3

(16)

(17)

La utilidad de optimizar el valor del parámetro p puede verse en el ejemplo presentado en la tabla 1. En ella el error estandar que se asigna a cada valor se ha determinado a partir del resultado de 10 estimaciones independientes de la integral, y el tiempo de cálculo se refiere a un PC486 a 33Mhz. Los resultados de la tabla, correspondientes a una integral , son representativos de los que se han obtenido para todas las integrales monoelectrónicas que se han calculado. El conjunto de puntos concreto que se emplee para calcular la integral influye, desde luego, en el resultado. Las variaciones encontradas resultan del orden del error estandar indicado en la última columna . Es interesante subrayar la buena coincidencia obtenida entre este error, calculado como si no se conociese el resultado exacto, y el error verdadero.

TABLA 1: Influencia de las condiciones del cálculo sobre el valor estimado para la integral . p α1 + α 2

numero de datos

tiempo de calculo

resultado del calculo

error estandar estimado

0.99599

0.0023

0.02

106

0.10

"

"

"

1.00169

0.0012

0.20(opt)

"

"

"

1.00101

0.00094

0.40

"

"

"

0.99964

0.0010

2.00

"

"

"

1.04436

0.129

opt

103

0.1 s.

0.973

0.027

opt

105

10.3 s.

1.0024

0.0029

opt

107

17 m. 6 s.

1.00039

0.00025

opt

108

2 h. 51 m.

1.000075

0.000081

opt

109

28 h. 28 m.

1.000034

0.000035

1 m. 43 s.

Otra cuestión a tener en cuenta es la influencia que pueda tener el tipo de integral calculada sobre la precisión del resultado. En la tabla 2 puede observarse que -como era de esperar- la clase de integral influye algo sobre la eficiencia del método, pero también se ve que las distintas integrales monoelectrónicas estudiadas -representativas de las que pueden aparecer en cálculos quimicocuánticos corrientes- se obtienen con errores de orden no muy distinto. Por consiguiente, al abordar los cálculos Roothaan, se ha descartado el empleo de diferente número de puntos para cada clase de integral.

TABLA 2: Influencia de la clase de integral. (Cálculos con 106 puntos, y p dado por la fórmula 17 )

Clase de integral

valor exacto

valor obtenido

error estimado

tiempo empleado

< 1s | 1s >

1.000000

1.00093

0.00076

1m. 43s.

< 1s | 2s >

0.866025

0.86689

0.00077

1m. 43s.

< 2s | 2s >

1.000000

1.00101

0.00094

1m. 43s.

< 2s | 2pz >

0.000000

-0.00019

0.0014

1m. 52s.

< 2pz | 2pz >

1.000000

1.00192

0.00138

2m. 2s.

< 2s| z |2pz >

2.886751

2.89272

0.0042

1m. 53s.

4 INTEGRALES BIELECTRONICAS Las integrales del tipo:

1 (µ, ν λ, σ)= ∫∫ fµ(1)fν(1) f (2)fσ(2)dV dV λ 1 2 r 12

(18)

se han calculado empleando una función densidad:

p p w(r r )= 1 2 exp(−p r − p r )senθ senθ 1 2 16π2 11 1 22 2

(19)

que corresponde a fórmulas de muestreo como (10), (11) y (12) aplicadas independientemente a cada electrón. Para orbitales de base de tipo Slater con potencia n y exponente α se ha tomado:

αµ + αν p = 1 nµ + nν + 3

(20)

α + ασ p = λ 2 n + nσ + 3 λ

(21)

Estas fórmulas están sugeridas por la (17) y su utilidad puede comprobarse en la tabla 3, en la cual se comparan los errores estandar de diversas integrales (µ,ν|λ,σ) obtenidas con p1 y p2 dados por (20) y (21), con los errores obtenidos utilizando el doble de estos valores, y la mitad. Puede verse que los valores dados por (20) y (21) son suficientemente buenos y, además, que cualquier intento de precisarlos mejor carecería de utilidad, pues la mejora hipotéticamente obtenida quedaría enmascarada por el "ruido" numérico que acompaña a los cálculos Monte Carlo.

TABLA 3: Comparación de los errores estandar obtenidos en diversas clases de integral (µ,ν|λ,σ), con 105 puntos y distintos valores de p1 y p2 .

INTEGRAL

p1 (20)

p2 (21)

σ

σ (0.5p

1

0.5p2)

(p

1

p2 )

σ (2p

1

2p2)

(1s,1s|1s,1s)

0.667

0.667

0.00288

0.00205

0.00416

(2s,2s|2s,2s)

0.400

0.400

0.00264

0.00184

0.00344

(3s,3s|3s,3s)

0.286

0.286

0.00248

0.00161

0.00298

(4s,4s|4s,4s)

0.222

0.222

0.00228

0.00143

0.00267

(2s,2s|4s,4s)

0.400

0.222

0.00208

0.00174

0.00273

(1s,2s|3s,4s)

0.500

0.250

0.00153

0.00138

0.00231

La tabla 4, por su parte, permite observar la influencia de diversas condiciones de cálculo sobre una integral representativa y comprobar, de paso, que los tiempos de cálculo no son mucho mayores que en el caso de las integrales monoelectrónicas.

TABLA 4: Influencia de las condiciones del cálculo sobre el valor obtenido para (2s,2s|3s,3s), con STOs de exponente unidad. (Su valor exacto es 0.293620). tiempo de cálculo

resultado

error estandar estimado

4 m. 39 s.

0.29464

0.00057

p1

p2

numero de datos

0.400

0.286

106

0.800

0.592

"

"

"

0.29507

0.00083

0.400

0.592

"

"

"

0.29400

0.00082

0.800

0.286

"

"

"

0.29292

0.00102

0.400

0.286

103

0.3 s.

0.26844

0.018

"

"

105

28 s.

0.29609

0.0018

"

"

107

46 m. 33 s.

0.29367

0.00022

"

"

108

7 h. 42 m.

0.29368

0.000060

5 APLICACION A SISTEMAS BIELECTRONICOS Se han llevado a cabo cálculos Roothaan con integración MC en unos cuantos sistemas bielectrónicos representativos, que permiten poner de manifiesto la utilidad y las limitaciones de esta metodología. En primer lugar se ha estudiado el estado fundamental del átomo de helio con base DZ [10]. Dado que los resultados de la integración MC dependen en algún grado de cual sea el conjunto de puntos seleccionado para la integración, resulta posible hacer una estimación del error estandard de los cálculos realizados: basta hacer estos un número suficiente de veces, con conjuntos de puntos independientes, y aplicar la definición de error estandard. Esto es lo que se ha hecho al final de las tablas 5 y 6. En la tabla 5 pueden verse los coeficientes C1 y C2 obtenidos para el orbital ocupado: 1/2

 3 α  φ=C  1  1 π     

1/2

 3 α  exp(−α r) + C  2  1 2 π     

exp(−α r) 2

(22)

la energía orbital y la energía total,cuando se emplean los exponentes de Clementi y Roetti [10]: α

1

= 1 . 45 α

2

= 2 . 91

(23)

y las integrales se calculan por Monte Carlo con solo 104 puntos. En la tabla aparecen los resultados obtenidos con diez conjuntos de puntos distintos, a partir de los cuales se calculan los errores estandar indicados al final de la tabla. La tabla 6 es análoga a la 5, con integraciones realizadas con 105 puntos. Puede verse que las energías Roothaan-Monte Carlo, incluso si se emplea un número de puntos tan pequeño como 104, tienen errores respecto a las energías exactas que son del mismo orden que los del propio método de Roothaan con integración analítica.

TABLA 5: Monte Carlo con 104 puntos. Base DZ con exponentes de Clementi y Roetti [10]

semilla

c1

c2

ε

E

1

0.8305

0.1931

-0.9234

-2.87906

2

0.8417

0.1844

-0.9296

-2.86305

3

0.8446

0.1799

-0.9166

-2.85234

4

0.8376

0.1876

-0.9030

-2.84019

5

0.8604

0.1609

-0.9189

-2.84686

6

0.8277

0.1995

-0.9306

-2.87666

7

0.8441

0.1811

-0.9183

-2.85296

8

0.8847

0.1340

-0.9125

-2.84902

9

0.8374

0.1871

-0.9134

-2.84983

10

0.8146

0.2154

-0.9285

-2.88939

ε = − 0.9194 ± 0.0087 E= − 2.860 ± 0.016 C1 = 0.842 ± 0.019 C2 = 0.182 ± 0.022

Valores exactos: ε = -0.91833 ;

E = -2.86167 ;

C1 = 0.8421 ;

C2 = 0.1827

TABLA 6: Monte Carlo con 105 puntos. Base DZ con exponentes de Clementi y Roetti [10]

semilla

c1

c2

ε

E

1

0.8447

0.1797

-0.9164

-2.8605

2

0.8368

0.1882

-0.9151

-2.8593

3

0.8422

0.1825

-0.9213

-2.8648

4

0.8381

0.1873

-0.9065

-2.8508

5

0.8446

0.1797

-0.9224

-2.8661

6

0.8366

0.1886

-0.9216

-2.8674

7

0.8413

0.1831

-0.9181

-2.8657

8

0.8482

0.1761

-0.9181

-2.8583

9

0.8452

0.1796

-0.9243

-2.8640

10

0.8450

0.1799

-0.9189

-2.8594

ε = − 0.9182 ± 0.0050 E = − 2.8616 ± 0.0050 C1 = 0.8423 ± 0.0040

C2 = 0.1825 ± 0.0043

Resultados similares a los obtenidos en el caso del helio con dos funciones de base, se obtienen al usar integración Monte Carlo en el caso de los iones bielectrónicos. Por ejemplo, en la tabla 7 puede verse el resultado de integrar con 104 y con 106 puntos, en diversos sistemas helioides.

TABLA 7: Monte Carlo con 104 y 106 puntos en sistemas helioides. Base DZ de exponentes α1= Z y α2= Z - 1. Sistema

MC 4 10 puntos

MC 6 10 puntos

Integración Analítica

Límite Hartree-Fock

He

-2.8655

-2.8586

-2.85847

-2.862

Li +

-7.2456

-7.2283

-7.23297

-7.237

Be2+

-13.6274

-13.5992

-13.6062

-13.615

B3+

-22.0134

-21.9696

-21.9798

-21.986

C4+

-32.4149

-32.3403

-32.3538

-32.361

A la vista de los resultados expuestos, puede afirmarse que la integración Monte Carlo es perfectamente aplicable a sistemas tratados con dos funciones de base. Algo parecido se obtiene al emplear tres o incluso cuatro funciones, pero la situación empeora al aumentar más el número de funciones de base. En primer lugar por el aumento del número de integrales a evaluar, pero también por la necesidad de determinar cada una de ellas con mayor precisión que si solo se empleasen dos o tres funciones de base. Un caso representativo de lo que queremos decir acontece en el cálculo del estado fundamental del átomo de helio con base quíntuple-zeta y exponentes de Clementi [11]: α =1.4300 α =2.4415 α =4.0996 1 2 3 α =6.4843 α =0.7978 4 5

(24)

Empleando integrales analíticas se obtienen los resultados de la primera entrada de la tabla 8. Las integrales MC con 104 puntos conducen a un calculo iterativo que no converge con ninguno de los 10 conjuntos de puntos ensayados. Lo mismo ocurre cuando se integra con 105 puntos. La convergencia en la resolución de las ecuaciones de Roothaan solo se alcanza a partir de integraciones con 106 puntos que, a razón de unos 4 minutos de tiempo de cálculo por cada integral (en un PC486 a 33Mhz) , implica unas 8 horas para realizar el cálculo completo. Si bien esto todavía constituye un esfuerzo aceptable, empieza a marcar los límites de aplicación de esta metodología en ordenadores de tipo personal

Resulta interesante comentar una característica de los resultados presentados en la tabla 8: Los obtenidos para la energía total o la energía orbital son mucho mas parecidos a los analíticos, que los obtenidos para los coeficientes del orbital ocupado. De hecho, estos últimos pueden ser muy distintos de los analíticos, como se ve en los de la penúltima columna de la tabla. En nuestra opinión, esto no se debe exclusivamente al empleo de integración Monte Carlo, sino también a una cuasidependencia lineal de las 5 funciones de la base de Clementi, que hace que el resultado final del cálculo Roothaan pueda ser expresado, con un grado de exactitud muy parecido, por combinaciones lineales distintas. Como comprobación, puede verse en la última parte de la tabla 8 que el solapamiento entre las tres funciones obtenidas por Monte Carlo con diferentes conjuntos de puntos y la exacta es bastante parecido. También puede verse que el valor medio de las potencias de la coordenada radial obtenido con la función exacta y las procedentes de Monte Carlo, son razonablemente próximos.

TABLA 8: Monte Carlo con tres conjuntos de 106 puntos. Base quíntuple-zeta con exponentes de Clementi [11]. Al final de la tabla se comparan las integrales de solapamiento función aproximada-función exacta, y los valores medios para n = -1, 1, 2 y 3. parámetro:

valor analítico MC (caso 1)

MC (caso 2)

MC (caso 3)

C1

0.78499

0.7949

1.1525

0.7041

C2

0.20287

0.1880

-0.1817

0.3065

C3

0.03692

0.0479

0.2535

-0.0314

C4

-0.00293

-0.0062

-0.0605

0.0161

C5

0.00326

0.0007

-0.1569

0.0312

ε

-0.91796

-0.9184

-0.9225

-0.9172

E

-2.86168

-2.86103

-2.86098

-2.86149

Sanalítico-caso i

1.00000

1.00009

0.99938

0.99987



1.68730

1.6877

1.6803

1.6863



0.92727

0.9279

0.9268

0.9258



1.18482

1.1864

1.1650

1.1838



1.94057

1.9435

1.8370

1.9546

6: ELECTRONES EN UNA CAJA ESFÉRICA Como acaba de verse, la mayor desventaja de emplear integración Monte Carlo al aplicar el método de Roothaan, se encuentra en la limitación del número de funciones de base que es posible manejar. Su punto fuerte, en cambio, es la posibilidad de elegir libremente el tipo de funciones de base. A título de ejemplo, hemos modificado el programa que realizaba los cálculos en sistemas helioides de manera que lleve a cabo cálculos Roothaan para un par de electrones en una caja esférica, de radio dado R0 , tomando como base funciones de tipo: 



2 −1  nπr  m f = r sen  ⋅ Y (θ, ϕ) n,l,m R R  l 0  0

(25)

y calculando las integrales (f1, f2 | f3, f4) por Monte Carlo con función densidad constante. Los resultados obtenidos con integraciones de 107 puntos (40 min. de cálculo por integral, en un PC486 a 33 Mhz.), para dos y cuatro funciones de base, pueden verse en las tablas 9 y 10. Observese como cambia la

TABLA 9:

Orbital de un par de electrones en una caja esférica, expresado con dos funciones de base de tipo (25) con l=0 y n=1 y 2. Las últimas entradas de la tabla corresponden a la energía despreciando la repulsión interelectrónica, E(0), y a la perturbación de primer orden E(1) = (1,1 | 1,1). R0

0.1

1

10

50

100

C1

1.00000

0.99977

0.98715

0.91535

0.84882

C2

-0.00221

-0.02133

-0.15981

-0.40267

-0.52864

ε

511.32

6.7004

0.21372

0.031669

0.014782

E

1004.81

11.6420

0.26685

0.034603

E(0)

986.96

9.8696

0.09870

0.015690 0.003948 0.987·10-3

E(0) +E(1)

1004.82

11.656

0.2773

0.03968

0.01885

importancia de emplear dos o cuatro funciones de base, según sea grande o pequeño el radio de la esfera: Para radios cortos los resultados de ambos tipos de cálculo son equivalentes. Incluso, si el radio es suficientemente pequeño, pueden obtenerse buenos resultados con una aproximación de perturbaciones de primer orden (final de la tabla 9), en la que solo se precisa calcular por Monte Carlo una integral.

TABLA 10: Orbital de un par de electrones en una caja esférica, expresado con cuatro funciones de base de tipo (25) con l=0 y n=1, 2, 3 y 4. En la última entrada de la tabla se reproduce el resultado obtenido con solo dos funciones de base. R0

0.1

1

10

50

100

C1

1.00000

0.99977

0.98661

0.90479

0.83171

-0.2206E-2 -0.02135

-0.16222

-0.40865

-0.51312

0.01683

C2 C3

0.149E-4

0.3842E-3

0.11617

0.20101

C4

-0.803E-5

-0.831E-4 -0.2079E-2 -0.02929

-0.06758

ε

511.32

6.7004

0.21326

0.030882

0.014140

E 4 funciones

1004.81

11.6419

0.26662

0.034084

0.015217

E 2 funciones

1004.81

11.6420

0.26685

0.034603

0.015690

Referencias: [1] M.H. Kalos y P.A. Whitlock: Monte Carlo Methods WileySons, N.Y. (1986); vol I, cap 4. [2] A. Lüchow y H. Kleindienst, Int. J. Quant. Chem. 45, 445 (1993) [3] M.H. Kalos y P.A. Whitlock, Loc. Cit. sec 2.6 [4] M.H. Kalos y P.A. Whitlock, Loc. Cit. sec 2.9 [5] M.H. Kalos y P.A. Whitlock, Loc. Cit. sec 4.1 [6] J.M. Hammersley y D.C. Handscomb: Monte Carlo Methods Chapman and Hall, , Londres (1964); cap 3 [7] W.H. Press, B.P. Flannery, S.A. Teukolsky y W.T. Vetterling: Numerical Recipes Cambridge University Press (1986); cap 7 [8] I.M. Sóbol: Método de Monte Carlo MIR, Moscú (1976) S. Ríos: Métodos Estadísticos Ediciones del Castillo, Madrid (1977), cap 16 [9] M. Fabbri y L.F. Perondi, Int. J. Quant. Chem. 26, 1029 (1984) [10] E. Clementi y C. Roetti: Atomic Data and Nuclear Data Tables Academic Press, N.Y. 14, 177 (1974) [11] E. Clementi, J. Chem. Phys. 38, 996 (1962)

Agradecimientos: Los autores desean hacer patente su agradecimiento a los profesores Luis Sesé Sánchez, de la Universidad Nacional de Educación a Distancia, Jaime Fernández Rico, de la Universidad Autónoma de Madrid y Antonino García Rendon, de la Universidad de Cádiz, por sus orientaciones en las primeras fases de este trabajo.