Python numpy.sort() Examples
The following are 30
code examples of numpy.sort().
You can vote up the ones you like or vote down the ones you don't like,
and go to the original project or source file by following the links above each example.
You may also want to check out all available functions/classes of the module
numpy
, or try the search function
.

Example #1
Source File: matcher_test.py From object_detector_app with MIT License | 6 votes |
def test_all_columns_accounted_for(self): # Note: deliberately setting to small number so not always # all possibilities appear (matched, unmatched, ignored) num_matches = 10 match_results = tf.random_uniform( [num_matches], minval=-2, maxval=5, dtype=tf.int32) match = matcher.Match(match_results) matched_column_indices = match.matched_column_indices() unmatched_column_indices = match.unmatched_column_indices() ignored_column_indices = match.ignored_column_indices() with self.test_session() as sess: matched, unmatched, ignored = sess.run([ matched_column_indices, unmatched_column_indices, ignored_column_indices ]) all_indices = np.hstack((matched, unmatched, ignored)) all_indices_sorted = np.sort(all_indices) self.assertAllEqual(all_indices_sorted, np.arange(num_matches, dtype=np.int32))
Example #2
Source File: gtf_utils.py From models with MIT License | 6 votes |
def add_exon(self, chrom, strand, start, stop): if strand != self.strand or chrom != self.chrom: print("The exon has different chrom or strand to the transcript.") return _exon = np.array([start, stop], "int").reshape(1, 2) self.exons = np.append(self.exons, _exon, axis=0) self.exons = np.sort(self.exons, axis=0) self.tranL += abs(int(stop) - int(start) + 1) self.exonNum += 1 self.seglen = np.zeros(self.exons.shape[0] * 2 - 1, "int") self.seglen[0] = self.exons[0, 1] - self.exons[0, 0] + 1 for i in range(1, self.exons.shape[0]): self.seglen[i * 2 - 1] = self.exons[i, 0] - self.exons[i - 1, 1] - 1 self.seglen[i * 2] = self.exons[i, 1] - self.exons[i, 0] + 1 if ["-", "-1", "0", 0, -1].count(self.strand) > 0: self.seglen = self.seglen[::-1]
Example #3
Source File: dataloader.py From models with MIT License | 6 votes |
def add_exon(self, chrom, strand, start, stop): if strand != self.strand or chrom != self.chrom: print("The exon has different chrom or strand to the transcript.") return _exon = np.array([start, stop], "int").reshape(1,2) self.exons = np.append(self.exons, _exon, axis=0) self.exons = np.sort(self.exons, axis=0) self.tranL += abs(int(stop) - int(start) + 1) self.exonNum += 1 self.seglen = np.zeros(self.exons.shape[0] * 2 - 1, "int") self.seglen[0] = self.exons[0,1]-self.exons[0,0] + 1 for i in range(1, self.exons.shape[0]): self.seglen[i*2-1] = self.exons[i,0]-self.exons[i-1,1] - 1 self.seglen[i*2] = self.exons[i,1]-self.exons[i,0] + 1 if ["-","-1","0",0,-1].count(self.strand) > 0: self.seglen = self.seglen[::-1]
Example #4
Source File: gtf_utils.py From models with MIT License | 6 votes |
def add_exon(self, chrom, strand, start, stop): if strand != self.strand or chrom != self.chrom: print("The exon has different chrom or strand to the transcript.") return _exon = np.array([start, stop], "int").reshape(1, 2) self.exons = np.append(self.exons, _exon, axis=0) self.exons = np.sort(self.exons, axis=0) self.tranL += abs(int(stop) - int(start) + 1) self.exonNum += 1 self.seglen = np.zeros(self.exons.shape[0] * 2 - 1, "int") self.seglen[0] = self.exons[0, 1] - self.exons[0, 0] + 1 for i in range(1, self.exons.shape[0]): self.seglen[i * 2 - 1] = self.exons[i, 0] - self.exons[i - 1, 1] - 1 self.seglen[i * 2] = self.exons[i, 1] - self.exons[i, 0] + 1 if ["-", "-1", "0", 0, -1].count(self.strand) > 0: self.seglen = self.seglen[::-1]
Example #5
Source File: dataloader_m.py From models with MIT License | 6 votes |
def map_cpg_tables(cpg_tables, chromo, chromo_pos): """Maps values from cpg_tables to `chromo_pos`. Positions in `cpg_tables` for `chromo` must be a subset of `chromo_pos`. Inserts `dat.CPG_NAN` for uncovered positions. """ chromo_pos.sort() mapped_tables = OrderedDict() for name, cpg_table in six.iteritems(cpg_tables): cpg_table = cpg_table.loc[cpg_table.chromo == chromo] cpg_table = cpg_table.sort_values('pos') mapped_table = map_values(cpg_table.value.values, cpg_table.pos.values, chromo_pos) assert len(mapped_table) == len(chromo_pos) mapped_tables[name] = mapped_table return mapped_tables
Example #6
Source File: test_PlanetPopulation.py From EXOSIMS with BSD 3-Clause "New" or "Revised" License | 6 votes |
def test_honor_arange(self): """ Tests that the input range for semi-major axis is properly set. """ exclude_setrange = ['EarthTwinHabZone1', 'EarthTwinHabZone2', 'JupiterTwin', 'AlbedoByRadiusDulzPlavchan', 'DulzPlavchan', 'EarthTwinHabZone1SDET', 'EarthTwinHabZone3', 'EarthTwinHabZoneSDET'] arangein = np.sort(np.random.rand(2)*10.0) for mod in self.allmods: if mod.__name__ not in exclude_setrange: with RedirectStreams(stdout=self.dev_null): obj = mod(arange = arangein,**self.spec) #test that the correct input arange is used self.assertTrue(obj.arange[0].value == arangein[0],'sma low bound set failed for %s'%mod.__name__) self.assertTrue(obj.arange[1].value == arangein[1],'sma high bound set failed for %s'%mod.__name__) # test that an incorrect input arange is corrected with RedirectStreams(stdout=self.dev_null): obj = mod(arange=arangein[::-1], **self.spec) self.assertTrue(obj.arange[0].value == arangein.min(), 'sma low bound set failed for %s' % mod.__name__) self.assertTrue(obj.arange[1].value == arangein.max(), 'sma high bound set failed for %s' % mod.__name__)
Example #7
Source File: matcher_test.py From DOTA_models with Apache License 2.0 | 6 votes |
def test_all_columns_accounted_for(self): # Note: deliberately setting to small number so not always # all possibilities appear (matched, unmatched, ignored) num_matches = 10 match_results = tf.random_uniform( [num_matches], minval=-2, maxval=5, dtype=tf.int32) match = matcher.Match(match_results) matched_column_indices = match.matched_column_indices() unmatched_column_indices = match.unmatched_column_indices() ignored_column_indices = match.ignored_column_indices() with self.test_session() as sess: matched, unmatched, ignored = sess.run([ matched_column_indices, unmatched_column_indices, ignored_column_indices ]) all_indices = np.hstack((matched, unmatched, ignored)) all_indices_sorted = np.sort(all_indices) self.assertAllEqual(all_indices_sorted, np.arange(num_matches, dtype=np.int32))
Example #8
Source File: box_list_ops_test.py From DOTA_models with Apache License 2.0 | 6 votes |
def test_convert_to_absolute_and_back(self): coordinates = np.random.uniform(size=(100, 4)) coordinates = np.sort(coordinates) coordinates[99, :] = [0, 0, 1, 1] img = tf.ones((128, 202, 202, 3)) boxlist = box_list.BoxList(tf.constant(coordinates, tf.float32)) boxlist = box_list_ops.to_absolute_coordinates(boxlist, tf.shape(img)[1], tf.shape(img)[2]) boxlist = box_list_ops.to_normalized_coordinates(boxlist, tf.shape(img)[1], tf.shape(img)[2]) with self.test_session() as sess: out = sess.run(boxlist.get()) self.assertAllClose(out, coordinates)
Example #9
Source File: box_list_ops_test.py From DOTA_models with Apache License 2.0 | 6 votes |
def test_convert_to_normalized_and_back(self): coordinates = np.random.uniform(size=(100, 4)) coordinates = np.round(np.sort(coordinates) * 200) coordinates[:, 2:4] += 1 coordinates[99, :] = [0, 0, 201, 201] img = tf.ones((128, 202, 202, 3)) boxlist = box_list.BoxList(tf.constant(coordinates, tf.float32)) boxlist = box_list_ops.to_normalized_coordinates(boxlist, tf.shape(img)[1], tf.shape(img)[2]) boxlist = box_list_ops.to_absolute_coordinates(boxlist, tf.shape(img)[1], tf.shape(img)[2]) with self.test_session() as sess: out = sess.run(boxlist.get()) self.assertAllClose(out, coordinates)
Example #10
Source File: gtf_utils.py From models with MIT License | 6 votes |
def add_exon(self, chrom, strand, start, stop): if strand != self.strand or chrom != self.chrom: print("The exon has different chrom or strand to the transcript.") return _exon = np.array([start, stop], "int").reshape(1, 2) self.exons = np.append(self.exons, _exon, axis=0) self.exons = np.sort(self.exons, axis=0) self.tranL += abs(int(stop) - int(start) + 1) self.exonNum += 1 self.seglen = np.zeros(self.exons.shape[0] * 2 - 1, "int") self.seglen[0] = self.exons[0, 1] - self.exons[0, 0] + 1 for i in range(1, self.exons.shape[0]): self.seglen[i * 2 - 1] = self.exons[i, 0] - self.exons[i - 1, 1] - 1 self.seglen[i * 2] = self.exons[i, 1] - self.exons[i, 0] + 1 if ["-", "-1", "0", 0, -1].count(self.strand) > 0: self.seglen = self.seglen[::-1]
Example #11
Source File: dataloader_m.py From models with MIT License | 6 votes |
def map_values(values, pos, target_pos, dtype=None, nan=dat.CPG_NAN): """Maps `values` array at positions `pos` to `target_pos`. Inserts `nan` for uncovered positions. """ assert len(values) == len(pos) assert np.all(pos == np.sort(pos)) assert np.all(target_pos == np.sort(target_pos)) values = values.ravel() pos = pos.ravel() target_pos = target_pos.ravel() idx = np.in1d(pos, target_pos) pos = pos[idx] values = values[idx] if not dtype: dtype = values.dtype target_values = np.empty(len(target_pos), dtype=dtype) target_values.fill(nan) idx = np.in1d(target_pos, pos).nonzero()[0] assert len(idx) == len(values) assert np.all(target_pos[idx] == pos) target_values[idx] = values return target_values
Example #12
Source File: tedana.py From me-ica with GNU Lesser General Public License v2.1 | 6 votes |
def scoreatpercentile(a, per, limit=(), interpolation_method='lower'): """ This function is grabbed from scipy """ values = np.sort(a, axis=0) if limit: values = values[(limit[0] <= values) & (values <= limit[1])] idx = per /100. * (values.shape[0] - 1) if (idx % 1 == 0): score = values[int(idx)] else: if interpolation_method == 'fraction': score = _interpolate(values[int(idx)], values[int(idx) + 1], idx % 1) elif interpolation_method == 'lower': score = values[int(np.floor(idx))] elif interpolation_method == 'higher': score = values[int(np.ceil(idx))] else: raise ValueError("interpolation_method can only be 'fraction', " \ "'lower' or 'higher'") return score
Example #13
Source File: qrdqn.py From tensorflow_RL with MIT License | 6 votes |
def train_model(self): minibatch = random.sample(self.memory, self.batch_size) state_stack = [mini[0] for mini in minibatch] next_state_stack = [mini[1] for mini in minibatch] action_stack = [mini[2] for mini in minibatch] reward_stack = [mini[3] for mini in minibatch] done_stack = [mini[4] for mini in minibatch] done_stack = [int(i) for i in done_stack] onehotaction = np.zeros([self.batch_size, self.output_size]) for i, j in zip(onehotaction, action_stack): i[j] = 1 action_stack = np.stack(onehotaction) Q_next_state = self.sess.run(self.target_network, feed_dict={self.targetNet.input: next_state_stack}) next_action = np.argmax(np.mean(Q_next_state, axis=2), axis=1) Q_next_state_next_action = [Q_next_state[i, action, :] for i, action in enumerate(next_action)] Q_next_state_next_action = np.sort(Q_next_state_next_action) T_theta = [np.ones(self.num_support) * reward if done else reward + self.gamma * Q for reward, Q, done in zip(reward_stack, Q_next_state_next_action, done_stack)] _, l = self.sess.run([self.train_op, self.loss], feed_dict={self.mainNet.input: state_stack, self.action: action_stack, self.Y: T_theta}) return l
Example #14
Source File: box_list_ops_test.py From object_detector_app with MIT License | 6 votes |
def test_convert_to_normalized_and_back(self): coordinates = np.random.uniform(size=(100, 4)) coordinates = np.round(np.sort(coordinates) * 200) coordinates[:, 2:4] += 1 coordinates[99, :] = [0, 0, 201, 201] img = tf.ones((128, 202, 202, 3)) boxlist = box_list.BoxList(tf.constant(coordinates, tf.float32)) boxlist = box_list_ops.to_normalized_coordinates(boxlist, tf.shape(img)[1], tf.shape(img)[2]) boxlist = box_list_ops.to_absolute_coordinates(boxlist, tf.shape(img)[1], tf.shape(img)[2]) with self.test_session() as sess: out = sess.run(boxlist.get()) self.assertAllClose(out, coordinates)
Example #15
Source File: box_list_ops_test.py From object_detector_app with MIT License | 6 votes |
def test_convert_to_absolute_and_back(self): coordinates = np.random.uniform(size=(100, 4)) coordinates = np.sort(coordinates) coordinates[99, :] = [0, 0, 1, 1] img = tf.ones((128, 202, 202, 3)) boxlist = box_list.BoxList(tf.constant(coordinates, tf.float32)) boxlist = box_list_ops.to_absolute_coordinates(boxlist, tf.shape(img)[1], tf.shape(img)[2]) boxlist = box_list_ops.to_normalized_coordinates(boxlist, tf.shape(img)[1], tf.shape(img)[2]) with self.test_session() as sess: out = sess.run(boxlist.get()) self.assertAllClose(out, coordinates)
Example #16
Source File: test_Completeness.py From EXOSIMS with BSD 3-Clause "New" or "Revised" License | 6 votes |
def test_revise_update(self): """ Ensure that target completeness update revisions are appropriately sized """ for mod in self.allmods: if 'revise_updates' not in mod.__dict__: continue with RedirectStreams(stdout=self.dev_null): obj = mod(**copy.deepcopy(self.spec)) obj.gen_update(self.TL) ind = np.sort(np.random.choice(self.TL.nStars,size=int(self.TL.nStars/2.0),replace=False)) obj.revise_updates(ind) self.assertTrue(obj.updates.shape == (len(ind),5),"Updates array improperly resized for %s"%mod.__name__)
Example #17
Source File: test_loss_functions.py From risk-slim with BSD 3-Clause "New" or "Revised" License | 6 votes |
def get_score_bounds_from_range(Z_min, Z_max, rho_lb, rho_ub, L0_max = None): "global variables: L0_reg_ind" edge_values = np.vstack([Z_min * rho_lb, Z_max * rho_lb, Z_min * rho_ub, Z_max * rho_ub]) if L0_max is None or L0_max == Z_min.shape[0]: s_min = np.sum(np.min(edge_values, axis = 0)) s_max = np.sum(np.max(edge_values, axis = 0)) else: min_values = np.min(edge_values, axis = 0) s_min_reg = np.sum(np.sort(min_values[L0_reg_ind])[0:L0_max]) s_min_no_reg = np.sum(min_values[~L0_reg_ind]) s_min = s_min_reg + s_min_no_reg max_values = np.max(edge_values, axis = 0) s_max_reg = np.sum(-np.sort(-max_values[L0_reg_ind])[0:L0_max]) s_max_no_reg = np.sum(max_values[~L0_reg_ind]) s_max = s_max_reg + s_max_no_reg return s_min, s_max #generate data
Example #18
Source File: setup_functions.py From risk-slim with BSD 3-Clause "New" or "Revised" License | 6 votes |
def get_score_bounds(Z_min, Z_max, rho_lb, rho_ub, L0_reg_ind = None, L0_max = None): edge_values = np.vstack([Z_min * rho_lb, Z_max * rho_lb, Z_min * rho_ub, Z_max * rho_ub]) if (L0_max is None) or (L0_reg_ind is None) or (L0_max == Z_min.shape[0]): s_min = np.sum(np.min(edge_values, axis=0)) s_max = np.sum(np.max(edge_values, axis=0)) else: min_values = np.min(edge_values, axis=0) s_min_reg = np.sum(np.sort(min_values[L0_reg_ind])[0:L0_max]) s_min_no_reg = np.sum(min_values[~L0_reg_ind]) s_min = s_min_reg + s_min_no_reg max_values = np.max(edge_values, axis=0) s_max_reg = np.sum(-np.sort(-max_values[L0_reg_ind])[0:L0_max]) s_max_no_reg = np.sum(max_values[~L0_reg_ind]) s_max = s_max_reg + s_max_no_reg return s_min, s_max
Example #19
Source File: reference_direction.py From pymoo with Apache License 2.0 | 6 votes |
def map_onto_unit_simplex(rnd, method): n_points, n_dim = rnd.shape if method == "sum": ret = rnd / rnd.sum(axis=1)[:, None] elif method == "kraemer": M = sys.maxsize rnd *= M rnd = rnd[:, :n_dim - 1] rnd = np.column_stack([np.zeros(n_points), rnd, np.full(n_points, M)]) rnd = np.sort(rnd, axis=1) ret = np.full((n_points, n_dim), np.nan) for i in range(1, n_dim + 1): ret[:, i - 1] = rnd[:, i] - rnd[:, i - 1] ret /= M else: raise Exception("Invalid unit simplex mapping!") return ret
Example #20
Source File: nanoplotter_main.py From NanoPlot with GNU General Public License v3.0 | 6 votes |
def yield_by_minimal_length_plot(array, name, path, title=None, color="#4CB391", figformat="png"): df = pd.DataFrame(data={"lengths": np.sort(array)[::-1]}) df["cumyield_gb"] = df["lengths"].cumsum() / 10**9 yield_by_length = Plot( path=path + "Yield_By_Length." + figformat, title="Yield by length") ax = sns.regplot( x='lengths', y="cumyield_gb", data=df, x_ci=None, fit_reg=False, color=color, scatter_kws={"s": 3}) ax.set( xlabel='Read length', ylabel='Cumulative yield for minimal length', title=title or yield_by_length.title) yield_by_length.fig = ax.get_figure() yield_by_length.save(format=figformat) plt.close("all") return yield_by_length
Example #21
Source File: test_order_crossover.py From pymoo with Apache License 2.0 | 6 votes |
def order_crossover_contributed_no_shift(x1, x2, seq=None): assert len(x1) == len(x2) if seq is not None: start, end = seq else: start, end = np.sort(np.random.choice(len(x1), 2, replace=False)) y1 = x1.copy() y2 = x2.copy() # build y1 and y2 segment1 = set(y1[start:end]) segment2 = set(y2[start:end]) I = np.concatenate((np.arange(0, start), np.arange(end, len(x1)))) # find elements in x2 that are not in segment1 y1[I] = [y for y in x2 if y not in segment1] # find elements in x1 that are not in segment2 y2[I] = [y for y in x1 if y not in segment2] return y1, y2
Example #22
Source File: signal_fixpeaks.py From NeuroKit with MIT License | 6 votes |
def _correct_misaligned(misaligned_idcs, peaks): corrected_peaks = peaks.copy() misaligned_idcs = np.array(misaligned_idcs) # Make sure to not generate negative indices, or indices that exceed # the total number of peaks. prev_peaks and next_peaks must have the # same number of elements. valid_idcs = np.logical_and( misaligned_idcs > 1, misaligned_idcs < len(corrected_peaks) - 1 # pylint: disable=E1111 ) misaligned_idcs = misaligned_idcs[valid_idcs] prev_peaks = corrected_peaks[[i - 1 for i in misaligned_idcs]] next_peaks = corrected_peaks[[i + 1 for i in misaligned_idcs]] half_ibi = (next_peaks - prev_peaks) / 2 peaks_interp = prev_peaks + half_ibi # Shift the R-peaks from the old to the new position. corrected_peaks = np.delete(corrected_peaks, misaligned_idcs) corrected_peaks = np.concatenate((corrected_peaks, peaks_interp)).astype(int) corrected_peaks.sort(kind="mergesort") return corrected_peaks
Example #23
Source File: post_proc.py From HorizonNet with MIT License | 6 votes |
def get_rot_rad(init_coorx, coory, z=50, coorW=1024, coorH=512, floorW=1024, floorH=512, tol=5): gpid = get_gpid(init_coorx, coorW) coor = np.hstack([np.arange(coorW)[:, None], coory[:, None]]) xy = np_coor2xy(coor, z, coorW, coorH, floorW, floorH) xy_cor = [] rot_rad_suggestions = [] for j in range(len(init_coorx)): pca = PCA(n_components=1) pca.fit(xy[gpid == j]) rot_rad_suggestions.append(_get_rot_rad(*pca.components_[0])) rot_rad_suggestions = np.sort(rot_rad_suggestions + [1e9]) rot_rad = np.mean(rot_rad_suggestions[:-1]) best_rot_rad_sz = -1 last_j = 0 for j in range(1, len(rot_rad_suggestions)): if rot_rad_suggestions[j] - rot_rad_suggestions[j-1] > tol: last_j = j elif j - last_j > best_rot_rad_sz: rot_rad = rot_rad_suggestions[last_j:j+1].mean() best_rot_rad_sz = j - last_j dx = int(round(rot_rad * 1024 / 360)) return dx, rot_rad
Example #24
Source File: post_proc.py From HorizonNet with MIT License | 6 votes |
def vote(vec, tol): vec = np.sort(vec) n = np.arange(len(vec))[::-1] n = n[:, None] - n[None, :] + 1.0 l = squareform(pdist(vec[:, None], 'minkowski', p=1) + 1e-9) invalid = (n < len(vec) * 0.4) | (l > tol) if (~invalid).sum() == 0 or len(vec) < tol: best_fit = np.median(vec) p_score = 0 else: l[invalid] = 1e5 n[invalid] = -1 score = n max_idx = score.argmax() max_row = max_idx // len(vec) max_col = max_idx % len(vec) assert max_col > max_row best_fit = vec[max_row:max_col+1].mean() p_score = (max_col - max_row + 1) / len(vec) l1_score = np.abs(vec - best_fit).mean() return best_fit, p_score, l1_score
Example #25
Source File: sample_onnx.py From iAI with MIT License | 6 votes |
def process_output(output, file_format, ref_file, topK): if file_format == "ascii": gold_ref = read_ascii_file(ref_file, output.size) res_vec = np.argsort(-output)[:topK] ref_vec = np.argsort(-gold_ref)[:topK] for k in range(0, topK): print(str(res_vec[k]) + " " + str(ref_vec[k])) elif file_format == "ppm": ref_vec = read_reference_file(ref_file) desc_output = np.argsort(-output) output_sorted = np.sort(output)[::-1] for k in range(0, topK): print("Class : "+str(desc_output[k] + 1)+" == "+ ref_vec[desc_output[k]]) else: print("Format Not Supported") sys.exit()
Example #26
Source File: 21-matrix_A_B.py From pyscf with Apache License 2.0 | 6 votes |
def diagonalize(a, b, nroots=4): a_aa, a_ab, a_bb = a b_aa, b_ab, b_bb = b nocc_a, nvir_a, nocc_b, nvir_b = a_ab.shape a_aa = a_aa.reshape((nocc_a*nvir_a,nocc_a*nvir_a)) a_ab = a_ab.reshape((nocc_a*nvir_a,nocc_b*nvir_b)) a_bb = a_bb.reshape((nocc_b*nvir_b,nocc_b*nvir_b)) b_aa = b_aa.reshape((nocc_a*nvir_a,nocc_a*nvir_a)) b_ab = b_ab.reshape((nocc_a*nvir_a,nocc_b*nvir_b)) b_bb = b_bb.reshape((nocc_b*nvir_b,nocc_b*nvir_b)) a = numpy.bmat([[ a_aa , a_ab], [ a_ab.T, a_bb]]) b = numpy.bmat([[ b_aa , b_ab], [ b_ab.T, b_bb]]) e = numpy.linalg.eig(numpy.bmat([[a , b ], [-b.conj(),-a.conj()]]))[0] lowest_e = numpy.sort(e[e.real > 0].real)[:nroots] return lowest_e
Example #27
Source File: test_tduks.py From pyscf with Apache License 2.0 | 6 votes |
def diagonalize(a, b, nroots=4): a_aa, a_ab, a_bb = a b_aa, b_ab, b_bb = b nocc_a, nvir_a, nocc_b, nvir_b = a_ab.shape a_aa = a_aa.reshape((nocc_a*nvir_a,nocc_a*nvir_a)) a_ab = a_ab.reshape((nocc_a*nvir_a,nocc_b*nvir_b)) a_bb = a_bb.reshape((nocc_b*nvir_b,nocc_b*nvir_b)) b_aa = b_aa.reshape((nocc_a*nvir_a,nocc_a*nvir_a)) b_ab = b_ab.reshape((nocc_a*nvir_a,nocc_b*nvir_b)) b_bb = b_bb.reshape((nocc_b*nvir_b,nocc_b*nvir_b)) a = numpy.bmat([[ a_aa , a_ab], [ a_ab.T, a_bb]]) b = numpy.bmat([[ b_aa , b_ab], [ b_ab.T, b_bb]]) e = numpy.linalg.eig(numpy.bmat([[a , b ], [-b.conj(),-a.conj()]]))[0] lowest_e = numpy.sort(e[e.real > 0].real)[:nroots] lowest_e = lowest_e[lowest_e > 1e-3] return lowest_e
Example #28
Source File: linalg_helper.py From pyscf with Apache License 2.0 | 6 votes |
def _sort_by_similarity(w, v, nroots, conv, vlast, emin=None, heff=None): if not any(conv) or vlast is None: return w[:nroots], v[:,:nroots] head, nroots = vlast.shape conv = numpy.asarray(conv[:nroots]) ovlp = vlast[:,conv].T.conj().dot(v[:head]) ovlp = numpy.einsum('ij,ij->j', ovlp, ovlp) nconv = numpy.count_nonzero(conv) nleft = nroots - nconv idx = ovlp.argsort() sorted_idx = numpy.zeros(nroots, dtype=int) sorted_idx[conv] = numpy.sort(idx[-nconv:]) sorted_idx[~conv] = numpy.sort(idx[:-nconv])[:nleft] e = w[sorted_idx] c = v[:,sorted_idx] return e, c
Example #29
Source File: linalg_helper.py From pyscf with Apache License 2.0 | 6 votes |
def pick_real_eigs(w, v, nroots, envs): '''This function searchs the real eigenvalues or eigenvalues with small imaginary component. ''' threshold = 1e-3 abs_imag = abs(w.imag) # Grab `nroots` number of e with small(est) imaginary components max_imag_tol = max(threshold, numpy.sort(abs_imag)[min(w.size,nroots)-1]) real_idx = numpy.where((abs_imag <= max_imag_tol))[0] nbelow_thresh = numpy.count_nonzero(abs_imag[real_idx] < threshold) if nbelow_thresh < nroots and w.size >= nroots: warnings.warn('Only %d eigenvalues (out of %3d requested roots) with imaginary part < %4.3g.\n' % (nbelow_thresh, min(w.size,nroots), threshold)) # Guess whether the matrix to diagonalize is real or complex if envs.get('dtype') == numpy.double: w, v, idx = _eigs_cmplx2real(w, v, real_idx, real_eigenvectors=True) else: w, v, idx = _eigs_cmplx2real(w, v, real_idx, real_eigenvectors=False) return w, v, idx
Example #30
Source File: addons.py From pyscf with Apache License 2.0 | 6 votes |
def dynamic_occ_(mf, tol=1e-3): assert(isinstance(mf, hf.RHF)) old_get_occ = mf.get_occ def get_occ(mo_energy, mo_coeff=None): mol = mf.mol nocc = mol.nelectron // 2 sort_mo_energy = numpy.sort(mo_energy) lumo = sort_mo_energy[nocc] if abs(sort_mo_energy[nocc-1] - lumo) < tol: mo_occ = numpy.zeros_like(mo_energy) mo_occ[mo_energy<lumo] = 2 lst = abs(mo_energy - lumo) < tol mo_occ[lst] = 0 logger.warn(mf, 'set charge = %d', mol.charge+int(lst.sum())*2) logger.info(mf, 'HOMO = %.12g LUMO = %.12g', sort_mo_energy[nocc-1], sort_mo_energy[nocc]) logger.debug(mf, ' mo_energy = %s', sort_mo_energy) else: mo_occ = old_get_occ(mo_energy, mo_coeff) return mo_occ mf.get_occ = get_occ return mf