Commit ea6e13c7 authored by Carl Poelking's avatar Carl Poelking
Browse files

Graph dMap support, DEBUG option.

parent 560f13c7
......@@ -202,8 +202,10 @@ void AtomicSpectrumFT::contractDeep() {
int lm = l*l+l+m;
phi_l_s1s2 += f1(0, lm)*std::conj(f2(0, lm));
}
GLOG() << l << ":" << l << std::endl;
GLOG() << "Save @ " << offset << std::endl;
#ifdef DBG
GLOG() << l << ":" << l << std::endl;
GLOG() << "Save @ " << offset << std::endl;
#endif
coeffs(offset, 0) = pow(this->getCenter()->getSigma(), 2*l)*inv_alpha_l*phi_l_s1s2.real();
offset += 1;
} // l
......@@ -224,8 +226,10 @@ void AtomicSpectrumFT::contractDeep() {
phi_l1_l1l1 += w111 * f1(0,lm1) * f2(0,lm2) * f2(0,lm3);
}}}
// Store l1l1l2_s1s1s2
GLOG() << l1 << "::" << l1 << ":" << l1 << " " << w111_sq << std::endl;
GLOG() << "Save @ " << offset << std::endl;
#ifdef DBG
GLOG() << l1 << "::" << l1 << ":" << l1 << " " << w111_sq << std::endl;
GLOG() << "Save @ " << offset << std::endl;
#endif
coeffs(offset, 0) =
sqrt(inv_alpha_l1)*sqrt(inv_alpha_l1)*sqrt(inv_alpha_l1)
* pow(this->getCenter()->getSigma(), l1+l1+l1+0.5) *
......@@ -260,8 +264,10 @@ void AtomicSpectrumFT::contractDeep() {
phi_l1_l1l2 += w112 * f1(0,lm1) * f2(0,lm2) * f2(0,lm3);
}}}
// Store l1l1l2_s1s1s2
GLOG() << l1 << ":" << l1 << "::" << l2 << " " << w112_sq << std::endl;
GLOG() << "Save @ " << offset << std::endl;
#ifdef DBG
GLOG() << l1 << ":" << l1 << "::" << l2 << " " << w112_sq << std::endl;
GLOG() << "Save @ " << offset << std::endl;
#endif
coeffs(offset, 0) =
sqrt(inv_alpha_l1)*sqrt(inv_alpha_l1)*sqrt(inv_alpha_l2)
* pow(this->getCenter()->getSigma(), l1+l1+l2+0.5) *
......@@ -269,8 +275,10 @@ void AtomicSpectrumFT::contractDeep() {
* phi_l1l1_l2.real();
offset += 1;
// Store l1l1l2_s1s2s2
GLOG() << l1 << "::" << l1 << ":" << l2 << " " << w112_sq << std::endl;
GLOG() << "Save @ " << offset << std::endl;
#ifdef DBG
GLOG() << l1 << "::" << l1 << ":" << l2 << " " << w112_sq << std::endl;
GLOG() << "Save @ " << offset << std::endl;
#endif
coeffs(offset, 0) =
sqrt(inv_alpha_l1)*sqrt(inv_alpha_l1)*sqrt(inv_alpha_l2)
* pow(this->getCenter()->getSigma(), l1+l1+l2+0.5) *
......@@ -303,8 +311,10 @@ void AtomicSpectrumFT::contractDeep() {
phi_l3_l1l2 += w123 * f1(0,lm3) * f2(0,lm1) * f2(0,lm2);
}}}
// Store (l1)(l2l3)
GLOG() << l1 << "::" << l2 << ":" << l3 << " " << w123_sq << std::endl;
GLOG() << "Save @ " << offset << std::endl;
#ifdef DBG
GLOG() << l1 << "::" << l2 << ":" << l3 << " " << w123_sq << std::endl;
GLOG() << "Save @ " << offset << std::endl;
#endif
coeffs(offset, 0) =
sqrt(inv_alpha_l1)*sqrt(inv_alpha_l2)*sqrt(inv_alpha_l3)
* pow(this->getCenter()->getSigma(), l1+l2+l3+0.5) *
......@@ -312,8 +322,10 @@ void AtomicSpectrumFT::contractDeep() {
* phi_l1_l2l3.real();
offset += 1;
// Store (l2)(l3l1)
GLOG() << l2 << "::" << l3 << ":" << l1 << " " << w123_sq << std::endl;
GLOG() << "Save @ " << offset << std::endl;
#ifdef DBG
GLOG() << l2 << "::" << l3 << ":" << l1 << " " << w123_sq << std::endl;
GLOG() << "Save @ " << offset << std::endl;
#endif
coeffs(offset, 0) =
sqrt(inv_alpha_l1)*sqrt(inv_alpha_l2)*sqrt(inv_alpha_l3)
* pow(this->getCenter()->getSigma(), l1+l2+l3+0.5) *
......@@ -321,8 +333,10 @@ void AtomicSpectrumFT::contractDeep() {
* phi_l2_l3l1.real();
offset += 1;
// Store (l3)(l1l2)
GLOG() << l3 << "::" << l1 << ":" << l2 << " " << w123_sq << std::endl;
GLOG() << "Save @ " << offset << std::endl;
#ifdef DBG
GLOG() << l3 << "::" << l1 << ":" << l2 << " " << w123_sq << std::endl;
GLOG() << "Save @ " << offset << std::endl;
#endif
coeffs(offset, 0) =
sqrt(inv_alpha_l1)*sqrt(inv_alpha_l2)*sqrt(inv_alpha_l3)
* pow(this->getCenter()->getSigma(), l1+l2+l3+0.5) *
......
......@@ -19,38 +19,109 @@ class Graph(object):
idx=-1,
label='',
feature_mat=None,
feature_mat_avg=None,
position_mat=None,
connectivity_mat=None,
vertex_info=[],
graph_info={}):
# Labels
self.idx = idx
self.label = label
self.graph_info = graph_info
# Vertex data: descriptors
self.P = feature_mat
self.P_avg = feature_mat_avg
self.P_type_str = str(type(self.P))
# Vertex data: positions, labels
self.R = position_mat
self.C = connectivity_mat
self.vertex_info = vertex_info
self.graph_info = graph_info
# Edge data
self.C = connectivity_mat
return
def save_to_h5(self, h5f, dtype='float32'):
group = h5f.create_group('%06d' % self.idx)
group.create_dataset('feature_mat', data=self.P, compression='gzip', dtype=dtype)
if self.P_type_str == "<class 'soap.soapy.kernel.DescriptorMapMatrix'>":
# Save list of descriptor maps
g0 = group.create_group('feature_dmap')
for dmap_idx, dmap in enumerate(self.P):
g1 = g0.create_group('%d' % dmap_idx)
for key in dmap:
g1.create_dataset(
key,
data=dmap[key],
compression='gzip',
dtype=dtype)
# Save averaged descriptor map
g0_avg = group.create_group('feature_dmap_avg')
for key in self.P_avg:
g0_avg.create_dataset(
key,
data=self.P_avg[key],
compression='gzip',
dtype=dtype)
elif self.P_type_str == "<type 'numpy.ndarray'>":
# Save numpy arrays
group.create_dataset(
'feature_mat',
data=self.P,
compression='gzip',
dtype=dtype)
group.create_dataset(
'feature_mat_avg',
data=self.P_avg,
compression='gzip',
dtype=dtype)
else: raise NotImplementedError(self.P_type_str)
group.create_dataset('position_mat', data=self.R)
group.create_dataset('connectivity_mat', data=self.C)
group.attrs['idx'] = self.idx
group.attrs['label'] = self.label
group.attrs['vertex_info'] = json.dumps(self.vertex_info)
group.attrs['graph_info'] = json.dumps(self.graph_info)
group.attrs['P_type_str'] = self.P_type_str
return
def load_from_h5(self, h5f):
self.idx = h5f.attrs['idx']
self.label = h5f.attrs['label']
self.vertex_info = json.loads(h5f.attrs['vertex_info'])
self.graph_info = json.loads(h5f.attrs['graph_info'])
self.P = h5f['feature_mat'].value
# Determine feature matrix type
if 'P_type_str'in h5f.attrs:
self.P_type_str = h5f.attrs['P_type_str']
else:
self.P_type_str = "<type 'numpy.ndarray'>"
if self.P_type_str == "<class 'soap.soapy.kernel.DescriptorMapMatrix'>":
# Load list of descriptor maps
self.P = soap.soapy.kernel.DescriptorMapMatrix()
g0 = h5f['feature_dmap']
for i in range(len(g0)):
Pi = soap.soapy.kernel.DescriptorMap()
g1 = g0['%d' % i]
for key in g1:
Pi[key] = g1[key].value
self.P.append(Pi)
# Load averaged descriptor map
self.P_avg = soap.soapy.kernel.DescriptorMap()
g0_avg = h5f['feature_dmap_avg']
for key in g0_avg:
self.P_avg[key] = g0_avg[key].value
elif self.P_type_str == "<type 'numpy.ndarray'>":
self.P = h5f['feature_mat'].value
self.P_avg = h5f['feature_mat_avg'].value
else: raise NotImplementedError(self.P_type_str)
self.R = h5f['position_mat'].value
self.C = h5f['connectivity_mat'].value
return self
def load_graphs(hdf5, n=-1):
h = h5py.File(hdf5, 'r')
gsec = h['graphs']
if n == None or n < 0: n = len(gsec)
graphs = [ Graph().load_from_h5(
gsec['%06d' % i]) for i in range(len(gsec)) if i < n ]
h.close()
return graphs
def calculate_graph_connectivity(graph, zero_diagonal):
"""
Suited for periodic-table elements:
......@@ -297,7 +368,7 @@ def mp_compute_graph(
vertex_info = frag_labels
else:
connectivity_mat = atom_bond_mat
vertex_info = atom_labels
vertex_info = type_vec
graph = Graph(
idx = config.info['idx'],
label = str(config.info['label']),
......@@ -332,7 +403,7 @@ def mp_compute_kernel_block_hdf5(block, kernel, log, dtype_result, h5_file):
#soap.soapy.util.MP_LOCK.release()
return mp_compute_kernel_block([g_rows, g_cols], kernel, log, dtype_result)
def mp_compute_kernel_block(block, kernel, log, dtype_result):
def mp_compute_kernel_block(block, kernel, log, dtype_result, symmetric=True):
gi_block = block[0]
gj_block = block[1]
if log:
......@@ -342,7 +413,7 @@ def mp_compute_kernel_block(block, kernel, log, dtype_result):
kmat = np.zeros((len(gi_block),len(gj_block)), dtype=dtype_result)
for i,gi in enumerate(gi_block):
for j,gj in enumerate(gj_block):
if gi.idx > gj.idx: pass
if symmetric and gi.idx > gj.idx: pass
else: kmat[i,j] = kernel.compute(gi,gj)
return kmat
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment