Fluc_randHOp.cpp File Reference

Classes
struct	statm_t

Functions
void	read_off_memory_status (statm_t &result)
int	main (int argc, char **argv)

Variables
namespace	filesystem = std::experimental::filesystem
statm_t	Usage

Function Documentation

main()

int main	(	int	argc,
		char **	argv
	)

                                {
  #define USAGE_commom "Usage: 1.This 2.n_max 3.n_min 4.division No. 5.Output dir "
  #define Nargs_common 5
  #define USAGE_opt "(Opt1).dE (Opt2).repMin (Opt3).repMax "
  #define Nargs_opt 3
  #include "PBC_TI/Framework/Fragments/setVariablesForEnsemble.cpp"
  #include "PBC_TI/Framework/Fragments/setVariablesForMCAverage.cpp"
 
  Integer_t dim_sub; //全ヒルベルト空間の次元
  Integer_t ndata, info, failed=0; // カウンタ
  double gE, EnergyRange, OpRange, OpMin;
  double Fluc, MaxFluc, MaxPos; // ETHの指標
  debug_print("# Successfully initialized.");
  //******************** Check for the directory structure ********************
  debug_print("# Checking for the directory structure.");
  #include "PBC_TI/Framework/Fragments/MeasureOfETH_CreateOutputFiles.cpp"
  std::ofstream OutSample;
  //******************** (END) Check for the directory structure ********************
 
  //***************************************************************************
  //******************** Allocation & Initialization **************************
  //***************************************************************************
#ifdef GPU
  magma_init();
  magma_queue_t queue = NULL;
  magma_int_t     dev = 0;
  magma_getdevice( &dev );
  magma_queue_create( dev, &queue );
#else
  void* GPUconf = nullptr;
#endif
 
  //******************** Translation invariance ********************
  debug_print("# Calculating translation-invariant sectors.");
  #include "PBC_TI/Framework/Fragments/TranslationInvariantSectors.cpp"
  for(Integer_t n = n_max;n >= n_min; --n) {
    OutMeasure[n] << "# (n=" << n << ") dim=" << Sector[n].dim() << "\n"
                  << "# 1.(Sample No.) 2.(Fluc_2) 3.(Fluc_infty) 4.(Pos_max_flux) 5.(ndata)\n"
                  << std::endl << std::scientific;
  }
  //******************** (END)Translation invariance ********************
 
  //********** Allocate CPU memories **********//
  debug_print("# Allocating CPU memories.");
  Integer_t const dim_max = Sector[n_max].dim();
  std::vector<Integer_t>   NdataInShell(dim_max);
  std::vector<double>         eigenEnergy(dim_max);
  std::vector<double>            EXPvalue(dim_max);
  std::vector<double>               MCAverage(dim_max);
#ifndef GPU
  matrix<Complex_t>    h_tot(dim_max, dim_max);
  matrix<Complex_t>  loc_tot(dim_max, dim_max);
  #define   dh_tot   h_tot
  #define dloc_tot loc_tot
#endif
  //********** (END) Allocate CPU memories **********//
 
#ifdef GPU
  //********** Allocate GPU memories **********//
  debug_print("# Allocating GPU memories.");
  constexpr Integer_t GPU_UNIT  = 32;
  Integer_t const LDT = magma_roundup(dim_max, GPU_UNIT);
  matrix_gpu<Complex_t>   dh_tot(LDT, dim_max);
  matrix_gpu<Complex_t> dloc_tot(LDT, dim_max);
  //********** (END) Allocate GPU memories **********//
 
  //********** Determine GPU configuration **********//
  #include "PBC_TI/Framework/Fragments/getAttributesOfMatrixElementsInSector.cpp"
  //********** (END) Determine GPU configuration **********//
#endif // #ifdef GPU
 
  double start, t_int, end, temp_t;
  double T_diag=0, T_post=0, T_pre=0;
  init_genrand(SEED);
  start = getETtime();
  for(Integer_t repetition = 0;repetition < repMin; ++repetition) {
    for(Integer_t n = n_min;n <= n_max; ++n) {
      RandomMatrix_GUE(h_tot,   Sector[n].dim());
      RandomMatrix_GUE(loc_tot, Sector[n].dim());
    }
  }
  end = getETtime();
  std::cout << "(init_genrand): time=" << std::fixed << (end-start) << std::endl;
  //***************************************************************************
  //******************** (END) Allocation & Initialization ********************
  //***************************************************************************
 
  start = getETtime();
  end   = start;
  for(Integer_t repetition = repMin;repetition <= repMax; ++repetition) {
    // 局所ハミルトニアンと局所物理量をランダムにとる ******************************//
  #if !defined(NDEBUG) && defined(__linux__)
    read_off_memory_status(Usage);
        std::cerr << "# (GPU, rep= " << std::setw(6) << repetition
        << ") size="       << Usage.size
        << ", resident="   << Usage.resident
        << ", share="      << Usage.share
        << ", text="       << Usage.text
        << ", lib="        << Usage.lib
        << ", data="       << Usage.data
        << ", dt="         << Usage.dt
        << std::endl;
  #endif
 
    if(repetition < 10000) {
      #include "PBC_TI/Framework/Fragments/createFilesForSampleData.cpp"
    }
    for(Integer_t n = n_max;n >= n_min; --n) {
      debug_print("# (rep,n)=(" << repetition << "," << n << ")");
      debug_print("# Constructing global matrices in the sector.");
      dim_sub = Sector[n].dim();
      temp_t  = getETtime();
      {
          RandomMatrix_GUE(h_tot,   dim_sub);
          RandomMatrix_GUE(loc_tot, dim_sub);
        #ifdef GPU
          magma_setmatrix(dim_sub, dim_sub, sizeof(Complex_t),   &*h_tot.begin(), dim_sub,   dh_tot.ptr(), dim_sub, GPUconf.queue());
          magma_setmatrix(dim_sub, dim_sub, sizeof(Complex_t), &*loc_tot.begin(), dim_sub, dloc_tot.ptr(), dim_sub, GPUconf.queue());
        #endif
      }
      T_pre += getETtime() -temp_t;
 
      debug_print("# Calculating eigenstate expectation values.");
      temp_t = getETtime();
        info = EigenExpValue(EXPvalue, eigenEnergy, dim_sub, dh_tot, dloc_tot, GPUconf);
          if(info != 0) {
            failed += 1;
            std::stringstream buff("");
            buff << "/Failed "  << PRECISION << "_N" << n << "_No" << repetition << ".txt";
            std::string filename(baseDirName); filename += buff.str();
            std::ofstream ErrFs;
            ErrFs.open(filename);
            if( !ErrFs.is_open() ) {
              std::cerr << "# Warning: Can't open a file " << filename << ". Continue..." << std::endl;
            } else {
              ErrFs << "# info=" << info << "\n\n"
                       "# N=" << n << ", repetition=" << repetition << "\n";
              if(dim_sub < 20) {
                ErrFs << "# Total Hamiltonian\n";
                dh_tot.print(  ErrFs,dim_sub,dim_sub,"All", GPUconf);
                ErrFs << "# Total operator\n";
                dloc_tot.print(ErrFs,dim_sub,dim_sub,"All", GPUconf);
              }
              ErrFs.close();
            }
            continue;
          }
        debug_print("# Calculating the spectral range of the observable.");
        OpRange = SpectralRange(OpMin, dim_sub, dloc_tot);
      T_diag += getETtime() -temp_t;
 
      debug_print("# Calculating measures of the ETH.");
      temp_t  = getETtime();
      {
        EnergyRange = eigenEnergy[dim_sub-1] -eigenEnergy[0];
        gE = eigenEnergy[0];
        for(Integer_t i = 0;i < dim_sub; ++i) {
          eigenEnergy[i]  = (eigenEnergy[i] -gE)/EnergyRange;
             EXPvalue[i] /= OpRange;
        }
 
        // Calculate Microcanonical averages with energy corresponding to each eigenstate.
        MCAverage.resize(dim_sub); NdataInShell.resize(dim_sub);
        std::fill(   MCAverage.begin(),    MCAverage.end(), 0.0);
        std::fill(NdataInShell.begin(), NdataInShell.end(), 0.0);
        for(Integer_t i = 0;i < dim_sub; ++i) {
          MCAverage[i] = MicroCanonicalAverage(NdataInShell[i], eigenEnergy[i], dE, dim_sub, eigenEnergy, EXPvalue, i);
        }
 
        ndata = MeasureOfETH(Fluc, MaxFluc, MaxPos, Emin, Emax, dim_sub, eigenEnergy, EXPvalue, MCAverage);
      }
      T_post += getETtime() -temp_t;
 
      OutMeasure[n] << std::setw(6)  << repetition << " " << std::showpos
                    << std::setw(13) << Fluc       << " "
                    << std::setw(13) << MaxFluc    << " "
                    << std::setw(13) << MaxPos     << " " << std::noshowpos << ndata << "\n";
      if(repetition < 10000) {
        for(Integer_t i = 0;i < dim_sub; ++i) {
          OutSample << std::setw(2)  << n << " " << std::showpos
                    << std::setw(13) << eigenEnergy[i]*EnergyRange+gE << " "
                    << std::setw(13) << eigenEnergy[i] << " "
                    << std::setw(13) << EXPvalue[i]    << "\n" << std::noshowpos;
        }
      }
    }
    if(repetition < 10000) OutSample.close();
    if(repetition%10 == 9) {
      t_int = end; end = getETtime();
      std::cerr << "(total=" << std::setw(6) << repetition+1
                << "): timeINT="    << std::setprecision(6)  << std::setw(8)  << (end-t_int)
                << ", timeTOT="     << std::setprecision(6)  << std::setw(8)  << (end-start)
                << ", T_construct=" << std::setprecision(6)  << std::setw(10) << T_pre  << "(" << std::setprecision(1) << 100*T_pre /(end-start) << "%)"
                << ", T_diag="      << std::setprecision(6)  << std::setw(8)  << T_diag << "(" << std::setprecision(1) << 100*T_diag/(end-start) << "%)"
                << ", T_process="   << std::setprecision(6)  << std::setw(8)  << T_post << "(" << std::setprecision(1) << 100*T_post/(end-start) << "%)"
                << std::endl;
      for(Integer_t n = n_min;n <= n_max; ++n) OutMeasure[n].flush();
    }
  }
 
  return 0;
}

read_off_memory_status()

void read_off_memory_status

(

statm_t &

result

)

  {
    const char* statm_path = "/proc/self/statm";
 
    FILE *f = fopen(statm_path,"r");
    if(!f){
      perror(statm_path);
      abort();
    }
    if(7 != fscanf(f,"%ld %ld %ld %ld %ld %ld %ld",
      &result.size,&result.resident,&result.share,&result.text,&result.lib,&result.data,&result.dt))
    {
      perror(statm_path);
      abort();
    }
    fclose(f);
  }

Variable Documentation

filesystem

namespace filesystem = std::experimental::filesystem

Usage

statm_t Usage