Skip to content

Commit d2c57bb

Browse files
committed
Combine DiagramEvaluator classes
1 parent 6c69690 commit d2c57bb

12 files changed

Lines changed: 196 additions & 268 deletions

c++/triqs_xca/block_sparse_backbone.cpp

Lines changed: 41 additions & 50 deletions
Original file line numberDiff line numberDiff line change
@@ -11,9 +11,9 @@
1111
#include <triqs_xca/block_sparse.hpp>
1212
#include <triqs_xca/block_sparse_backbone.hpp>
1313

14-
DiagramBlockSparseEvaluator::DiagramBlockSparseEvaluator(double beta, imtime_ops &itops, nda::array_const_view<dcomplex, 3> hyb,
15-
nda::array_const_view<dcomplex, 3> hyb_refl, nda::vector_const_view<double> hyb_poles,
16-
BlockDiagOpFun &Gt, BlockOpSymQuartet &Fq)
14+
DiagramEvaluator::DiagramEvaluator(double beta, imtime_ops &itops, nda::array_const_view<dcomplex, 3> hyb,
15+
nda::array_const_view<dcomplex, 3> hyb_refl, nda::vector_const_view<double> hyb_poles, BlockDiagOpFun &Gt,
16+
BlockOpSymQuartet &Fq)
1717
: beta(beta),
1818
r(itops.rank()),
1919
n(hyb.extent(1)),
@@ -31,13 +31,25 @@ DiagramBlockSparseEvaluator::DiagramBlockSparseEvaluator(double beta, imtime_ops
3131

3232
// allocate arrays
3333
T = nda::zeros<dcomplex>(r, Nmax, Nmax);
34+
U = nda::zeros<dcomplex>(r, Nmax, Nmax);
3435
GKt = nda::zeros<dcomplex>(r, Nmax, Nmax);
3536
Tkaps = nda::zeros<dcomplex>(n, r, Nmax, Nmax);
3637
Tmu = nda::zeros<dcomplex>(r, Nmax, Nmax);
3738
}
3839

39-
void DiagramBlockSparseEvaluator::multiply_vertex_block(Backbone &backbone, int v_ix, nda::vector_const_view<int> ind_path,
40-
nda::vector_const_view<int> block_dims) {
40+
void DiagramEvaluator::reset() {
41+
T = 0;
42+
U = 0;
43+
GKt = 0;
44+
Tkaps = 0;
45+
Tmu = 0;
46+
for (int i = 0; i < Sigma.get_num_block_cols(); i++) {
47+
Sigma.set_block(i, nda::zeros<dcomplex>(r, Sigma.get_block_size(i), Sigma.get_block_size(i)));
48+
}
49+
}
50+
51+
void DiagramEvaluator::multiply_vertex_block(Backbone &backbone, int v_ix, nda::vector_const_view<int> ind_path,
52+
nda::vector_const_view<int> block_dims) {
4153
int o_ix = backbone.get_vertex_orb(v_ix); // orbital_index
4254
// split backbone orbital index into symmetry set index and orbital index within the symmetry set
4355
// i.e. have mapping between backbone orbital index and symmetry set index
@@ -83,8 +95,8 @@ void DiagramBlockSparseEvaluator::multiply_vertex_block(Backbone &backbone, int
8395
}
8496
}
8597

86-
void DiagramBlockSparseEvaluator::compose_with_edge_block(Backbone &backbone, int e_ix, nda::vector_const_view<int> ind_path,
87-
nda::vector_const_view<int> block_dims) {
98+
void DiagramEvaluator::compose_with_edge_block(Backbone &backbone, int e_ix, nda::vector_const_view<int> ind_path,
99+
nda::vector_const_view<int> block_dims) {
88100

89101
int n_col_r = e_ix < backbone.get_topology(0, 1) ? block_dims(1) : block_dims(0);
90102
int b_ix = ind_path(e_ix); // block index for the edge e_ix
@@ -104,8 +116,8 @@ void DiagramBlockSparseEvaluator::compose_with_edge_block(Backbone &backbone, in
104116
itops.vals2coefs(T(_, range(0, block_dims(e_ix + 1)), range(0, n_col_r))), TIME_ORDERED);
105117
}
106118

107-
void DiagramBlockSparseEvaluator::multiply_zero_vertex_block(Backbone &backbone, bool is_forward, int b_ix_0, int p_kap, int p_mu,
108-
nda::vector_const_view<int> ind_path, nda::vector_const_view<int> block_dims) {
119+
void DiagramEvaluator::multiply_zero_vertex_block(Backbone &backbone, bool is_forward, int b_ix_0, int p_kap, int p_mu,
120+
nda::vector_const_view<int> ind_path, nda::vector_const_view<int> block_dims) {
109121

110122
int b_ix_mu = ind_path(backbone.get_topology(0, 1) - 1); // block index for F_mu
111123
if (is_forward) {
@@ -155,7 +167,7 @@ void DiagramBlockSparseEvaluator::multiply_zero_vertex_block(Backbone &backbone,
155167
}
156168
}
157169

158-
void DiagramBlockSparseEvaluator::multiply_prefactor(Backbone &backbone, nda::vector_const_view<int> block_dims) {
170+
void DiagramEvaluator::multiply_prefactor(Backbone &backbone, nda::vector_const_view<int> block_dims) {
159171
int m = backbone.m;
160172
// Multiply by prefactor
161173
for (int m_ix = 0; m_ix < m - 1; m_ix++) {
@@ -169,20 +181,10 @@ void DiagramBlockSparseEvaluator::multiply_prefactor(Backbone &backbone, nda::ve
169181
}
170182
}
171183

172-
void DiagramBlockSparseEvaluator::reset() {
173-
T = 0;
174-
GKt = 0;
175-
Tkaps = 0;
176-
Tmu = 0;
177-
for (int i = 0; i < Sigma.get_num_block_cols(); i++) {
178-
Sigma.set_block(i, nda::zeros<dcomplex>(r, Sigma.get_block_size(i), Sigma.get_block_size(i)));
179-
}
180-
}
181-
182-
void DiagramBlockSparseEvaluator::eval_diagram_block_sparse(Backbone &backbone) {
184+
void DiagramEvaluator::eval_self_energy(Backbone &backbone) {
183185

184186
int m = backbone.m;
185-
nda::vector<int> ind_path(2 * m - 1); // tracks block indices of factors for computing a particular block of the self-energy
187+
nda::vector<int> ind_path(2 * m - 1); // tracks block indices of factors for computing a particular block of the self-energy
186188
nda::vector<int> block_dims(2 * m + 1); // tracks the dimensions of the blocks in these factors
187189

188190
// loop over all flat indices
@@ -268,7 +270,7 @@ void DiagramBlockSparseEvaluator::eval_diagram_block_sparse(Backbone &backbone)
268270
if (fork_all_nonzero) {
269271
// evaluate the diagram with these directions, poles, and orbital indices
270272
// b_ix is the block index for the first edge
271-
eval_backbone_fixed_indices_block_sparse(backbone, b_ix, p_kap, p_mu, ind_path, block_dims);
273+
eval_self_energy_fixed_indices(backbone, b_ix, p_kap, p_mu, ind_path, block_dims);
272274
}
273275
}
274276
}
@@ -279,9 +281,8 @@ void DiagramBlockSparseEvaluator::eval_diagram_block_sparse(Backbone &backbone)
279281
Sigma.set_zero_block_indices(); // set zero_block_indices according to current blocks
280282
}
281283

282-
void DiagramBlockSparseEvaluator::eval_backbone_fixed_indices_block_sparse(Backbone &backbone, int b_ix, int p_kap, int p_mu,
283-
nda::vector_const_view<int> ind_path,
284-
nda::vector_const_view<int> block_dims) {
284+
void DiagramEvaluator::eval_self_energy_fixed_indices(Backbone &backbone, int b_ix, int p_kap, int p_mu, nda::vector_const_view<int> ind_path,
285+
nda::vector_const_view<int> block_dims) {
285286

286287
int m = backbone.m;
287288

@@ -307,16 +308,8 @@ void DiagramBlockSparseEvaluator::eval_backbone_fixed_indices_block_sparse(Backb
307308
Sigma.add_block(b_ix, T(_, range(0, block_dims(2 * m)), range(0, block_dims(0))));
308309
}
309310

310-
CorrelatorDiagramBlockSparseEvaluator::CorrelatorDiagramBlockSparseEvaluator(double beta, imtime_ops &itops, nda::array_const_view<dcomplex, 3> hyb,
311-
nda::array_const_view<dcomplex, 3> hyb_refl,
312-
nda::vector_const_view<double> hyb_poles, BlockDiagOpFun &Gt,
313-
BlockOpSymQuartet &Fq)
314-
: DiagramBlockSparseEvaluator(beta, itops, hyb, hyb_refl, hyb_poles, Gt, Fq) {
315-
U = nda::zeros<dcomplex>(r, Nmax, Nmax);
316-
}
317-
318-
void CorrelatorDiagramBlockSparseEvaluator::multiply_vertex_corr_left_block(Backbone &backbone, int v_ix, nda::vector_const_view<int> ind_path,
319-
nda::vector_const_view<int> block_dims) {
311+
void DiagramEvaluator::multiply_vertex_corr_block(Backbone &backbone, int v_ix, nda::vector_const_view<int> ind_path,
312+
nda::vector_const_view<int> block_dims) {
320313

321314
int o_ix = backbone.get_vertex_orb(v_ix); // orbital_index
322315
// split backbone orbital index into symmetry set index and orbital index within the symmetry set
@@ -363,8 +356,8 @@ void CorrelatorDiagramBlockSparseEvaluator::multiply_vertex_corr_left_block(Back
363356
}
364357
}
365358

366-
void CorrelatorDiagramBlockSparseEvaluator::compose_with_edge_corr_left_block(Backbone &backbone, int e_ix, nda::vector_const_view<int> ind_path,
367-
nda::vector_const_view<int> block_dims) {
359+
void DiagramEvaluator::compose_with_edge_corr_block(Backbone &backbone, int e_ix, nda::vector_const_view<int> ind_path,
360+
nda::vector_const_view<int> block_dims) {
368361

369362
int n_col_r = e_ix < backbone.get_topology(0, 1) ? block_dims(1) : block_dims(0);
370363
int b_ix = ind_path(e_ix); // block index for the edge e_ix
@@ -384,14 +377,13 @@ void CorrelatorDiagramBlockSparseEvaluator::compose_with_edge_corr_left_block(Ba
384377
itops.vals2coefs(U(_, range(0, block_dims(e_ix + 1)), range(0, n_col_r))), TIME_ORDERED);
385378
}
386379

387-
nda::array<dcomplex, 3> CorrelatorDiagramBlockSparseEvaluator::eval_diagram_block_sparse(Backbone &backbone, std::vector<BlockOp> mu_ops,
388-
std::vector<BlockOp> kap_ops) {
380+
nda::array<dcomplex, 3> DiagramEvaluator::eval_correlator(Backbone &backbone, std::vector<BlockOp> mu_ops, std::vector<BlockOp> kap_ops) {
389381
int m = backbone.m;
390-
nda::vector<int> ind_path(2 * m); // tracks block indices of factors for computing a particular block's contribution to the correlator
382+
nda::vector<int> ind_path(2 * m); // tracks block indices of factors for computing a particular block's contribution to the correlator
391383
nda::vector<int> block_dims(2 * m + 1); // tracks the dimensions of the blocks in these factors
392-
int f_ix_max = static_cast<int>(backbone.fb_ix_max * backbone.o_ix_max * pow(hyb_poles.size(), m - 1));
384+
int f_ix_max = static_cast<int>(backbone.fb_ix_max * backbone.o_ix_max * pow(hyb_poles.size(), m - 1));
393385
nda::array<dcomplex, 3> correlator = nda::zeros<dcomplex>(r, mu_ops.size(), kap_ops.size());
394-
nda::array<dcomplex, 3> Tmuop = nda::zeros<dcomplex>(r, Nmax, Nmax);
386+
nda::array<dcomplex, 3> Tmuop = nda::zeros<dcomplex>(r, Nmax, Nmax);
395387
// loop over all flat indices
396388
for (int f_ix = 0; f_ix < f_ix_max; ++f_ix) {
397389
backbone.set_flat_index(f_ix, hyb_poles); // set directions, pole indices, and orbital indices from a single integer index
@@ -403,13 +395,12 @@ nda::array<dcomplex, 3> CorrelatorDiagramBlockSparseEvaluator::eval_diagram_bloc
403395
-------- -------- -------- -------- -------- -------- -------- --------
404396
*/
405397
}
406-
398+
// TODO: implement the full correlator evaluation by looping over all block indices and accumulating contributions
407399
return correlator;
408400
}
409401

410-
void CorrelatorDiagramBlockSparseEvaluator::eval_backbone_fixed_indices_block_sparse(Backbone &backbone, int b_ix, int p_kap, int p_mu,
411-
nda::vector_const_view<int> ind_path,
412-
nda::vector_const_view<int> block_dims) {
402+
void DiagramEvaluator::eval_correlator_fixed_indices(Backbone &backbone, int b_ix, int p_kap, int p_mu, nda::vector_const_view<int> ind_path,
403+
nda::vector_const_view<int> block_dims) {
413404

414405
int m = backbone.m;
415406

@@ -437,11 +428,11 @@ void CorrelatorDiagramBlockSparseEvaluator::eval_backbone_fixed_indices_block_sp
437428
}
438429
}
439430
for (int v = backbone.get_topology(0, 1) + 1; v < 2 * m - 1; ++v) {
440-
multiply_vertex_corr_left_block(backbone, v, ind_path, block_dims);
441-
compose_with_edge_corr_left_block(backbone, v, ind_path, block_dims);
431+
multiply_vertex_corr_block(backbone, v, ind_path, block_dims);
432+
compose_with_edge_corr_block(backbone, v, ind_path, block_dims);
442433
}
443434
// multiply by the last vertex
444-
multiply_vertex_corr_left_block(backbone, 2 * m - 1, ind_path, block_dims);
435+
multiply_vertex_corr_block(backbone, 2 * m - 1, ind_path, block_dims);
445436
// convolve with last edge
446437
GKt(_, range(0, block_dims(2 * m)), range(0, block_dims(2 * m))) = Gt.get_block(ind_path(2 * m - 1));
447438
int bv = backbone.get_vertex_Ksign(2 * m - 1); // sign on K

0 commit comments

Comments
 (0)