Work : boucle ER, début reflexion boys
This commit is contained in:
parent
82540d1a7c
commit
1268019bc2
317
Work.ipynb
317
Work.ipynb
@ -890,6 +890,80 @@
|
||||
"## Rotation de deux orbitales (Boys)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"On part de la matrice des OMs donné par $HF$ -> m_C.On calcul $A_{ij}^{x,y,z}$, et $B_{ij}^{x,y,z}$, tel que :\n",
|
||||
"\n",
|
||||
"$A_{ij}^x = \\langle \\phi_i | x^2 | \\phi_j \\rangle$ \n",
|
||||
"\n",
|
||||
"$= \\sum_a \\sum_b c_{ai} c_{bj} \\langle \\chi_a | x^2 | \\chi_b \\rangle$\n",
|
||||
"\n",
|
||||
"$B_{ij}^x = \\langle \\phi_i | x | \\phi_j \\rangle$ \n",
|
||||
"\n",
|
||||
"$= \\sum_a \\sum_b c_{ai} c_{bj} \\langle \\chi_a | x | \\chi_b \\rangle$\n",
|
||||
"\n",
|
||||
"On forme ainsi des matrices (n par n) $A_{ij}^x$, $A_{ij}^y$, $A_{ij}^z$, $B_{ij}^x$, $B_{ij}^y$ et $B_{ij}^z$ que l'on peut sommer pour obtenir :\n",
|
||||
"\n",
|
||||
"$A_{ij} =A_{ij}^x + A_{ij}^y + A_{ij}^z$ et $B_{ij} =B_{ij}^x + B_{ij}^y + B_{ij}^z$\n",
|
||||
"\n",
|
||||
"$D$ est défini comme :\n",
|
||||
"\n",
|
||||
"$D(\\theta) = D(0) + \\frac{1}{4} [(1-cos4\\theta)\\beta_{ij}+sin4\\theta \\gamma_{ij}] $\n",
|
||||
"\n",
|
||||
"Où\n",
|
||||
"\n",
|
||||
"$\\beta_{ij}= (B^x_{ii}-B^x_{jj})^2 - 4 {(B^x_{ij})}^2 + [...y...] + [...z...]$\n",
|
||||
"\n",
|
||||
"$=(B_{ii}-B_{jj})^2 - 4 {(B_{ij})}^2 $\n",
|
||||
"\n",
|
||||
"Et \n",
|
||||
"\n",
|
||||
"$\\gamma_{ij}= 4 B^x_{ij} (B^{x}_{ii}-B^x_{jj}) + [...y...] + [...z...]$\n",
|
||||
"\n",
|
||||
"$=4 B_{ij} (B_{ii}-B_{jj})$\n",
|
||||
"\n",
|
||||
"Avec\n",
|
||||
"\n",
|
||||
"$D(0)= D^{x}(0) + D^{y}(0) + D^{z}(0)$\n",
|
||||
"\n",
|
||||
"$=A_{ii}^x + A_{jj}^x - (\\tilde B_{ii}^x)^2 - (\\tilde B_{jj}^x)^2$\n",
|
||||
"$+A_{ii}^y + A_{jj}^y - (\\tilde B_{ii}^y)^2 - (\\tilde B_{jj}^y)^2$\n",
|
||||
"$+A_{ii}^z + A_{jj}^z - (\\tilde B_{ii}^z)^2 - (\\tilde B_{jj}^z)^2$\n",
|
||||
"\n",
|
||||
"$= A_{ii} + A_{jj} - \\tilde B_{ii}^2 - \\tilde B_{jj}^2$\n",
|
||||
"\n",
|
||||
"Avec\n",
|
||||
"\n",
|
||||
"${\\tilde B}^{x2}_{ii}= (cos^2\\theta B^x_{ii} + sin^2\\theta B^x_{jj} - sin2\\theta B^x_{ij})^2$\n",
|
||||
"\n",
|
||||
"${\\tilde B}^{x2}_{jj}= (sin^2\\theta B^x_{ii} + cos^2\\theta B^x_{jj} + sin2\\theta B^x_{ij})^2$\n",
|
||||
"\n",
|
||||
"Et donc :\n",
|
||||
"\n",
|
||||
"${\\tilde B}^{2}_{ii}= (cos^2\\theta B_{ii} + sin^2\\theta B_{jj} - sin2\\theta B_{ij})^2$\n",
|
||||
"\n",
|
||||
"${\\tilde B}^{2}_{jj}= (sin^2\\theta B_{ii} + cos^2\\theta B_{jj} + sin2\\theta B_{ij})^2$\n",
|
||||
"\n",
|
||||
"Ainsi, on a : \n",
|
||||
"\n",
|
||||
"$D(\\theta) = A_{ii} + A_{jj} - (cos^2\\theta B_{ii} + sin^2\\theta B_{jj} - sin2\\theta B_{ij})^2 - (sin^2\\theta B_{ii} + cos^2\\theta B_{jj} + sin2\\theta B_{ij})^2$\n",
|
||||
"\n",
|
||||
"Le but est de minimiser la matrice $D(\\theta)$ il faut donc faire une rotation des orbitales i et j donnants le plus grand élément de la matrice $D(\\theta)$\n",
|
||||
"\n",
|
||||
"Si je comprends bien, on doit prendre le plus grand élément $D(i,j)$ de la matrice $D(\\theta)$ pour $\\theta = 0$ au départ et d'appliquer des rotation en faisant varier $\\theta$ pour minimiser $D(i,j)$ puis recommencer avec le nouveau $D(i,j)max$\n",
|
||||
"\n",
|
||||
"\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": []
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
@ -928,8 +1002,6 @@
|
||||
"$$\n",
|
||||
"{\\cal B}_x(\\theta)= A^x_{11}+A^x_{22} \n",
|
||||
"- [(cos^4\\theta+ sin^4\\theta)({B^x_{11}}^2+ {B^x_{22}}^2 )\n",
|
||||
"$$\n",
|
||||
"$$\n",
|
||||
"+ 2 sin^2 2\\theta {B^x_{12}}^2\n",
|
||||
"+ 2 sin 2\\theta cos 2\\theta (({B^x_{22}} -{B^x_{11}} ) {B^x_{12}}]\n",
|
||||
"$$\n",
|
||||
@ -941,8 +1013,6 @@
|
||||
"$$\n",
|
||||
"{\\cal B}_x(\\theta)= A^x_{11}+A^x_{22} \n",
|
||||
"- [ \\frac{1}{4} ( 3 + cos4\\theta)({B^x_{11}}^2+ {B^x_{22}}^2 )\n",
|
||||
"$$\n",
|
||||
"$$\n",
|
||||
"+ (1 -cos 4\\theta) {B^x_{12}}^2\n",
|
||||
"+ sin 4\\theta (({B^x_{22}} -{B^x_{11}} ) {B^x_{12}}]\n",
|
||||
"$$\n",
|
||||
@ -1073,7 +1143,30 @@
|
||||
"\n",
|
||||
"$= \\left(\\sum_a c_{ai} \\sum_b c_{bi} - \\sum_a c_{aj} \\sum_b c_{bj} \\right) [ \\chi_a \\chi_b| \\phi_i \\phi_j ] $\n",
|
||||
"\n",
|
||||
"$= \\left(\\sum_a c_{ai} \\sum_b c_{bi} - \\sum_a c_{aj} \\sum_b c_{bj} \\right) \\sum_e \\sum_f c_{ei} c_{fj} [ \\chi_a \\chi_b| \\chi_e \\chi_f ] $"
|
||||
"$= \\left(\\sum_a c_{ai} \\sum_b c_{bi} - \\sum_a c_{aj} \\sum_b c_{bj} \\right) \\sum_e \\sum_f c_{ei} c_{fj} [ \\chi_a \\chi_b| \\chi_e \\chi_f ] $\n",
|
||||
"\n",
|
||||
"On extrait de la matrice des OMs les orbitales i et j (vecteurs colonnes) pour former une nouvelle matrice m_Ksi de dimension n * 2. \n",
|
||||
"\n",
|
||||
"Pour effectuer la rotation des orbitales i et j on utilise la matrice de rotation m_R pour la rotation de 2 orbitales, qui est définie comme :\n",
|
||||
"\n",
|
||||
"m_R =\n",
|
||||
"\n",
|
||||
" ( cos(alpha) -sin(alpha) )\n",
|
||||
"\n",
|
||||
" ( sin(alpha) cos(alpha) )\n",
|
||||
" \n",
|
||||
"On applique m_R à m_Ksi : m_R m_Ksi = m_Ksi_thilde \n",
|
||||
"\n",
|
||||
"On obtient m_Ksi_tilde la matrice contenant les deux nouvelles OMs i~ et j~ obtenues par rotation de i et j.\n",
|
||||
"\n",
|
||||
"On réinjecte ces deux nouvelles orbitales i~ et j~ à la place des anciennes orbitales i et j dans la matrice des OMs m_Phi, ce qui nous donne une nouvelle matrice des OMs m_Phi_tilde. Pour cela on créer des matrices intérmédiaires:\n",
|
||||
" - une matrice ( m_Psi ) n par n avec tous ses éléments nuls sauf les colonnes qui contiennent les OMs i et j \n",
|
||||
" - une matrice ( m_Psi_tilde ) n par n avec tous ses éléments nuls sauf les colonnes qui contiennent les OMs i~ et j~\n",
|
||||
" - une matrice ( m_interm ) n par n où l'on a soustrait m_Psi à m_Phi pour créer des 0 sur les colonnes des OMs i et j \n",
|
||||
"\n",
|
||||
"Ainsi, on soustrait à la matrice des OMs Phi la matrice m_ksi pour supprimer les OMs i et j de celle ci puis on additionne cette nouvelle matrice à la matrice m_ksi_tilde pour créer la nouvelle matrice des OMs m_Phi_tilde.\n",
|
||||
"\n",
|
||||
"On repart de cette nouvelle matrice m_Phi_tilde et on cherche la paire d'orbitale (k,l) ayant le plus grand angle de rotation alpha. Et on procède comme nous l'avons précedemment de manière intérative. Le but étant de maximiser D, c'est à dire maximiser le terme de répulsion.\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
@ -1086,13 +1179,18 @@
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"metadata": {
|
||||
"jupyter": {
|
||||
"outputs_hidden": true
|
||||
}
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"(* Localisation de Edminstion *)\n",
|
||||
"let ee_ints = AOBasis.ee_ints ao_basis;;\n",
|
||||
"let m_C = MOBasis.mo_coef mo_basis;;\n",
|
||||
"let n_ao = Mat.dim1 m_C ;;\n",
|
||||
"let n_mo = n_ao - 1;;\n",
|
||||
"let sum a = \n",
|
||||
" Array.fold_left (fun accu x -> accu +. x) 0. a\n",
|
||||
";;\n",
|
||||
@ -1119,18 +1217,22 @@
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"metadata": {
|
||||
"jupyter": {
|
||||
"outputs_hidden": true
|
||||
}
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"(* Calcul de toutes les intégrales B_12 -> Matrice *)\n",
|
||||
"let m_b12 = Mat.init_cols n_ao n_ao (fun i j -> \n",
|
||||
"let m_b12 = Mat.init_cols n_mo n_mo (fun i j -> \n",
|
||||
"integral_general (fun a b e f i j ->\n",
|
||||
" ( m_C.{a,i} *. m_C.{b,i} -. m_C.{a,j} *. m_C.{b,j} ) *. m_C.{e,i} *. m_C.{f,j}\n",
|
||||
" ) i j\n",
|
||||
" )\n",
|
||||
"\n",
|
||||
"(* Calcul de toutes les intégrales A_12 -> Matrice *)\n",
|
||||
"let m_a12 = Mat.init_cols n_ao n_ao (fun i j ->\n",
|
||||
"let m_a12 = Mat.init_cols n_mo n_mo (fun i j ->\n",
|
||||
"integral_general (fun a b e f i j -> m_C.{a,i} *. m_C.{b,i} *. m_C.{e,i} *. m_C.{f,j} \n",
|
||||
" -. 0.25 *. ( m_C.{e,i} *. m_C.{f,i} -. m_C.{e,j} *. m_C.{f,j} ) \n",
|
||||
" *. ( m_C.{a,i} *. m_C.{b,i} -. m_C.{a,j} *. m_C.{b,j} )\n",
|
||||
@ -1138,8 +1240,8 @@
|
||||
" )\n",
|
||||
"\n",
|
||||
"(* Calcul de tous les alpha -> Matrice *)\n",
|
||||
"let m_alpha = Mat.init_cols n_ao n_ao ( fun i j ->\n",
|
||||
" if i<>j then asin(m_b12.{i,j} /. sqrt((m_a12.{i,j})**2 + (m_b12.{i,j})**2 )\n",
|
||||
"let m_alpha = Mat.init_cols n_mo n_mo ( fun i j ->\n",
|
||||
" asin(m_b12.{i,j} /. sqrt((m_a12.{i,j}**2.) +. (m_b12.{i,j}**2.) )\n",
|
||||
" )\n",
|
||||
" )"
|
||||
]
|
||||
@ -1147,7 +1249,11 @@
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"metadata": {
|
||||
"jupyter": {
|
||||
"outputs_hidden": true
|
||||
}
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"(* Extraction du plus grand angle alpha*)\n",
|
||||
@ -1155,8 +1261,11 @@
|
||||
"|> iamax;;\n",
|
||||
"\n",
|
||||
"(* Indices du plus grand angle alpha *)\n",
|
||||
"let indice_i = ij_max_alpha mod n_ao;;\n",
|
||||
"let indice_j = ij_max_alpha / n_ao + 1;;"
|
||||
"let indice_i = ij_max_alpha mod n_mo;;\n",
|
||||
"let indice_j = ij_max_alpha / n_mo + 1;;\n",
|
||||
"\n",
|
||||
"(* Valeur du plus grand angle alpha*)\n",
|
||||
"let alpha = m_alpha.{indice_i,indice_j};;"
|
||||
]
|
||||
},
|
||||
{
|
||||
@ -1169,15 +1278,19 @@
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"metadata": {
|
||||
"jupyter": {
|
||||
"outputs_hidden": true
|
||||
}
|
||||
},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"let m_Phi = m_C\n",
|
||||
"\n",
|
||||
"(* Matrice de rotation 2 par 2 *)\n",
|
||||
"let m_R = Mat.init_cols 2 2 (fun i j -> if i=j then cos a\n",
|
||||
" else if i>j then sin a \n",
|
||||
" else -. sin a );;\n",
|
||||
"let m_R = Mat.init_cols 2 2 (fun i j -> if i=j then cos alpha\n",
|
||||
" else if i>j then sin alpha \n",
|
||||
" else -. sin alpha );;\n",
|
||||
"\n",
|
||||
"(* Fonction d'extraction des 2 vecteurs propres i et j de la matrice des OMs pour les mettres dans la matrice Ksi (n par 2)\n",
|
||||
"pour appliquer R afin d'effectuer la rotation des orbitales *) (* {1,2} -> 1ere ligne, 2e colonne *)\n",
|
||||
@ -1185,9 +1298,9 @@
|
||||
"let q = indice_j;;\n",
|
||||
"let m_Ksi = Mat.init_cols n 2 (fun i j -> if j=1 then m_Phi.{i,p} else m_Phi.{i,q} );;\n",
|
||||
"\n",
|
||||
"(* Fonction de calcul de ksi~ (matrice n par 2), nouvelle matrice par application de la matrice de rotation dans laquelle on obtient\n",
|
||||
"les deux orbitales que l'on va réinjecter dans la matrice Phi*)\n",
|
||||
"let m_Ksi_tilde = gemm ksi m_R;;"
|
||||
"(* Fonction de calcul de ksi~ (matrice n par 2), nouvelle matrice par application de la matrice de rotation dans laquelle on\n",
|
||||
"obtient les deux orbitales que l'on va réinjecter dans la matrice Phi*)\n",
|
||||
"let m_Ksi_tilde = gemm m_Ksi m_R;;"
|
||||
]
|
||||
},
|
||||
{
|
||||
@ -1209,32 +1322,42 @@
|
||||
"celle contenant i~ et j~ *)\n",
|
||||
"\n",
|
||||
"(* Matrice intérmédiare pour l'injection de ksi~ (i~ et j~) dans la matrice Phi *)\n",
|
||||
"let m_Psi_tilde = Mat.init_cols n n (fun i j -> if j=q then m_Ksi_tilde.{i,p}\n",
|
||||
"let m_Psi_tilde = Mat.init_cols n_mo n_mo (fun i j -> if j=q then m_Ksi_tilde.{i,p}\n",
|
||||
" else if j=p then m_Ksi_tilde.{i,q}\n",
|
||||
" else 0.);;\n",
|
||||
" \n",
|
||||
"(* Matrice intermédiaire pour supprimer ksi (i et j) dans la matrice Phi *) \n",
|
||||
"let m_Psi = Mat.init_cols n n (fun i j -> if j=q then m_Ksi.{i,p}\n",
|
||||
"let m_Psi = Mat.init_cols n_mo n_mo (fun i j -> if j=q then m_Ksi.{i,p}\n",
|
||||
" else if j=p then m_Ksi.{i,q}\n",
|
||||
" else 0.);;\n",
|
||||
" \n",
|
||||
"(* Matrice intérmédiaire où les orbitales i et j ont été supprimées et remplacées par des 0, par soustraction de la matrice Phi\n",
|
||||
"par la matrice *)\n",
|
||||
"let m_interm = Mat.sub m_C M_Psi;;\n",
|
||||
"let m_interm = Mat.sub m_C m_Psi;;\n",
|
||||
"\n",
|
||||
"(* Construction de la nouvelle matrice d'OMs phi~ *)\n",
|
||||
"let m_Phi_tilde = Mat.add m_Psi_tilde m_interm;; "
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"# TESTS"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": []
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"(* Fonction pour itérer *) (*cc -> critère de convergence*)\n",
|
||||
"\n",
|
||||
"while m_alpha.{indice_i,indice_j} > cc do m_C = m_Phi_tilde\n"
|
||||
"## Réécriture pour itérerer ER"
|
||||
]
|
||||
},
|
||||
{
|
||||
@ -1243,6 +1366,132 @@
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"(* Fonction pour itérer , cc critère de convergence*) \n",
|
||||
"(* Localisation de Edminstion *)\n",
|
||||
"let ee_ints = AOBasis.ee_ints ao_basis;;\n",
|
||||
"let m_C = MOBasis.mo_coef mo_basis;;\n",
|
||||
"let n_ao = Mat.dim1 m_C ;;\n",
|
||||
"let n_mo = n_ao - 1;;\n",
|
||||
"let sum a = \n",
|
||||
" Array.fold_left (fun accu x -> accu +. x) 0. a\n",
|
||||
";;\n",
|
||||
"\n",
|
||||
"(*Calcul des intégrales*) \n",
|
||||
"let integral_general g i j =\n",
|
||||
"Array.map (fun a ->\n",
|
||||
" let v = \n",
|
||||
" Array.map (fun b ->\n",
|
||||
" let u = \n",
|
||||
" Array.map (fun e ->\n",
|
||||
" let t = Array.map (fun f ->\n",
|
||||
" (g a b e f i j) *. ERI.get_phys ee_ints a e b f\n",
|
||||
" ) (Util.array_range 1 n_ao)\n",
|
||||
" in sum t\n",
|
||||
" ) (Util.array_range 1 n_ao)\n",
|
||||
" in sum u\n",
|
||||
" ) (Util.array_range 1 n_ao)\n",
|
||||
" in sum v\n",
|
||||
") (Util.array_range 1 n_ao)\n",
|
||||
"|> sum \n",
|
||||
"\n",
|
||||
"\n",
|
||||
"let m_Phi = m_C;;\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"(* Calcul de la nouvelle matrice Phi après rotations *)\n",
|
||||
"let rec new_Phi m_Phi n = \n",
|
||||
"\n",
|
||||
" (* Calcul de tous les alpha -> Matrice *)\n",
|
||||
" let m_alpha m_C =\n",
|
||||
"\n",
|
||||
" (* Calcul de toutes les intégrales B_12 -> Matrice *)\n",
|
||||
" let m_b12 = Mat.init_cols n_mo n_mo (fun i j -> \n",
|
||||
" integral_general (fun a b e f i j ->\n",
|
||||
" ( m_C.{a,i} *. m_C.{b,i} -. m_C.{a,j} *. m_C.{b,j} ) *. m_C.{e,i} *. m_C.{f,j}\n",
|
||||
" ) i j\n",
|
||||
" )\n",
|
||||
" in\n",
|
||||
" (* Calcul de toutes les intégrales A_12 -> Matrice *)\n",
|
||||
" let m_a12 = Mat.init_cols n_mo n_mo (fun i j ->\n",
|
||||
" integral_general (fun a b e f i j -> m_C.{a,i} *. m_C.{b,i} *. m_C.{e,i} *. m_C.{f,j} \n",
|
||||
" -. 0.25 *. ( m_C.{e,i} *. m_C.{f,i} -. m_C.{e,j} *. m_C.{f,j} ) \n",
|
||||
" *. ( m_C.{a,i} *. m_C.{b,i} -. m_C.{a,j} *. m_C.{b,j} )\n",
|
||||
" ) i j\n",
|
||||
" )\n",
|
||||
" in\n",
|
||||
" Mat.init_cols n_mo n_mo ( fun i j ->\n",
|
||||
" asin(m_b12.{i,j} /. sqrt((m_a12.{i,j}**2.) +. (m_b12.{i,j}**2.) )\n",
|
||||
" )\n",
|
||||
" )\n",
|
||||
"\n",
|
||||
" in\n",
|
||||
" (* Matrice de rotation 2 par 2 *)\n",
|
||||
" let m_R m_alpha = \n",
|
||||
"\n",
|
||||
" (* Valeur du plus grand angle alpha*)\n",
|
||||
" let alpha = \n",
|
||||
"\n",
|
||||
" (* Extraction du plus grand angle alpha*)\n",
|
||||
" let ij_max_alpha = Mat.as_vec m_alpha\n",
|
||||
" |> iamax\n",
|
||||
" in\n",
|
||||
"\n",
|
||||
" (* Indices du plus grand angle alpha *)\n",
|
||||
" let indice_i = ij_max_alpha mod n_mo\n",
|
||||
" in\n",
|
||||
" let indice_j = ij_max_alpha / n_mo + 1\n",
|
||||
"\n",
|
||||
" in\n",
|
||||
" m_alpha.{indice_i,indice_j}\n",
|
||||
" in\n",
|
||||
" Printf.printf \"%f\\n%!\" alpha;\n",
|
||||
"\n",
|
||||
" in\n",
|
||||
" Mat.init_cols 2 2 (fun i j -> if i=j then cos alpha\n",
|
||||
" else if i>j then sin alpha \n",
|
||||
" else -. sin alpha )\n",
|
||||
"\n",
|
||||
" in\n",
|
||||
" (* Fonction de calcul de ksi~ (matrice n par 2), nouvelle matrice par application de la matrice de rotation dans laquelle \n",
|
||||
" on obtient les deux orbitales que l'on va réinjecter dans la matrice Phi*)\n",
|
||||
" let m_Ksi_tilde m_R = \n",
|
||||
"\n",
|
||||
" (* Fonction d'extraction des 2 vecteurs propres i et j de la matrice des OMs pour les mettres dans la matrice Ksi \n",
|
||||
" (n par 2) pour appliquer R afin d'effectuer la rotation des orbitales *) (* {1,2} -> 1ere ligne, 2e colonne *)\n",
|
||||
" let p = indice_i\n",
|
||||
" in\n",
|
||||
" let q = indice_j\n",
|
||||
" in\n",
|
||||
" let m_Ksi = Mat.init_cols n 2 (fun i j -> if j=1 then m_Phi.{i,p} else m_Phi.{i,q} )\n",
|
||||
" in\n",
|
||||
" gemm m_Ksi m_R\n",
|
||||
"\n",
|
||||
" in\n",
|
||||
" (* Construction de la nouvelle matrice d'OMs phi~ *)\n",
|
||||
" let m_Phi_tilde m_Ksi m_Ksi_tilde = \n",
|
||||
"\n",
|
||||
" (* Matrice intérmédiaire où les orbitales i et j ont été supprimées et remplacées par des 0, par soustraction de la \n",
|
||||
" matrice Phi par la matrice *)\n",
|
||||
" let m_interm = \n",
|
||||
"\n",
|
||||
" (* Matrice intérmédiare pour l'injection de ksi~ (i~ et j~) dans la matrice Phi *)\n",
|
||||
" let m_Psi_tilde = Mat.init_cols n_mo n_mo (fun i j -> if j=q then m_Ksi_tilde.{i,p}\n",
|
||||
" else if j=p then m_Ksi_tilde.{i,q}\n",
|
||||
" else 0.)\n",
|
||||
" in\n",
|
||||
" (* Matrice intermédiaire pour supprimer ksi (i et j) dans la matrice Phi *) \n",
|
||||
" let m_Psi = Mat.init_cols n_mo n_mo (fun i j -> if j=q then m_Ksi.{i,p}\n",
|
||||
" else if j=p then m_Ksi.{i,q}\n",
|
||||
" else 0.)\n",
|
||||
"\n",
|
||||
" in\n",
|
||||
" Mat.sub m_C m_Psi\n",
|
||||
" in\n",
|
||||
" Mat.add m_Psi_tilde m_interm\n",
|
||||
"in\n",
|
||||
"if n > 0 then new_Phi m_Phi_tilde (n-1);;\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"\n"
|
||||
]
|
||||
},
|
||||
@ -1253,6 +1502,20 @@
|
||||
"outputs": [],
|
||||
"source": []
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## Boys"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": []
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
|
||||
Loading…
Reference in New Issue
Block a user