2121
2222import os , sys , math
2323import numpy as np
24+ from .mapping import sd2indx
2425
2526# Fisrt, we add the location of the library to test to the PYTHON path
2627if sys .platform == "cygwin" :
2930 from liblibra_core import *
3031
3132
32-
3333def ovlp_arb (_SD1 , _SD2 , S , active_space = None , verbose = False ):
3434 """Compute the overlap of two generic SDs: <SD1|SD2>
35-
35+ This functions translates the explicit use of two Python `for` loops
36+ into a more efficient way using numpy
3637 Args:
3738
3839 _SD1 ( lists of ints ): first SD, such that:
@@ -41,7 +42,7 @@ def ovlp_arb(_SD1, _SD2, S, active_space=None, verbose=False):
4142 _SD2 ( lists of ints ): second SD, such that:
4243 SeeAlso: ```inp``` in the ```sd2indx(inp)``` function
4344
44- S ( CMATRIX (N,N) ): is the matrix in the space of 1-el orbital
45+ S ( numpy array (N,N) ): is the matrix in the space of 1-el orbital
4546
4647 active_space ( list of ints ): indices of the orbitals (starting from 1) to
4748 include into consideration. If None - all the orbitals will be used [ default: None ]
@@ -77,7 +78,95 @@ def ovlp_arb(_SD1, _SD2, S, active_space=None, verbose=False):
7778 # Apply the determinant formula
7879 det_size = len (SD1 )
7980 s = np .zeros ( (det_size , det_size ), dtype = np .float64 );
81+ """ Commented on 7/22/2025
82+ # ============================== The numpy version of the above double-for-loops
83+ # Find the tensor product of the SD1 and SD2
84+ SD_tensor_product = np.tensordot(SD1, SD2, axes=0)
85+ # Next, find the sign of these elementwise multiplications
86+ SD_tensor_product_sign = np.sign(SD_tensor_product)
87+ # Now, find where we have alpha-beta indices so that the
88+ negative_sign_indices = np.where(SD_tensor_product_sign < 0)
89+ s[negative_sign_indices] = 0
90+ """
91+ # Let's build the matrix related to sd1 and sd2 from the KS orbitals
92+ # For this, we reuire to turn each element into matrix indices
93+ sd1 = sd2indx (SD1 , nbasis , False , 0 ) # adopted from the `mapping` module
94+ sd2 = sd2indx (SD2 , nbasis , False , 0 )
95+ # What about beta indices?! We should add `nbasis/2` to them. This is added on 7/22/2025
96+ # With this simple approach, there is no need for tensor product or the `if else` clause brought
97+ # previously since the elements corresponding to the alpha-beta orbitals are already
98+ # zero in the two-spinor format of the KS matrices :)
99+ beta_indices = np .where (np .array (SD1 ) < 0 )
100+ sd1 = np .array (sd1 )
101+ sd1 [beta_indices ] += int (nbasis / 2 )
102+
103+ beta_indices = np .where (np .array (SD2 ) < 0 )
104+ sd2 = np .array (sd2 )
105+ sd2 [beta_indices ] += int (nbasis / 2 )
106+
107+
108+ # Making the numpy array
109+ s = S [sd1 ,:][:,sd2 ]
80110
111+ if verbose == True :
112+ print (s )
113+ res = np .linalg .det (s ) #* phase
114+
115+ return res
116+
117+
118+
119+ def ovlp_arb_2 (_SD1 , _SD2 , S , active_space = None , verbose = False ):
120+ """Compute the overlap of two generic SDs: <SD1|SD2>
121+ The difference of this function with the above one is that it uses explicitly
122+ two Pythonic `for` loops which are not optimized for very large sets
123+
124+ Args:
125+
126+ _SD1 ( lists of ints ): first SD, such that:
127+ SeeAlso: ```inp``` in the ```sd2indx(inp)``` function
128+
129+ _SD2 ( lists of ints ): second SD, such that:
130+ SeeAlso: ```inp``` in the ```sd2indx(inp)``` function
131+
132+ S ( numpy array(N,N) ): is the matrix in the space of 1-el orbital
133+
134+ active_space ( list of ints ): indices of the orbitals (starting from 1) to
135+ include into consideration. If None - all the orbitals will be used [ default: None ]
136+
137+ verbose ( Bool ): whether to print some extra info [ default: False]
138+
139+ Returns:
140+ complex: the overlap of the two determinants <SD1|SD2>
141+
142+ """
143+
144+ nbasis = S .shape [0 ] #num_of_rows
145+
146+ SD1 , SD2 = [], []
147+
148+ if active_space == None :
149+ SD1 , SD2 = _SD1 , _SD2
150+ else :
151+ for a in _SD1 :
152+ if abs (a ) in active_space :
153+ SD1 .append (a )
154+ for a in _SD2 :
155+ if abs (a ) in active_space :
156+ SD2 .append (a )
157+
158+ res = 0.0 + 0j
159+ if len (SD1 ) != len (SD2 ):
160+ print ("Trying to compute an overlap of the SDs with different number of electrons" )
161+ print ("Exiting..." )
162+ sys .exit (0 )
163+
164+
165+ # Apply the determinant formula
166+ det_size = len (SD1 )
167+ s = np .zeros ( (det_size , det_size ), dtype = np .float64 );
168+ # ============================== The explicit version with for loops
169+ # Note: It seems that this explicit vesion doesn't consider the unrestricted spin calculations
81170 # Forming the overlap of the SDs
82171 for indx_I , I in enumerate (SD1 ):
83172 for indx_J , J in enumerate (SD2 ):
0 commit comments