Skip to content

Commit

Permalink
formmating
Browse files Browse the repository at this point in the history
  • Loading branch information
jasonkaye committed Sep 21, 2023
1 parent 2b4440a commit fb25b3d
Show file tree
Hide file tree
Showing 6 changed files with 24 additions and 27 deletions.
2 changes: 1 addition & 1 deletion c++/cppdlr/dlr_build.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -203,7 +203,7 @@ namespace cppdlr {

int nom = om.size();

if (statistic==Fermion){
if (statistic == Fermion) {

// 2*n+1 should go from -2*nmax+1 to 2*nmax-1, so n goes from -nmax to nmax-1
auto kmat = nda::matrix<dcomplex>(2 * nmax, nom);
Expand Down
4 changes: 2 additions & 2 deletions c++/cppdlr/dlr_build.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -79,12 +79,12 @@ namespace cppdlr {

/**
* @copydoc build_k_it(nda::vector_const_view<double>, nda::vector_const_view<double>)
*/
*/
nda::vector<double> build_k_it(double t, nda::vector_const_view<double> om);

/**
* @copydoc build_k_it(nda::vector_const_view<double>, nda::vector_const_view<double>)
*/
*/
nda::vector<double> build_k_it(nda::vector_const_view<double> t, double om);

/**
Expand Down
11 changes: 5 additions & 6 deletions c++/cppdlr/dlr_imfreq.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -47,19 +47,19 @@ namespace cppdlr {

} else { // Symmetrized bosonic case: enforce n = 0 as DLR imaginary frequency node

auto kvec0 = nda::vector<dcomplex>(kmat(nmax,range(r))); // K at n = 0: K(0, om)
auto kvec0 = nda::vector<dcomplex>(kmat(nmax, range(r))); // K at n = 0: K(0, om)

// Remove row of K matrix corresponding to n = 0 and column corresponding
// to omega=0
kmat(range(nmax,2*nmax),range(r)) = kmat(range(nmax+1,2*nmax+1),range(r));
kmat(range(nmax, 2 * nmax), range(r)) = kmat(range(nmax + 1, 2 * nmax + 1), range(r));

// Pivoted Gram-Schmidt on rows of modified K matrix, augmented by vector
// K(0, om) which was removed, to obtain DLR imaginary frequency nodes
auto [q, norms, piv] = pivrgs_sym(kmat(range(2*nmax),range(r)), kvec0, 1e-100);
auto [q, norms, piv] = pivrgs_sym(kmat(range(2 * nmax), range(r)), kvec0, 1e-100);

std::sort(piv.begin(), piv.end()); // Sort pivots in ascending order
for (int i = 0; i < (r - 1) / 2; ++i) { // Fill in negative im freq nodes
dlr_if(i) = piv(i + 1)-1 - nmax; // Shift by 1 to obtain pivots of original matrix, not augmented matrix
dlr_if(i) = piv(i + 1) - 1 - nmax; // Shift by 1 to obtain pivots of original matrix, not augmented matrix
}
dlr_if((r - 1) / 2) = 0; // i*nu_n = 0 is always chosen
for (int i = (r - 1) / 2 + 1; i < r; ++i) { // Fill in time nodes on (1/2,1]
Expand All @@ -72,9 +72,8 @@ namespace cppdlr {
}
for (int j = 0; j < r; ++j) { cf2if((r - 1) / 2, j) = kvec0(j); } // Entries corresponding to i*nu_n = 0
for (int i = (r - 1) / 2 + 1; i < r; ++i) { // Entries corresponding to positive im freq nodes
for (int j = 0; j < r; ++j) { cf2if(i, j) = kmat(piv(i)-1, j); }
for (int j = 0; j < r; ++j) { cf2if(i, j) = kmat(piv(i) - 1, j); }
}

}

// Prepare imaginary time values to coefficients transformation by computing
Expand Down
2 changes: 1 addition & 1 deletion test/c++/imfreq_ops.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -234,7 +234,7 @@ TEST(imfreq_ops, interp_matrix_sym) {
*/
TEST(imfreq_ops, interp_matrix_sym_bos) {

double lambda = 1000; // DLR cutoff
double lambda = 1000; // DLR cutoff
double eps = 1e-10; // DLR tolerance
auto statistic = Boson; // Bosonic Green's function

Expand Down
28 changes: 13 additions & 15 deletions test/c++/symcompare_if.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -45,15 +45,14 @@ nda::matrix<dcomplex> gfun(int norb, double beta, int n, statistic_t statistic)
return g;
}


/**
* @brief Compare symmetrized and unsymmetrized DLR for interpolation and
* evaluation of matrix-valued Green's function in imaginary time
*/
int main() {

double lambda = 1000; // DLR cutoff
double eps = 1e-10; // DLR tolerance
double lambda = 1000; // DLR cutoff
double eps = 1e-10; // DLR tolerance
auto statistic = Boson; // Fermionic Green's statistics

double beta = 1000; // Inverse temperature
Expand All @@ -69,15 +68,14 @@ int main() {
int rsym = dlr_rf_sym.size();

// Get DLR imaginary frequency object
auto ifops = imfreq_ops(lambda, dlr_rf, statistic);
auto ifops = imfreq_ops(lambda, dlr_rf, statistic);
auto ifops_sym = imfreq_ops(lambda, dlr_rf_sym, statistic, SYM);

// Sample Green's function G at DLR imaginary frequency nodes
auto const &dlr_if = ifops.get_ifnodes();
auto const &dlr_if = ifops.get_ifnodes();
auto const &dlr_if_sym = ifops_sym.get_ifnodes();


auto g = nda::array<dcomplex, 3>(r, norb, norb);
auto g = nda::array<dcomplex, 3>(r, norb, norb);
auto g_sym = nda::array<dcomplex, 3>(rsym, norb, norb);
for (int i = 0; i < r; ++i) { g(i, _, _) = gfun(norb, beta, dlr_if(i), statistic); }
for (int i = 0; i < rsym; ++i) { g_sym(i, _, _) = gfun(norb, beta, dlr_if_sym(i), statistic); }
Expand All @@ -87,18 +85,18 @@ int main() {
auto gc_sym = ifops_sym.vals2coefs(beta, g_sym);

// Compute L infinity error
auto gtru = nda::matrix<dcomplex>(norb, norb);
auto gtst = nda::matrix<dcomplex>(norb, norb);
auto gtst_sym = nda::matrix<dcomplex>(norb, norb);
auto gtru = nda::matrix<dcomplex>(norb, norb);
auto gtst = nda::matrix<dcomplex>(norb, norb);
auto gtst_sym = nda::matrix<dcomplex>(norb, norb);
double err = 0, err_sym = 0;
for (int n = -nmaxtst; n < nmaxtst; ++n) {
gtru = gfun(norb, beta, n, statistic);
gtst = ifops.coefs2eval(beta, gc, n);

gtst = ifops.coefs2eval(beta, gc, n);
gtst_sym = ifops_sym.coefs2eval(beta, gc_sym, n);
err = std::max(err, max_element(abs(gtru - gtst)));
err_sym = std::max(err_sym, max_element(abs(gtru - gtst_sym)));

err = std::max(err, max_element(abs(gtru - gtst)));
err_sym = std::max(err_sym, max_element(abs(gtru - gtst_sym)));
}

// Print results
Expand Down
4 changes: 2 additions & 2 deletions test/c++/symcompare_it.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -49,8 +49,8 @@ nda::matrix<double> gfun(int norb, double beta, double t) {
*/
int main() {

double lambda = 1000; // DLR cutoff
double eps = 1e-10; // DLR tolerance
double lambda = 1000; // DLR cutoff
double eps = 1e-10; // DLR tolerance
auto statistic = Boson; // Green's function statistics

double beta = 1000; // Inverse temperature
Expand Down

0 comments on commit fb25b3d

Please sign in to comment.