Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Adding RMSD #38

Open
wants to merge 4 commits into
base: master
Choose a base branch
from
Open

Adding RMSD #38

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
e7e2fe6
Adding RMSD
AlexandreDeCamargo Apr 9, 2024
da68838
Merge branch 'theochem:master' into alex_PyCI
AlexandreDeCamargo May 1, 2024
29738a6
Merge branch 'theochem:master' into alex_PyCI
AlexandreDeCamargo May 3, 2024
9c819e0
Adding Neumman and Dkystra algorithms.
AlexandreDeCamargo May 3, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
111 changes: 96 additions & 15 deletions notebooks/fanci_tutorial.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -51,7 +51,7 @@
},
{
"cell_type": "code",
"execution_count": 1,
"execution_count": 68,
"metadata": {},
"outputs": [],
"source": [
Expand All @@ -63,9 +63,12 @@
"from pyci.test import datafile\n",
"\n",
"# System information\n",
"filename = datafile(\"lih_sto6g.fcidump\")\n",
"# filename = datafile(\"lih_sto6g.fcidump\")\n",
"filename = datafile(\"/workspaces/PyCI/pyci/test/data/h4_2_0.0_STO-6G.fcidump\")\n",
"\n",
"ham = pyci.hamiltonian(filename)\n",
"e_dict = {}"
"e_dict = {}\n",
"f_dict = {}"
]
},
{
Expand All @@ -79,7 +82,7 @@
},
{
"cell_type": "code",
"execution_count": 2,
"execution_count": 59,
"metadata": {},
"outputs": [],
"source": [
Expand All @@ -101,7 +104,7 @@
},
{
"cell_type": "code",
"execution_count": 3,
"execution_count": 60,
"metadata": {},
"outputs": [],
"source": [
Expand All @@ -124,13 +127,28 @@
},
{
"cell_type": "code",
"execution_count": 4,
"execution_count": 61,
"metadata": {},
"outputs": [],
"source": [
"# Optimize wavefunction\n",
"ap1rog_results = ap1rog.optimize(ap1_params, use_jac=True)\n",
"e_dict[\"AP1roG\"] = ap1rog_results.x[-1]"
"e_dict[\"AP1roG\"] = (ap1rog_results.x[-1])\n",
"\n",
"# f_dict['AP1roG'] = (np.std(ap1rog_results.fun)) #Calculating the std of the residuals \n",
"\n",
"#Calculation of root mean square deviation (RMSD) \n",
"\"\"\" \n",
"RMSD = \\sqrt ( \\sum_i^n (residuals)^2 / (n - 2)) \n",
"n = size of sample \n",
"resilduals = Vector of residuals at the solution.\n",
"\"\"\"\n",
"f_dict['AP1roG'] = (np.sqrt(np.sum((ap1rog_results.fun)**2)/(len(ap1rog_results.fun)-2)))\n",
"\n",
"# print(f_dict['AP1roG'])\n",
"\n",
"# print(ap1rog_results) \n",
"\n"
]
},
{
Expand All @@ -146,7 +164,7 @@
},
{
"cell_type": "code",
"execution_count": 5,
"execution_count": 62,
"metadata": {},
"outputs": [],
"source": [
Expand All @@ -160,7 +178,22 @@
"\n",
"# Optimize wavefunction\n",
"pccds_results = pccds.optimize(pccds_params, use_jac=False)\n",
"e_dict[\"pCCD+S\"] = pccds_results.x[-1]"
"e_dict[\"pCCD+S\"] = pccds_results.x[-1]\n",
"# f_dict[\"pCCD+S\"] = (np.std(ap1rog_results.fun))\n",
"\n",
"#Calculation of root mean square deviation (RMSD) \n",
"\"\"\" \n",
"RMSD = \\sqrt ( \\sum_i^n (residuals)^2 / (n - 2)) \n",
"n = size of sample \n",
"resilduals = Vector of residuals at the solution.\n",
"\"\"\"\n",
"\n",
"f_dict['pCCD+S'] = (np.sqrt(np.sum((pccds_results.fun)**2)/(len(pccds_results.fun)-2)))\n",
"\n",
"\n",
"# print(pccds_results)\n",
"# print(f_dict['pCCD+S'])\n",
"\n"
]
},
{
Expand All @@ -174,17 +207,17 @@
},
{
"cell_type": "code",
"execution_count": 6,
"execution_count": 63,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
" METHOD, ENERGY\n",
" HF, -8.947289175e+00\n",
" AP1roG, -8.963531034e+00\n",
" pCCD+S, -8.963613544e+00\n"
" HF, -3.875079433e+00\n",
" AP1roG, -3.955402403e+00\n",
" pCCD+S, -3.816787845e+00\n"
]
}
],
Expand All @@ -195,6 +228,54 @@
" print(f\"{name:>8s}, {energy:>16.9e}\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Calculate the standard deviation "
]
},
{
"cell_type": "code",
"execution_count": 66,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"2.4683556322835916e-14\n",
"1.9667674453325034e-09\n"
]
}
],
"source": [
"print(f_dict[\"AP1roG\"])\n",
"print(f_dict[\"pCCD+S\"])\n"
]
},
{
"cell_type": "code",
"execution_count": 65,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
" METHOD, RMSD\n",
" AP1roG, 2.468355632e-14\n",
" pCCD+S, 1.966767445e-09\n"
]
}
],
"source": [
"# Print RMSD from various methods\n",
"print(f\"{'METHOD':>8s}, {'RMSD':>16s}\")\n",
"for name, fun in f_dict.items():\n",
" print(f\"{name:>8s}, {fun:>16.9e}\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
Expand Down Expand Up @@ -242,7 +323,7 @@
"output_type": "stream",
"text": [
"Number of electrons = 4.1\n",
"Number of pairs = 6.1\n"
"Number of pairs = 6.2\n"
]
}
],
Expand All @@ -269,7 +350,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.6"
"version": "3.10.13"
}
},
"nbformat": 4,
Expand Down
65 changes: 65 additions & 0 deletions pyci/rdm/algorithms.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,65 @@
import numpy as np
from abc import ABCMeta, abstractmethod
class Projection(metaclass=ABCMeta):

@abstractmethod
def __init__(self, initial_guess,constraints):
self.initial_guess = initial_guess
self.constraints = constraints

@abstractmethod
def optimize(self):
pass

class Dykstra(metaclass=ABCMeta):
def __init__(self, initial_guess:float, constraints:list, alpha:int =1, max_iterations:int =100, eps:int =1e-6):
"""
Dykstra's algorithm for projection onto convex sets.

Parameters
----------
initial_guess : float
Initializing with guessed density matrix \Gamma_0.
constraints : list
A list of projection operators onto all $J$ of the constraints as a single set, $\{\mathcal{P}_j\}_{j=0}^{J-1}$.
alpha : int, optional
Tunning parameter, by default 1
max_iterations : int, optional
Number of maximum iterations, by default 100
eps : int, optional
Tolerance, by default 1e-6
"""
super().__init__(initial_guess, constraints)
self.alpha = alpha
self.max_iterations = max_iterations
self.eps = eps
def optimize(self):
X = [np.zeros(self.initial_guess) for i in range(len(self.constraints))]
D = np.copy(self.initial_guess)
norm = []
for i in range (self.max_iterations):
for j, projection in enumerate(self.constraints):
C = D - X[j]
D = projection(C)
X[j] = D - C
norm.append(np.linalg.norm(D - self.initial_guess))
is_stop = self.alpha * abs(norm[i] - norm[i - 1]) + (1 - self.alpha) * norm[i] < self.eps
if is_stop:
break
class Neumann(metaclass=ABCMeta):
def __init__(self, initial_guess, constraints, alpha, max_iterations=100, eps=1e-6):
super().__init__(initial_guess, constraints)
self.alpha = alpha
self.max_iterations = max_iterations
self.eps = eps
def optimize(self):
D = np.copy(self.initial_guess)
norm = []
for i in range (self.max_iterations):
for projection in enumerate(self.constraints):
D = projection(D)
norm.append(np.linalg.norm(D - self.initial_guess))
is_stop = self.alpha * abs(norm[i] - norm[i - 1]) + (1 - self.alpha) * norm[i] < self.eps
if is_stop:
break

115 changes: 115 additions & 0 deletions pyci/test/data/h4_2_0.0_STO-6G.fcidump
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,115 @@
&FCI NORB= 4,NELEC= 4,MS2=0,
ORBSYM=1,1,1,1,
ISYM=1,
&END
0.5068898550551018 1 1 1 1
0.003794527559320846 1 1 2 1
0.4098635263371225 1 1 2 2
1.025916213093225e-09 1 1 3 1
-8.035957177465747e-10 1 1 3 2
0.445100515768005 1 1 3 3
0.09719871204645167 1 1 4 1
0.05096181960669316 1 1 4 2
-2.511664404236669e-09 1 1 4 3
0.5183452230109675 1 1 4 4
0.003794527559320839 2 1 1 1
0.08813159971644449 2 1 2 1
0.07782243723669213 2 1 2 2
-5.653522719417658e-10 2 1 3 1
-4.283259379123106e-09 2 1 3 2
-0.08764576192895299 2 1 3 3
0.02826446097635334 2 1 4 1
-0.08001497578900814 2 1 4 2
7.330390461302017e-11 2 1 4 3
0.05589849747028398 2 1 4 4
0.4098635263371225 2 2 1 1
0.07782243723669217 2 2 2 1
0.521293056868916 2 2 2 2
-5.046337625547181e-09 2 2 3 1
1.404773519730895e-11 2 2 3 2
0.3561543663891216 2 2 3 3
-0.001427950323337701 2 2 4 1
-0.09114783741363394 2 2 4 2
3.876907371525018e-09 2 2 4 3
0.4713483506720944 2 2 4 4
1.025916174929309e-09 3 1 1 1
-5.653522754112128e-10 3 1 2 1
-5.046337653302757e-09 3 1 2 2
0.1293028829330219 3 1 3 1
-0.1072466481774506 3 1 3 2
4.05264671821115e-09 3 1 3 3
-1.581443448417374e-09 3 1 4 1
7.217758335453794e-11 3 1 4 2
-0.0846859118356684 3 1 4 3
-1.893284105880522e-09 3 1 4 4
-8.035956726437643e-10 3 2 1 1
-4.283259382592552e-09 3 2 2 1
1.404779764735409e-11 3 2 2 2
-0.1072466481774506 3 2 3 1
0.09458993631478847 3 2 3 2
1.40467162901281e-09 3 2 3 3
7.391295561709477e-11 3 2 4 1
4.802033475720391e-09 3 2 4 2
0.07788008022122561 3 2 4 3
-1.244524983068374e-09 3 2 4 4
0.445100515768005 3 3 1 1
-0.08764576192895296 3 3 2 1
0.3561543663891216 3 3 2 2
4.052646773722302e-09 3 3 3 1
1.404671601257235e-09 3 3 3 2
0.5359374140556439 3 3 3 3
-0.006072330211121087 3 3 4 1
0.1003613149338521 3 3 4 2
-3.081953497963141e-09 3 3 4 3
0.4053549128731515 3 3 4 4
0.09719871204645164 4 1 1 1
0.02826446097635335 4 1 2 1
-0.001427950323337657 4 1 2 2
-1.581443447550013e-09 4 1 3 1
7.391295756865868e-11 4 1 3 2
-0.006072330211121068 4 1 3 3
0.1186180412469477 4 1 4 1
0.03612665777791915 4 1 4 2
-1.119209735306326e-09 4 1 4 3
0.122452034039052 4 1 4 4
0.05096181960669319 4 2 1 1
-0.08001497578900814 4 2 2 1
-0.09114783741363391 4 2 2 2
7.217762845734832e-11 4 2 3 1
4.802033409800899e-09 4 2 3 2
0.1003613149338522 4 2 3 3
0.03612665777791913 4 2 4 1
0.1272246486516095 4 2 4 2
-1.005679005350313e-09 4 2 4 3
0.004098379798186072 4 2 4 4
-2.511664473625608e-09 4 3 1 1
7.330390461302017e-11 4 3 2 1
3.876907420097275e-09 4 3 2 2
-0.08468591183566843 4 3 3 1
0.07788008022122564 4 3 3 2
-3.081953497963141e-09 4 3 3 3
-1.119209738775773e-09 4 3 4 1
-1.005678960247502e-09 4 3 4 2
0.06715293748280152 4 3 4 3
-8.629244367674982e-10 4 3 4 4
0.5183452230109674 4 4 1 1
0.05589849747028395 4 4 2 1
0.4713483506720945 4 4 2 2
-1.89328408506384e-09 4 4 3 1
-1.244525038579525e-09 4 4 3 2
0.4053549128731514 4 4 3 3
0.1224520340390518 4 4 4 1
0.004098379798186058 4 4 4 2
-8.629244610536269e-10 4 4 4 3
0.5934456368866865 4 4 4 4
-1.80985796887928 1 1 0 0
-0.08161693543269405 2 1 0 0
-1.373368656657632 2 2 0 0
4.184586471194904e-09 3 1 0 0
5.930225890093923e-10 3 2 0 0
-1.366056716783876 3 3 0 0
-0.174357790674463 4 1 0 0
0.01748866944921174 4 2 0 0
5.075449520787198e-10 4 3 0 0
-1.115558063357834 4 4 0 0
2.168608539661471 0 0 0 0
Loading