Unverified Commit 8473fa9c authored by Peter Doak's avatar Peter Doak Committed by GitHub
Browse files

Merge pull request #176 from gbalduzz/new_models 

New models
parents 0afc61c0 1316ebee

cmake/dca_config.cmake

+28 −12
Original line number Diff line number Diff line 
################################################################################

# Author: Urs R. Haehner (haehneru@itp.phys.ethz.ch)

#         Giovanni Badlduzzi (gbalduzz@itp.phys.ethz.ch)

#

# Build options for DCA++.
@@ -53,8 +54,8 @@ endif() 
# TODO: Add more point groups and lattices.



# Point group

set(DCA_POINT_GROUP "D4" CACHE STRING "Point group symmetry, options are: C6 | D4.")

set_property(CACHE DCA_POINT_GROUP PROPERTY STRINGS C6 D4)

set(DCA_POINT_GROUP "D4" CACHE STRING "Point group symmetry, options are: C6 | D4 | no_symmetry<2>.")

set_property(CACHE DCA_POINT_GROUP PROPERTY STRINGS C6 D4 no_symmetry<2>)



if (DCA_POINT_GROUP STREQUAL "C6")

  set(DCA_POINT_GROUP_INCLUDE
@@ -63,14 +64,18 @@ if (DCA_POINT_GROUP STREQUAL "C6") 
elseif (DCA_POINT_GROUP STREQUAL "D4")

  set(DCA_POINT_GROUP_INCLUDE "dca/phys/domains/cluster/symmetries/point_groups/2d/2d_square.hpp")



elseif (DCA_POINT_GROUP STREQUAL "no_symmetry<2>")

  set(DCA_POINT_GROUP_INCLUDE "dca/phys/domains/cluster/symmetries/point_groups/no_symmetry.hpp")



else()

  message(FATAL_ERROR "Please set DCA_POINT_GROUP to a valid option: C6 | D4.")

endif()



# Lattice type

set(DCA_LATTICE "square" CACHE STRING

    "Lattice type, options are: bilayer | square | triangular | twoband_chain | singleband_chain.")

set_property(CACHE DCA_LATTICE PROPERTY STRINGS bilayer square triangular twoband_chain singleband_chain)

set(DCA_LATTICE "square" CACHE STRING "Lattice type, options are: bilayer | square | triangular |

    hund | twoband_Cu | threeband | FeAs.")

set_property(CACHE DCA_LATTICE PROPERTY STRINGS bilayer square triangular hund twoband_Cu threeband

             FeAs)



if (DCA_LATTICE STREQUAL "bilayer")

  set(DCA_LATTICE_TYPE dca::phys::models::bilayer_lattice<PointGroup>)
@@ -86,19 +91,30 @@ elseif (DCA_LATTICE STREQUAL "triangular") 
  set(DCA_LATTICE_TYPE dca::phys::models::triangular_lattice<PointGroup>)

  set(DCA_LATTICE_INCLUDE

    "dca/phys/models/analytic_hamiltonians/triangular_lattice.hpp")

elseif (DCA_LATTICE STREQUAL "hund")

  set(DCA_LATTICE_TYPE dca::phys::models::HundLattice<PointGroup>)



elseif (DCA_LATTICE STREQUAL "threeband")

  set(DCA_LATTICE_TYPE dca::phys::models::ThreebandHubbard<PointGroup>)

  set(DCA_LATTICE_INCLUDE

    "dca/phys/models/analytic_hamiltonians/threeband_hubbard.hpp")



elseif (DCA_LATTICE STREQUAL "twoband_chain")

  set(DCA_LATTICE_TYPE dca::phys::models::twoband_chain<dca::phys::domains::no_symmetry<1>>)

  set(DCA_LATTICE_INCLUDE

      "dca/phys/models/analytic_hamiltonians/twoband_chain.hpp")



elseif (DCA_LATTICE STREQUAL "singleband_chain")

  set(DCA_LATTICE_TYPE dca::phys::models::singleband_chain<dca::phys::domains::no_symmetry<1>>)

      "dca/phys/models/analytic_hamiltonians/hund_lattice.hpp")

elseif (DCA_LATTICE STREQUAL "FeAs")

  set(DCA_LATTICE_TYPE dca::phys::models::FeAsLattice<PointGroup>)

  set(DCA_LATTICE_INCLUDE

      "dca/phys/models/analytic_hamiltonians/singleband_chain.hpp")

      "dca/phys/models/analytic_hamiltonians/fe_as_lattice.hpp")

elseif (DCA_LATTICE STREQUAL "twoband_Cu")

  set(DCA_LATTICE_TYPE dca::phys::models::TwoBandCu<PointGroup>)

  set(DCA_LATTICE_INCLUDE

      "dca/phys/models/analytic_hamiltonians/twoband_Cu.hpp")



else()

  message(FATAL_ERROR

          "Please set DCA_LATTICE to a valid option: bilayer | square | triangular | twoband_chain | singleband_chain.")

  message(FATAL_ERROR "Please set DCA_LATTICE to a valid option: bilayer | square | triangular |

          hund | twoband_Cu | threeband | FeAs.")

endif()



# Model type

include/dca/function/domains/dmn_0.hpp

+3 −3
Original line number Diff line number Diff line
@@ -37,7 +37,7 @@ public: 
    return parameters::get_size();

  }



  static std::vector<element_type>& get_elements() {

  static const std::vector<element_type>& get_elements() {

    return parameters::get_elements();

  }
@@ -75,7 +75,7 @@ void dmn_0<parameters>::reset() { 
  dmn_0::initialize();

}



}  // func

}  // dca

}  // namespace func

}  // namespace dca



#endif  // DCA_FUNCTION_DOMAINS_DMN_0_HPP

include/dca/phys/dca_algorithms/compute_band_structure.hpp

+1 −1
Original line number Diff line number Diff line
@@ -83,7 +83,7 @@ void compute_band_structure::execute(const ParametersType& parameters, 
  // Compute H(k).

  func::function<std::complex<double>, func::dmn_variadic<nu, nu, k_domain_cut_dmn_type>> H_k(

      "H_k");

  ParametersType::lattice_type::initialize_H_0(parameters, H_k);

  ParametersType::lattice_type::initializeH0(parameters, H_k);



  // Compute the bands.

  dca::linalg::Vector<double, dca::linalg::CPU> L_vec(nu::dmn_size());

include/dca/phys/dca_analysis/bse_solver/bse_cluster_solver.hpp

+10 −8
Original line number Diff line number Diff line
@@ -38,14 +38,16 @@ public: 
  using profiler_t = typename ParametersType::profiler_type;

  using concurrency_t = typename ParametersType::concurrency_type;



  using Lattice = typename ParametersType::lattice_type;



  using w = func::dmn_0<domains::frequency_domain>;

  using w_VERTEX = func::dmn_0<domains::vertex_frequency_domain<domains::COMPACT>>;

  using WVertexDmn = func::dmn_0<domains::vertex_frequency_domain<domains::COMPACT>>;

  using b = func::dmn_0<domains::electron_band_domain>;

  using k_DCA =

      func::dmn_0<domains::cluster_domain<double, ParametersType::lattice_type::DIMENSION, domains::CLUSTER,

                                          domains::MOMENTUM_SPACE, domains::BRILLOUIN_ZONE>>;



  using cluster_eigenvector_dmn_t = func::dmn_variadic<b, b, k_DCA, w_VERTEX>;

  using cluster_eigenvector_dmn_t = func::dmn_variadic<b, b, k_DCA, WVertexDmn>;

  using DCA_matrix_dmn_t = func::dmn_variadic<cluster_eigenvector_dmn_t, cluster_eigenvector_dmn_t>;



  BseClusterSolver(ParametersType& parameters, DcaDataType& data);
@@ -82,7 +84,7 @@ private: 
  diagrammatic_symmetries<ParametersType> diagrammatic_symmetries_obj;



  func::function<std::complex<ScalarType>, DCA_matrix_dmn_t> Gamma_cluster;

  func::function<std::complex<double>, func::dmn_variadic<b, b, b, b, k_DCA, w_VERTEX>> G_II_0_function;

  func::function<std::complex<double>, func::dmn_variadic<b, b, b, b, k_DCA, WVertexDmn>> G_II_0_function;

};



template <typename ParametersType, typename DcaDataType, typename ScalarType>
@@ -129,9 +131,9 @@ void BseClusterSolver<ParametersType, DcaDataType, ScalarType>::apply_symmetries 


  profiler_t prof(__FUNCTION__, __FILE__, __LINE__);



  symmetrize::execute(data_.Sigma, data_.H_symmetry);

  symmetrize::execute<Lattice>(data_.Sigma, data_.H_symmetry);



  symmetrize::execute(data_.G_k_w, data_.H_symmetry);

  symmetrize::execute<Lattice>(data_.G_k_w, data_.H_symmetry);

}



template <typename ParametersType, typename DcaDataType, typename ScalarType>
@@ -209,10 +211,10 @@ void BseClusterSolver<ParametersType, DcaDataType, ScalarType>::load_G_II_0( 


  G_II_0 = 0.;



  func::dmn_variadic<k_DCA, w_VERTEX> k_w_dmn;

  func::dmn_variadic<k_DCA, WVertexDmn> k_w_dmn;



  int W = parameters.get_sp_fermionic_frequencies();  // TODO: Replace using w::dmn_size().

  int W_vertex = w_VERTEX::dmn_size() / 2;

  int W_vertex = WVertexDmn::dmn_size() / 2;

  int q = parameters.get_four_point_momentum_transfer_index();



  int w_nu = parameters.get_four_point_frequency_transfer();
@@ -285,7 +287,7 @@ void BseClusterSolver<ParametersType, DcaDataType, ScalarType>::load_G_II_0_func 
  if (concurrency.id() == concurrency.first())

    std::cout << "\t" << __FUNCTION__ << "\n\n";



  for (int w_ind = 0; w_ind < w_VERTEX::dmn_size(); w_ind++)

  for (int w_ind = 0; w_ind < WVertexDmn::dmn_size(); w_ind++)

    for (int K_ind = 0; K_ind < k_DCA::dmn_size(); K_ind++)



      for (int m2 = 0; m2 < b::dmn_size(); m2++)

include/dca/phys/dca_analysis/bse_solver/bse_lattice_solver.hpp

+17 −17
Original line number Diff line number Diff line
@@ -59,7 +59,7 @@ public: 
  const static int num_harmonics = 3;

  using CubicHarmonicsDmn = func::dmn_0<func::dmn<num_harmonics, int>>;



  using w_VERTEX = func::dmn_0<domains::vertex_frequency_domain<domains::COMPACT>>;

  using WVertexDmn = func::dmn_0<domains::vertex_frequency_domain<domains::COMPACT>>;

  using b = func::dmn_0<domains::electron_band_domain>;

  using b_b = func::dmn_variadic<b, b>;
@@ -81,10 +81,10 @@ public: 


  using chi_vector_dmn_t = func::dmn_variadic<b, b, crystal_harmonics_expansion_dmn_t>;



  using LatticeEigenvectorDmn = func::dmn_variadic<b, b, k_HOST_VERTEX, w_VERTEX>;

  using LatticeEigenvectorDmn = func::dmn_variadic<b, b, k_HOST_VERTEX, WVertexDmn>;

  using crystal_eigenvector_dmn_t =

      func::dmn_variadic<b, b, crystal_harmonics_expansion_dmn_t, w_VERTEX>;

  using CubicHarmonicsEigenvectorDmn = func::dmn_variadic<b, b, CubicHarmonicsDmn, w_VERTEX>;

      func::dmn_variadic<b, b, crystal_harmonics_expansion_dmn_t, WVertexDmn>;

  using CubicHarmonicsEigenvectorDmn = func::dmn_variadic<b, b, CubicHarmonicsDmn, WVertexDmn>;



  using HOST_matrix_dmn_t = func::dmn_variadic<LatticeEigenvectorDmn, LatticeEigenvectorDmn>;
@@ -143,7 +143,7 @@ private: 
  DcaDataType& MOMS;



  func::function<std::complex<ScalarType>, HOST_matrix_dmn_t> Gamma_lattice;

  func::function<std::complex<ScalarType>, func::dmn_variadic<b_b, b_b, k_HOST_VERTEX, w_VERTEX>> chi_0_lattice;

  func::function<std::complex<ScalarType>, func::dmn_variadic<b_b, b_b, k_HOST_VERTEX, WVertexDmn>> chi_0_lattice;

  // Matrix in \vec{k} and \omega_n with the diagonal = chi_0_lattice.

  func::function<std::complex<ScalarType>, HOST_matrix_dmn_t> chi_0_lattice_matrix;
@@ -250,7 +250,7 @@ void BseLatticeSolver<ParametersType, DcaDataType, ScalarType>::initialize() { 
            std::exp(I * math::util::innerProduct(r_vec, k_HOST_VERTEX::get_elements()[k_ind])) /

            std::sqrt(double(k_HOST_VERTEX::dmn_size()));



      for (int w_ind = 0; w_ind < w_VERTEX::dmn_size(); w_ind++)

      for (int w_ind = 0; w_ind < WVertexDmn::dmn_size(); w_ind++)

        for (int k_ind = 0; k_ind < k_HOST_VERTEX::dmn_size(); k_ind++)

          for (int m_ind = 0; m_ind < b::dmn_size(); m_ind++)

            for (int n_ind = 0; n_ind < b::dmn_size(); n_ind++)
@@ -346,7 +346,7 @@ void BseLatticeSolver<ParametersType, DcaDataType, ScalarType>::computeChi0Latti 
  // Renormalize and set diagonal \chi_0 matrix.

  const ScalarType renorm = 1. / (parameters.get_beta() * k_HOST_VERTEX::dmn_size());



  for (int w_ind = 0; w_ind < w_VERTEX::dmn_size(); w_ind++)

  for (int w_ind = 0; w_ind < WVertexDmn::dmn_size(); w_ind++)

    for (int K_ind = 0; K_ind < k_HOST_VERTEX::dmn_size(); K_ind++)

      for (int m2 = 0; m2 < b::dmn_size(); m2++)

        for (int n2 = 0; n2 < b::dmn_size(); n2++)
@@ -739,12 +739,12 @@ void BseLatticeSolver<ParametersType, DcaDataType, ScalarType>::compute_folded_s 


    chi_q = 0.;



    for (int w2 = 0; w2 < w_VERTEX::dmn_size(); w2++)

    for (int w2 = 0; w2 < WVertexDmn::dmn_size(); w2++)

      for (int K2 = 0; K2 < M; K2++)

        for (int m2 = 0; m2 < b::dmn_size(); m2++)

          for (int n2 = 0; n2 < b::dmn_size(); n2++)



            for (int w1 = 0; w1 < w_VERTEX::dmn_size(); w1++)

            for (int w1 = 0; w1 < WVertexDmn::dmn_size(); w1++)

              for (int K1 = 0; K1 < M; K1++)

                for (int m1 = 0; m1 < b::dmn_size(); m1++)

                  for (int n1 = 0; n1 < b::dmn_size(); n1++)
@@ -809,9 +809,9 @@ void BseLatticeSolver<ParametersType, DcaDataType, ScalarType>::printOnShell() { 
        std::complex<ScalarType> norm = 0;



        for (int j = 0; j < k_HOST_VERTEX::dmn_size(); j++) {

          scal_prod += conj(psi_k(j, l)) * leading_eigenvectors(i, 0, 0, j, w_VERTEX::dmn_size() / 2);

          norm += conj(leading_eigenvectors(i, 0, 0, j, w_VERTEX::dmn_size() / 2)) *

                  leading_eigenvectors(i, 0, 0, j, w_VERTEX::dmn_size() / 2);

          scal_prod += conj(psi_k(j, l)) * leading_eigenvectors(i, 0, 0, j, WVertexDmn::dmn_size() / 2);

          norm += conj(leading_eigenvectors(i, 0, 0, j, WVertexDmn::dmn_size() / 2)) *

                  leading_eigenvectors(i, 0, 0, j, WVertexDmn::dmn_size() / 2);

        }

        assert(std::abs(norm.imag()) < 1.e-3);

        const std::complex<ScalarType> result = scal_prod / std::sqrt(norm.real() + 1.e-16);
@@ -838,9 +838,9 @@ void BseLatticeSolver<ParametersType, DcaDataType, ScalarType>::printOnShell() { 


      for (int i = 0; i < num_evals; ++i) {

        std::cout << "\t["

                  << leading_symmetry_decomposition(i, 0, 0, harmonic, w_VERTEX::dmn_size() / 2).real()

                  << leading_symmetry_decomposition(i, 0, 0, harmonic, WVertexDmn::dmn_size() / 2).real()

                  << ", "

                  << leading_symmetry_decomposition(i, 0, 0, harmonic, w_VERTEX::dmn_size() / 2).imag()

                  << leading_symmetry_decomposition(i, 0, 0, harmonic, WVertexDmn::dmn_size() / 2).imag()

                  << "]";

      }

      std::cout << "\n";
@@ -868,13 +868,13 @@ void BseLatticeSolver<ParametersType, DcaDataType, ScalarType>::symmetrizeLeadin 
                                           std::sin((2. * M_PI * l) / num_phases));

      ScalarType diff = 0;



      for (int w = 0; w < w_VERTEX::dmn_size() / 2; w++)

      for (int w = 0; w < WVertexDmn::dmn_size() / 2; w++)

        for (int k = 0; k < k_HOST_VERTEX::dmn_size(); k++)

          for (int b2 = 0; b2 < b::dmn_size(); b2++)

            for (int b1 = 0; b1 < b::dmn_size(); b1++)

              diff += std::abs(

                  phase * leading_eigenvectors(i, b1, b2, k, w) -

                  conj(phase * leading_eigenvectors(i, b1, b2, k, w_VERTEX::dmn_size() - 1 - w)));

                  conj(phase * leading_eigenvectors(i, b1, b2, k, WVertexDmn::dmn_size() - 1 - w)));



      if (diff < diff_min) {

        diff_min = diff;
@@ -898,7 +898,7 @@ void BseLatticeSolver<ParametersType, DcaDataType, ScalarType>::characterizeLead 
    std::cout << "\n" << __FUNCTION__ << std::endl;



  for (int ev = 0; ev < num_evals; ev++) {

    for (int w = 0; w < w_VERTEX::dmn_size(); w++) {

    for (int w = 0; w < WVertexDmn::dmn_size(); w++) {

      for (int b2 = 0; b2 < b::dmn_size(); b2++) {

        for (int b1 = 0; b1 < b::dmn_size(); b1++) {

          for (int harmonic = 0; harmonic < num_harmonics; harmonic++) {