@@ -3295,80 +3295,81 @@ def symmetric_low_rank_spectrum_update(update, orig_e, U, proj_U, force_rank_inc
32953295    #return the new eigenvalues 
32963296    return  new_evals , True 
32973297
3298- #Note: This function won't work for our purposes because of the assumptions 
3299- #about the rank of the update on the nullspace of the matrix we're updating, 
3300- #but keeping this here commented for future reference. 
3301- #Function for doing fast calculation of the updated inverse trace: 
3302- #def riedel_style_inverse_trace(update, orig_e, U, proj_U, force_rank_increase=True): 
3303- #    """ 
3304- #    input: 
3305- #     
3306- #    update : ndarray 
3307- #        symmetric low-rank update to perform. 
3308- #        This is the first half the symmetric rank decomposition s.t. 
3309- #        update@update.T= the full update matrix. 
3310- #     
3311- #    orig_e : ndarray 
3312- #        Spectrum of the original matrix. This is a 1-D array. 
3313- #         
3314- #    proj_U : ndarray 
3315- #        Projector onto the complement of the column space of the 
3316- #        original matrix's eigenvectors. 
3317- #         
3318- #    output: 
3319- #     
3320- #    trace : float 
3321- #        Value of the trace of the updated psuedoinverse matrix. 
3322- #     
3323- #    updated_rank : int 
3324- #        total rank of the updated matrix. 
3325- #         
3326- #    rank_increase_flag : bool 
3327- #        a flag that is returned to indicate is a candidate germ failed to amplify additional parameters.  
3328- #        This indicates things short circuited and so the scoring function should skip this germ. 
3329- #    """ 
3330- #     
3331- #    #First we need to for the matrix P, whose column space 
3332- #    #forms an orthonormal basis for the component of update 
3333- #    #that is in the complement of U. 
3334- #     
3335- #    proj_update= proj_U@update 
3336- #     
3337- #    #Next take the RRQR decomposition of this matrix: 
3338- #    q_update, r_update, _ = _sla.qr(proj_update, mode='economic', pivoting=True) 
3339- #     
3340- #    #Construct P by taking the columns of q_update corresponding to non-zero values of r_A on the diagonal. 
3341- #    nonzero_indices_update= _np.nonzero(_np.diag(r_update)>1e-10) #HARDCODED (threshold is hardcoded) 
3342- #     
3343- #    #if the rank doesn't increase then we can't use the Riedel approach. 
3344- #    #Abort early and return a flag to indicate the rank did not increase. 
3345- #    if len(nonzero_indices_update[0])==0 and force_rank_increase: 
3346- #        return None, None, False 
3347- #     
3348- #    P= q_update[: , nonzero_indices_update[0]] 
3349- #     
3350- #    updated_rank= len(orig_e)+ len(nonzero_indices_update[0]) 
3351- #     
3352- #    #Now form the matrix R_update which is given by P.T @ proj_update. 
3353- #    R_update= P.T@proj_update 
3354- #     
3355- #    #R_update gets concatenated with U.T@update to form 
3356- #    #a block column matrixblock_column= np.concatenate([U.T@update, R_update], axis=0)     
3357- #     
3358- #    Uta= U.T@update 
3359- #     
3360- #    try: 
3361- #        RRRDinv= R_update@_np.linalg.inv(R_update.T@R_update)  
3362- #    except _np.linalg.LinAlgError as err: 
3363- #        print('Numpy thinks this matrix is singular, condition number is: ', _np.linalg.cond(R_update.T@R_update)) 
3364- #        print((R_update.T@R_update).shape) 
3365- #        raise err 
3366- #    pinv_orig_e_mat= _np.diag(1/orig_e) 
3367- #     
3368- #    trace= _np.sum(1/orig_e) + _np.trace( RRRDinv@(_np.eye(Uta.shape[1]) + Uta.T@pinv_orig_e_mat@Uta)@RRRDinv.T ) 
3369- #     
3370- #    return trace, updated_rank, True 
3298+ # Note: Th function below won't work for our purposes because of the assumptions 
3299+ # about the rank of the update on the nullspace of the matrix we're updating, 
3300+ # but keeping this here commented for future reference. 
3301+ ''' 
3302+ def riedel_style_inverse_trace(update, orig_e, U, proj_U, force_rank_increase=True): 
3303+     """ 
3304+     input: 
3305+      
3306+     update : ndarray 
3307+         symmetric low-rank update to perform. 
3308+         This is the first half the symmetric rank decomposition s.t. 
3309+         update@update.T= the full update matrix. 
3310+      
3311+     orig_e : ndarray 
3312+         Spectrum of the original matrix. This is a 1-D array. 
3313+          
3314+     proj_U : ndarray 
3315+         Projector onto the complement of the column space of the 
3316+         original matrix's eigenvectors. 
3317+          
3318+     output: 
3319+      
3320+     trace : float 
3321+         Value of the trace of the updated psuedoinverse matrix. 
3322+      
3323+     updated_rank : int 
3324+         total rank of the updated matrix. 
3325+          
3326+     rank_increase_flag : bool 
3327+         a flag that is returned to indicate is a candidate germ failed to amplify additional parameters.  
3328+         This indicates things short circuited and so the scoring function should skip this germ. 
3329+     """ 
33713330     
3331+     #First we need to for the matrix P, whose column space 
3332+     #forms an orthonormal basis for the component of update 
3333+     #that is in the complement of U. 
3334+      
3335+     proj_update= proj_U@update 
3336+      
3337+     #Next take the RRQR decomposition of this matrix: 
3338+     q_update, r_update, _ = _sla.qr(proj_update, mode='economic', pivoting=True) 
3339+      
3340+     #Construct P by taking the columns of q_update corresponding to non-zero values of r_A on the diagonal. 
3341+     nonzero_indices_update= _np.nonzero(_np.diag(r_update)>1e-10) #HARDCODED (threshold is hardcoded) 
3342+      
3343+     #if the rank doesn't increase then we can't use the Riedel approach. 
3344+     #Abort early and return a flag to indicate the rank did not increase. 
3345+     if len(nonzero_indices_update[0])==0 and force_rank_increase: 
3346+         return None, None, False 
3347+      
3348+     P= q_update[: , nonzero_indices_update[0]] 
3349+      
3350+     updated_rank= len(orig_e)+ len(nonzero_indices_update[0]) 
3351+      
3352+     #Now form the matrix R_update which is given by P.T @ proj_update. 
3353+     R_update= P.T@proj_update 
3354+      
3355+     #R_update gets concatenated with U.T@update to form 
3356+     #a block column matrixblock_column= np.concatenate([U.T@update, R_update], axis=0)     
3357+      
3358+     Uta= U.T@update 
3359+      
3360+     try: 
3361+         RRRDinv= R_update@_np.linalg.inv(R_update.T@R_update)  
3362+     except _np.linalg.LinAlgError as err: 
3363+         print('Numpy thinks this matrix is singular, condition number is: ', _np.linalg.cond(R_update.T@R_update)) 
3364+         print((R_update.T@R_update).shape) 
3365+         raise err 
3366+     pinv_orig_e_mat= _np.diag(1/orig_e) 
3367+      
3368+     trace= _np.sum(1/orig_e) + _np.trace( RRRDinv@(_np.eye(Uta.shape[1]) + Uta.T@pinv_orig_e_mat@Uta)@RRRDinv.T ) 
3369+      
3370+     return trace, updated_rank, True 
3371+ ''' 
3372+ 
33723373def  minamide_style_inverse_trace (update , orig_e , U , proj_U , force_rank_increase = False ):
33733374    """ 
33743375    This function performs a low-rank update to the components of 
0 commit comments