mirror of
https://github.com/paboyle/Grid.git
synced 2025-06-07 09:57:07 +01:00
Compare commits
4 Commits
b517e88db3
...
dc546aaa4b
Author | SHA1 | Date | |
---|---|---|---|
|
dc546aaa4b | ||
|
5364d580c9 | ||
|
2a9a6347e3 | ||
|
cfdb56f314 |
@ -252,6 +252,11 @@ void WilsonFlow<Gimpl>::smear(GaugeField& out, const GaugeField& in) const{
|
||||
|
||||
out = in;
|
||||
RealD taus = 0.;
|
||||
|
||||
// Perform initial t=0 measurements
|
||||
for(auto const &meas : this->functions)
|
||||
meas.second(0,taus,out);
|
||||
|
||||
for (unsigned int step = 1; step <= Nstep; step++) { //step indicates the number of smearing steps applied at the time of measurement
|
||||
auto start = std::chrono::high_resolution_clock::now();
|
||||
evolve_step(out, taus);
|
||||
@ -336,6 +341,11 @@ void WilsonFlowAdaptive<Gimpl>::smear(GaugeField& out, const GaugeField& in) con
|
||||
RealD taus = 0.;
|
||||
RealD eps = init_epsilon;
|
||||
unsigned int step = 0;
|
||||
|
||||
// Perform initial t=0 measurements
|
||||
for(auto const &meas : this->functions)
|
||||
meas.second(step,taus,out);
|
||||
|
||||
do{
|
||||
int step_success = evolve_step_adaptive(out, taus, eps);
|
||||
step += step_success; //step will not be incremented if the integration step fails
|
||||
|
@ -66,6 +66,7 @@ namespace Grid{
|
||||
};
|
||||
}
|
||||
|
||||
|
||||
template <class T> void writeFile(T& in, std::string const fname){
|
||||
#ifdef HAVE_LIME
|
||||
// Ref: https://github.com/paboyle/Grid/blob/feature/scidac-wp1/tests/debug/Test_general_coarse_hdcg_phys48.cc#L111
|
||||
@ -107,8 +108,18 @@ int main(int argc, char **argv) {
|
||||
|
||||
for (int conf = CPar.StartConfiguration; conf <= CPar.EndConfiguration; conf+= CPar.Skip){
|
||||
|
||||
#if 0
|
||||
CPNersc.CheckpointRestore(conf, Umu, sRNG, pRNG);
|
||||
#else
|
||||
// Don't require Grid format RNGs
|
||||
FieldMetaData header;
|
||||
std::string file, filesmr;
|
||||
file = CPar.conf_path + "/" + CPar.conf_prefix + "." + std::to_string(conf);
|
||||
filesmr = CPar.conf_path + "/" + CPar.conf_smr_prefix + "." + std::to_string(conf);
|
||||
|
||||
NerscIO::readConfiguration(Umu,header,file);
|
||||
#endif
|
||||
|
||||
std::cout << std::setprecision(15);
|
||||
std::cout << GridLogMessage << "Initial plaquette: "<< WilsonLoops<PeriodicGimplR>::avgPlaquette(Umu) << std::endl;
|
||||
|
||||
@ -116,6 +127,7 @@ int main(int argc, char **argv) {
|
||||
std::string file_post = CPar.conf_prefix + "." + std::to_string(conf);
|
||||
|
||||
WilsonFlow<PeriodicGimplR> WF(WFPar.step_size,WFPar.steps,WFPar.meas_interval);
|
||||
|
||||
WF.addMeasurement(WFPar.meas_interval_density, [&file_pre,&file_post,&conf](int step, RealD t, const typename PeriodicGimplR::GaugeField &U){
|
||||
|
||||
typedef typename PeriodicGimplR::GaugeLinkField GaugeMat;
|
||||
@ -165,33 +177,48 @@ int main(int argc, char **argv) {
|
||||
//double coeff = 2.0 / (1.0 * Nd * (Nd - 1)) / 3.0;
|
||||
//Plq = coeff * Plq;
|
||||
|
||||
int tau = std::round(t);
|
||||
std::string efile = file_pre + "E_dnsty_" + std::to_string(tau) + "_" + file_post;
|
||||
writeFile(R,efile);
|
||||
std::string tfile = file_pre + "Top_dnsty_" + std::to_string(tau) + "_" + file_post;
|
||||
writeFile(qfield,tfile);
|
||||
|
||||
RealD WFlow_TC5Li = WilsonLoops<PeriodicGimplR>::TopologicalCharge5Li(U);
|
||||
|
||||
int tau = std::round(t);
|
||||
|
||||
std::string efile = file_pre + "E_dnsty_" + std::to_string(tau) + "_" + file_post;
|
||||
// writeFile(R,efile);
|
||||
|
||||
std::string tfile = file_pre + "Top_dnsty_" + std::to_string(tau) + "_" + file_post;
|
||||
// writeFile(qfield,tfile);
|
||||
|
||||
std::string ufile = file_pre + "U_" + std::to_string(tau) + "_" + file_post;
|
||||
{
|
||||
PeriodicGimplR::GaugeField Ucopy = U;
|
||||
NerscIO::writeConfiguration(Ucopy,ufile);
|
||||
}
|
||||
|
||||
RealD E = real(sum(R))/ RealD(U.Grid()->gSites());
|
||||
RealD T = real( sum(qfield) );
|
||||
Coordinate scoor; for (int mu=0; mu < Nd; mu++) scoor[mu] = 0;
|
||||
RealD E0 = real(peekSite(R,scoor));
|
||||
RealD T0 = real(peekSite(qfield,scoor));
|
||||
std::cout << GridLogMessage << "[WilsonFlow] Saved energy density (clover) & topo. charge density: " << conf << " " << step << " " << tau << " "
|
||||
<< "(E_avg,T_sum) " << E << " " << T << " (E, T at origin) " << E0 << " " << T0 << std::endl;
|
||||
<< "(E_avg,T_sum) " << E << " " << T << " (E, T at origin) " << E0 << " " << T0 << " Q5Li "<< WFlow_TC5Li << std::endl;
|
||||
|
||||
});
|
||||
|
||||
int t=WFPar.maxTau;
|
||||
WF.smear(Uflow, Umu);
|
||||
|
||||
NerscIO::writeConfiguration(Uflow,filesmr);
|
||||
|
||||
|
||||
RealD WFlow_plaq = WilsonLoops<PeriodicGimplR>::avgPlaquette(Uflow);
|
||||
RealD WFlow_TC = WilsonLoops<PeriodicGimplR>::TopologicalCharge(Uflow);
|
||||
RealD WFlow_TC5Li = WilsonLoops<PeriodicGimplR>::TopologicalCharge5Li(Uflow);
|
||||
RealD WFlow_T0 = WF.energyDensityPlaquette(t,Uflow); // t
|
||||
RealD WFlow_EC = WF.energyDensityCloverleaf(t,Uflow);
|
||||
std::cout << GridLogMessage << "Plaquette "<< conf << " " << WFlow_plaq << std::endl;
|
||||
std::cout << GridLogMessage << "T0 "<< conf << " " << WFlow_T0 << std::endl;
|
||||
std::cout << GridLogMessage << "TC0 "<< conf << " " << WFlow_EC << std::endl;
|
||||
std::cout << GridLogMessage << "TopologicalCharge "<< conf << " " << WFlow_TC << std::endl;
|
||||
std::cout << GridLogMessage << "Plaquette "<< conf << " " << WFlow_plaq << std::endl;
|
||||
std::cout << GridLogMessage << "T0 "<< conf << " " << WFlow_T0 << std::endl;
|
||||
std::cout << GridLogMessage << "TC0 "<< conf << " " << WFlow_EC << std::endl;
|
||||
std::cout << GridLogMessage << "TopologicalCharge "<< conf << " " << WFlow_TC << std::endl;
|
||||
std::cout << GridLogMessage << "TopologicalCharge5Li "<< conf << " " << WFlow_TC5Li<< std::endl;
|
||||
|
||||
std::cout<< GridLogMessage << " Admissibility check:\n";
|
||||
const double sp_adm = 0.067; // admissible threshold
|
||||
|
@ -1,3 +1,12 @@
|
||||
|
||||
spack load c-lime
|
||||
spack load fftw
|
||||
spack load hdf5+cxx
|
||||
|
||||
export FFTW=`spack find --paths fftw | grep ^fftw | awk '{print $2}' `
|
||||
export HDF5=`spack find --paths hdf5+cxx | grep ^hdf5 | awk '{print $2}' `
|
||||
export CLIME=`spack find --paths c-lime | grep ^c-lime | awk '{print $2}' `
|
||||
|
||||
../../configure \
|
||||
--enable-comms=mpi-auto \
|
||||
--enable-unified=yes \
|
||||
@ -5,12 +14,16 @@
|
||||
--enable-shm-fast-path=shmopen \
|
||||
--enable-accelerator=none \
|
||||
--enable-simd=AVX512 \
|
||||
--disable-accelerator-cshift \
|
||||
--with-lime=$CLIME \
|
||||
--with-hdf5=$HDF5 \
|
||||
--with-fftw=$FFTW \
|
||||
--disable-fermion-reps \
|
||||
--disable-gparity \
|
||||
CXX=clang++ \
|
||||
MPICXX=mpicxx \
|
||||
CXXFLAGS="-std=c++17"
|
||||
LIBS=-llime \
|
||||
LDFLAGS=-L$CLIME/lib/ \
|
||||
CXXFLAGS="-std=c++17 -fPIE"
|
||||
|
||||
|
||||
|
||||
|
@ -1,4 +1,5 @@
|
||||
source $HOME/spack/share/spack/setup-env.sh
|
||||
spack load llvm@17.0.4
|
||||
export LD_LIBRARY_PATH=/direct/sdcc+u/paboyle/spack/opt/spack/linux-almalinux8-icelake/gcc-8.5.0/llvm-17.0.4-laufdrcip63ivkadmtgoepwmj3dtztdu/lib:$LD_LIBRARY_PATH
|
||||
module load openmpi
|
||||
module load openmpi/4.1.8
|
||||
spack load c-lime
|
||||
|
@ -31,12 +31,23 @@ directory
|
||||
|
||||
using namespace std;
|
||||
using namespace Grid;
|
||||
;
|
||||
|
||||
//typedef WilsonFermionD FermionOp;
|
||||
typedef DomainWallFermionD FermionOp;
|
||||
typedef typename DomainWallFermionD::FermionField FermionField;
|
||||
|
||||
template <class T> void writeFile(T& in, std::string const fname){
|
||||
#ifdef HAVE_LIME
|
||||
// Ref: https://github.com/paboyle/Grid/blob/feature/scidac-wp1/tests/debug/Test_general_coarse_hdcg_phys48.cc#L111
|
||||
std::cout << Grid::GridLogMessage << "Writes to: " << fname << std::endl;
|
||||
Grid::emptyUserRecord record;
|
||||
Grid::ScidacWriter WR(in.Grid()->IsBoss());
|
||||
WR.open(fname);
|
||||
WR.writeScidacFieldRecord(in,record,0);
|
||||
WR.close();
|
||||
#endif
|
||||
// What is the appropriate way to throw error?
|
||||
}
|
||||
|
||||
RealD AllZero(RealD x) { return 0.; }
|
||||
|
||||
@ -121,7 +132,7 @@ int main(int argc, char** argv) {
|
||||
|
||||
int Ls=16;
|
||||
RealD M5=1.8;
|
||||
RealD mass = -1.0;
|
||||
RealD mass = 0.01;
|
||||
|
||||
mass=LanParams.mass;
|
||||
Ls=LanParams.Ls;
|
||||
@ -159,10 +170,10 @@ int main(int argc, char** argv) {
|
||||
U[mu] = PeekIndex<LorentzIndex>(Umu, mu);
|
||||
}
|
||||
*/
|
||||
|
||||
int Nstop = 10;
|
||||
int Nk = 20;
|
||||
int Nstop = Nk;
|
||||
int Np = 80;
|
||||
|
||||
Nstop=LanParams.Nstop;
|
||||
Nk=LanParams.Nk;
|
||||
Np=LanParams.Np;
|
||||
@ -201,8 +212,12 @@ int main(int argc, char** argv) {
|
||||
int Nconv;
|
||||
IRL.calc(eval, evec, src, Nconv);
|
||||
|
||||
std::cout << mass <<" : " << eval << std::endl;
|
||||
|
||||
std::cout << mass <<" : " << eval << std::endl;
|
||||
std::cout << " #evecs " << evec.size() << std::endl;
|
||||
std::cout << " Nconv " << Nconv << std::endl;
|
||||
std::cout << " Nm " << Nm << std::endl;
|
||||
if ( Nconv > evec.size() ) Nconv = evec.size();
|
||||
|
||||
#if 0
|
||||
Gamma g5(Gamma::Algebra::Gamma5) ;
|
||||
ComplexD dot;
|
||||
@ -232,6 +247,7 @@ int main(int argc, char** argv) {
|
||||
vector<LatticeFermion> finalevec(Nconv, FGrid);
|
||||
vector<RealD> eMe(Nconv), eMMe(Nconv);
|
||||
for(int i = 0; i < Nconv; i++){
|
||||
cout << "calculate the matrix element["<<i<<"]" << endl;
|
||||
G5R5Herm.HermOpAndNorm(evec[i], G5R5Mevec[i], eMe[i], eMMe[i]);
|
||||
}
|
||||
cout << "Re<evec, G5R5M(evec)>: " << endl;
|
||||
@ -326,13 +342,27 @@ int main(int argc, char** argv) {
|
||||
axpby_ssp(G5evec[i], -1., finalevec[i], 0., G5evec[i], j, j);
|
||||
}
|
||||
}
|
||||
|
||||
for(int i = 0; i < Nconv; i++){
|
||||
chiral_matrix_real[i].resize(Nconv);
|
||||
chiral_matrix[i].resize(Nconv);
|
||||
|
||||
std::string evfile("./evec_density");
|
||||
evfile = evfile+"_"+std::to_string(i);
|
||||
auto evdensity = localInnerProduct(finalevec[i],finalevec[i] );
|
||||
writeFile(evdensity,evfile);
|
||||
|
||||
for(int j = 0; j < Nconv; j++){
|
||||
chiral_matrix[i][j] = innerProduct(finalevec[i], G5evec[j]);
|
||||
std::cout <<" chiral_matrix_real signed "<<i<<" "<<j<<" "<< chiral_matrix_real[i][j] << std::endl;
|
||||
chiral_matrix_real[i][j] = abs(chiral_matrix[i][j]);
|
||||
std::cout <<" chiral_matrix_real "<<i<<" "<<j<<" "<< chiral_matrix_real[i][j] << std::endl;
|
||||
if ( chiral_matrix_real[i][j] > 0.8 ) {
|
||||
auto g5density = localInnerProduct(finalevec[i], G5evec[j]);
|
||||
std::string chfile("./chiral_density_");
|
||||
chfile = chfile +std::to_string(i)+"_"+std::to_string(j);
|
||||
writeFile(g5density,chfile);
|
||||
}
|
||||
}
|
||||
}
|
||||
for(int i = 0; i < Nconv; i++){
|
||||
@ -341,6 +371,43 @@ int main(int argc, char** argv) {
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
FILE *fp = fopen("lego-plot.py","w"); assert(fp!=NULL);
|
||||
#define PYTHON_LINE(A) fprintf(fp,A"\n");
|
||||
PYTHON_LINE("import matplotlib.pyplot as plt");
|
||||
PYTHON_LINE("import numpy as np");
|
||||
PYTHON_LINE("");
|
||||
PYTHON_LINE("fig = plt.figure()");
|
||||
PYTHON_LINE("ax = fig.add_subplot(projection='3d')");
|
||||
PYTHON_LINE("");
|
||||
PYTHON_LINE("x, y = np.random.rand(2, 100) * 4");
|
||||
PYTHON_LINE("hist, xedges, yedges = np.histogram2d(x, y, bins=10, range=[[0, 9], [0, 9]])");
|
||||
PYTHON_LINE("");
|
||||
PYTHON_LINE("# Construct arrays for the anchor positions of the 16 bars");
|
||||
PYTHON_LINE("xpos, ypos = np.meshgrid(xedges[:-1] + 0.25, yedges[:-1] + 0.25, indexing=\"ij\")");
|
||||
PYTHON_LINE("xpos = xpos.ravel()");
|
||||
PYTHON_LINE("ypos = ypos.ravel()");
|
||||
PYTHON_LINE("zpos = 0");
|
||||
PYTHON_LINE("");
|
||||
PYTHON_LINE("# Construct arrays with the dimensions for the 16 bars.");
|
||||
PYTHON_LINE("dx = dy = 0.5 * np.ones_like(zpos)");
|
||||
PYTHON_LINE("dz = np.array([");
|
||||
for(int i = 0; i < Nconv; i++){
|
||||
fprintf(fp,"\t[ ");
|
||||
for(int j = 0; j < Nconv; j++){
|
||||
fprintf(fp,"%lf ",chiral_matrix_real[i][j]);
|
||||
if(j<Nconv-1) fprintf(fp,",");
|
||||
else fprintf(fp," ");
|
||||
}
|
||||
fprintf(fp,"]");
|
||||
if(i<Nconv-1) fprintf(fp,",\n");
|
||||
else fprintf(fp,"\n");
|
||||
}
|
||||
|
||||
PYTHON_LINE("\t])");
|
||||
PYTHON_LINE("dz = dz.ravel()");
|
||||
PYTHON_LINE("ax.bar3d(xpos, ypos, zpos, dx, dy, dz, zsort='average')");
|
||||
PYTHON_LINE("plt.show()");
|
||||
fclose(fp);
|
||||
|
||||
Grid_finalize();
|
||||
}
|
||||
|
Loading…
x
Reference in New Issue
Block a user