From 961954547aa54dd375af85b3a9687430fbb6c5ad Mon Sep 17 00:00:00 2001 From: yann Date: Mon, 18 May 2020 16:04:18 +0200 Subject: [PATCH] work : calcul indice j avec i --- FullCI.ipynb | 2 +- Work.ipynb | 158 +++++++++++++++++++++------------------------------ 2 files changed, 67 insertions(+), 93 deletions(-) diff --git a/FullCI.ipynb b/FullCI.ipynb index 32623fd..95b2de2 100644 --- a/FullCI.ipynb +++ b/FullCI.ipynb @@ -125,7 +125,7 @@ " \n", "let hf = HartreeFock.make simulation\n", "\n", - "let mo_basis = MOBasis.of_hartree_fock hf" + "let mo_basis = MOBasis.of_hartree_fock hf;;" ] }, { diff --git a/Work.ipynb b/Work.ipynb index e43cb6d..cdba7e2 100644 --- a/Work.ipynb +++ b/Work.ipynb @@ -1299,11 +1299,6 @@ " \n", "let sum a = \n", " Array.fold_left (fun accu x -> accu +. x) 0. a\n", - " \n", - "type alphaij = {\n", - " alpha_max : float;\n", - " indice_ii : int;\n", - " indice_jj : int;};;\n", " \n", "let pi = 3.14159265358979323846264338;;\n" ] @@ -1528,14 +1523,8 @@ "\n", " let alpha_boys , d_boys = f_alpha_boys m_C\n", " in\n", - " (*let d_boys = d_boys m_C\n", - " in*)\n", " let alpha_er , d_er = f_alpha m_C\n", " in\n", - " (*let alpha_er = mat_alpha m_C\n", - " in\n", - " let d_er = s_D m_C\n", - " in*)\n", " let alpha methode =\n", " match methode with \n", " | \"Boys\"\n", @@ -1584,6 +1573,11 @@ "metadata": {}, "outputs": [], "source": [ + "type alphaij = {\n", + " alpha_max : float;\n", + " indice_ii : int;\n", + " indice_jj : int;};;\n", + "\n", "(* Détermination alpha_max et ses indices i et j.\n", "Si alpha max > pi/2 on soustrait pi/2 à la matrice des alphas de manière récursive *)\n", "let rec new_m_alpha m_alpha m_C n_rec_alpha=\n", @@ -1614,23 +1608,13 @@ " |> iamax\n", " in\n", " \n", - " (* indice i du alpha max *)\n", - " let indice_ii = \n", + " (* indice i et j du alpha max *)\n", + " let indice_ii, indice_jj = \n", " let max = max_element3 alpha_m \n", " in\n", - " \n", - " (*Printf.printf \"%i\\n%!\" max;*)\n", - " \n", - " (max - 1) mod n_mo +1 \n", + " (max - 1) mod n_mo +1, (max - 1) / n_mo +1\n", " in\n", - " \n", - " (* indice j du alpha max *)\n", - " let indice_jj = \n", - " let max = max_element3 alpha_m \n", - " in\n", - " (max - 1) / n_mo +1\n", - " in\n", - " \n", + " \n", " (* Valeur du alpha max*)\n", " let alpha alpha_m = \n", " let i = indice_ii \n", @@ -1756,43 +1740,7 @@ "metadata": {}, "outputs": [], "source": [ - "(* Pour la réinjection on créer des matrices intérmédiares, une matrice nulle partout sauf sur \n", - "les colonnes de i et j et de i~ et j~. On fait la différence de la première matrice avec la matrice\n", - "des OMs Phi afin de substituer les colonnes de i et j par des zéro et ensuite sommer cette matrice avec \n", - "celle contenant i~ et j~ *)\n", - " \n", - "(* Matrice intérmédiare pour l'injection de ksi~ (i~ et j~) dans la matrice Phi *)\n", - "(*\n", - "let f_Psi_tilde m_Ksi_tilde indice_i indice_j m_C=\n", - " let n_mo = Mat.dim2 m_C\n", - " in\n", - " let n_ao = Mat.dim1 m_C\n", - "in\n", - "Mat.init_cols n_ao n_mo (fun i j -> \n", - " if j=indice_i \n", - " then m_Ksi_tilde.{i,1}\n", - " else if j=indice_j \n", - " then m_Ksi_tilde.{i,2}\n", - " else 0.)\n", - " \n", - ";;\n", - "*)\n", - "(* Matrice intermédiaire pour supprimer ksi (i et j) dans la matrice Phi *) \n", - "(*\n", - "let f_Psi m_Ksi indice_i indice_j m_C = \n", - " let n_mo = Mat.dim2 m_C\n", - " in\n", - " let n_ao = Mat.dim1 m_C\n", - "in\n", - "Mat.init_cols n_ao n_mo (fun i j -> \n", - " if j=indice_i \n", - " then m_Ksi.{i,1}\n", - " else if j=indice_j \n", - " then m_Ksi.{i,2}\n", - " else 0.)\n", - ";;\n", - "*)\n", - "\n", + "(* Fonction pour la création de matrice intermédiaire pour supprimer i j et ajouter i~ et j~ *)\n", "let f_k mat indice_i indice_j m_C = \n", "let n_mo = Mat.dim2 m_C\n", " in\n", @@ -1805,10 +1753,8 @@ " then mat.{i,2}\n", " else 0.)\n", "(*************************)\n", - "\n", "(*\n", "let m_Psi = f_Psi m_Ksi indice_i indice_j;; \n", - "\n", "let m_Psi_tilde = f_Psi_tilde m_Ksi_tilde indice_i indice_j;; \n", "*)" ] @@ -1824,12 +1770,9 @@ "let f_interm m_C m_Psi = Mat.sub m_C m_Psi;;\n", "\n", "(*************************)\n", - "\n", "(*\n", "let m_interm = f_interm m_C m_Psi;; \n", - "\n", "let new_m_C m_C= Mat.add m_Psi_tilde m_interm;;\n", - "\n", "let m_new_m_C = new_m_C m_C;;\n", "*)" ] @@ -1893,8 +1836,8 @@ " (* Matrice des alphas *)\n", " let m_alpha = alphad.m_alpha\n", " in\n", - " (*let norme_alpha = norme m_alpha\n", - " in*)\n", + " let norme_alpha = norme m_alpha\n", + " in\n", " \n", " (*Printf.printf \"%f\\n%!\" norme_alpha;*)\n", " \n", @@ -1957,12 +1900,12 @@ " (*Util.debug_matrix \"m_interm\" m_interm;*)\n", " \n", " (* Matrice après rotation *)\n", - " ( Mat.add m_Psi_tilde m_interm, critere_D)\n", + " ( Mat.add m_Psi_tilde m_interm, critere_D, norme_alpha)\n", " in\n", - " let m_new_m_C , critere_D = new_m_C m_C methode loc_deloc \n", + " let m_new_m_C , critere_D, norme_alpha = new_m_C m_C methode loc_deloc \n", " in\n", - " let diff = prev_critere_D -. critere_D +. 1.\n", - " \n", + " let diff = prev_critere_D -. critere_D \n", + "\n", "in\n", "\n", "(*Util.debug_matrix \"new_alpha_m\" (f_alpha m_C);*)\n", @@ -1971,7 +1914,6 @@ "if diff**2. < cc**2.\n", " then m_new_m_C\n", " else\n", - "\n", "localisation m_new_m_C methode loc_deloc epsilon (n-1) critere_D cc;;" ] }, @@ -2066,32 +2008,64 @@ "metadata": {}, "outputs": [], "source": [ + "let new_m_boys = localisation m_C \"boys\" \"loc\" 1. 100 0. 10e-7;;\n", + "let new_m_er = localisation m_C \"ER\" \"loc\" 1. 100 0. 10e-7;;\n", + "\n", "let loc_m_occ_boys = localisation m_occ \"boys\" \"loc\" 1. 100 0. 10e-7;;\n", "let loc_m_occ_er = localisation m_occ \"ER\" \"loc\" 1. 100 0. 10e-7;;\n", "let loc_m_vir_boys = localisation m_vir \"boys\" \"loc\" 1. 100 0. 10e-7;;\n", - "let loc_m_vir_er = localisation m_vir \"ER\" \"loc\" 1. 100 0. 10e-7;;" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ + "let loc_m_vir_er = localisation m_vir \"ER\" \"loc\" 1. 100 0. 10e-7;;\n", + "\n", "let m_assemble_boys = assemble loc_m_occ_boys loc_m_vir_boys;;\n", - "let m_assemble_er = assemble loc_m_occ_er loc_m_vir_er;;" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ + "let m_assemble_er = assemble loc_m_occ_er loc_m_vir_er;;\n", + "\n", "let loc_assemble_boys = localisation m_assemble_boys \"boys\" \"loc\" 1. 100 0. 10e-7;;\n", "let loc_assemble_er = localisation m_assemble_er \"er\" \"loc\" 1. 100 0. 10e-7;;" ] }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "let nuclei =\n", + " Nuclei.of_xyz_string xyz_string\n", + " \n", + "let basis = \n", + " Basis.of_nuclei_and_basis_string nuclei basis_string\n", + " \n", + "let simulation = \n", + " Simulation.make ~charge:0 ~multiplicity:1 ~nuclei basis\n", + " \n", + "let ao_basis = \n", + " Simulation.ao_basis simulation\n", + " \n", + "let hf = HartreeFock.make simulation\n", + "\n", + "let mo_basis = MOBasis.of_hartree_fock hf" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "let ci =\n", + " DeterminantSpace.fci_of_mo_basis mo_basis ~frozen_core:false\n", + " |> CI.make\n", + " \n", + "let ci_coef, ci_energy = Lazy.force ci.eigensystem " + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, { "cell_type": "code", "execution_count": null,