ABACUS/src/EXECS/Heis_Fourier_to_Sqt.cc

117 řádky
3.5 KiB
C++

/**********************************************************
This software is part of J.-S. Caux's ABACUS library.
Copyright (c) J.-S. Caux.
-----------------------------------------------------------
File: Heis_Fourier_to_Sqt.cc
Purpose: Fourier transform to q and t-dependent correlator for Heis
***********************************************************/
#include "ABACUS.h"
using namespace std;
using namespace ABACUS;
int main(int argc, char* argv[])
{
if (argc != 10) { // Print out instructions
cout << "Usage of Heis_Fourier_to_Sqt executable: " << endl << endl;
cout << "Provide arguments using one of the following options:" << endl << endl;
cout << "whichDSF Delta N M iKmin iKmax devmax Npts_t tmax" << endl << endl;
}
else {
int index = 1;
char whichDSF = *argv[index++];
DP Delta = atof(argv[index++]);
int N = atoi(argv[index++]);
int M = atoi(argv[index++]);
int iKmin = atoi(argv[index++]);
int iKmax = atoi(argv[index++]);
DP devmax = atof(argv[index++]);
int Npts_t = atoi(argv[index++]);
DP tmax = atof(argv[index++]);
stringstream filenameprefix;
Data_File_Name (filenameprefix, whichDSF, Delta, N, M, iKmin, iKmax, 0.0, 0, "");
string prefix = filenameprefix.str();
stringstream RAW_stringstream; string RAW_string;
RAW_stringstream << prefix << ".raw";
RAW_string = RAW_stringstream.str(); const char* RAW_Cstr = RAW_string.c_str();
ifstream RAW_infile;
RAW_infile.open(RAW_Cstr);
if (RAW_infile.fail()) {
cout << RAW_Cstr << endl;
ABACUSerror("Could not open RAW_infile... ");
}
RecMat<complex<DP> > FT(Npts_t + 1, iKmax - iKmin + 1);
DP omega;
int iK;
DP FF;
DP dev;
string label;
// Now define time coordinates: between 0 and tmax
Vect_DP tlattice(Npts_t + 1);
for (int i = 0; i <= Npts_t; ++i) tlattice[i] = i * tmax/Npts_t;
while (RAW_infile.peek() != EOF) {
RAW_infile >> omega >> iK >> FF >> dev >> label;
for (int it = 0; it < Npts_t; ++it)
FT[it][iK] += FF * FF * exp(II * omega * tlattice[it]);
}
RAW_infile.close();
// Output to files:
Write_K_File (N, iKmin, iKmax);
Write_Time_File (Npts_t + 1, 0.0, tmax);
stringstream FTre_stringstream; string FTre_string;
FTre_stringstream << prefix << "_Sqt_tmin_" << 0 << "_tmax_" << tmax << "_Nt_" << Npts_t << ".dat_re";
FTre_string = FTre_stringstream.str(); const char* FTre_Cstr = FTre_string.c_str();
ofstream FTre_outfile;
FTre_outfile.open(FTre_Cstr);
if (FTre_outfile.fail()) ABACUSerror("Could not open FTre_outfile... ");
stringstream FTim_stringstream; string FTim_string;
FTim_stringstream << prefix << "_Sqt_tmin_" << 0 << "_tmax_" << tmax << "_Nt_" << Npts_t << ".dat_im";
FTim_string = FTim_stringstream.str(); const char* FTim_Cstr = FTim_string.c_str();
ofstream FTim_outfile;
FTim_outfile.open(FTim_Cstr);
if (FTim_outfile.fail()) ABACUSerror("Could not open FTim_outfile... ");
for (int iK = iKmin; iK <= iKmax; ++iK) {
FTre_outfile << setprecision(16);
FTim_outfile << setprecision(16);
if (iK > iKmin) {
FTre_outfile << endl;
FTim_outfile << endl;
}
FTre_outfile << real(FT[0][iK]);
FTim_outfile << real(FT[0][iK]);
for (int it = 1; it <= Npts_t; ++it) {
FTre_outfile << "\t" << real(FT[it][iK]);
FTim_outfile << "\t" << imag(FT[it][iK]);
}
}
FTre_outfile.close();
FTim_outfile.close();
}
return(0);
}