Skip to content

Generate HPhi chain momentum TransSym input #59

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
tmisawa merged 2 commits into from
Jul 7, 2026
Merged

Generate HPhi chain momentum TransSym input #59

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
fae83da
Generate HPhi chain momentum TransSym input
tmisawa Jul 7, 2026
a7e84c5
Reject unsupported HPhi momentum standard terms
tmisawa Jul 7, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
109 changes: 109 additions & 0 deletions src/StdFace_main.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1078,6 +1078,7 @@ static void StdFace_ResetVals(struct StdIntList *StdI) {
StdI->exct = StdI->NaN_i;
StdI->LanczosEps = StdI->NaN_i;
StdI->LanczosTarget = StdI->NaN_i;
StdI->MomentumIndex = StdI->NaN_i;
StdI->NumAve = StdI->NaN_i;
StdI->ExpecInterval = StdI->NaN_i;
StdI->dt = NaN_d;
Expand Down Expand Up @@ -1432,6 +1433,109 @@ static void PrintTrans(struct StdIntList *StdI){
fclose(fp);
fprintf(stdout, " trans.def is written.\n");
}/*static void PrintTrans*/
#if defined(_HPhi)
static int StdFace_UsesMomentumSymmetry(struct StdIntList *StdI)
{
return StdI->MomentumIndex != StdI->NaN_i;
}

static int StdFace_IsChainLattice(struct StdIntList *StdI)
{
return strcmp(StdI->lattice, "chain") == 0 ||
strcmp(StdI->lattice, "chainlattice") == 0;
}

static int StdFace_HasNonZeroRealTerms(double *terms, int nterms, double eps)
{
int iterm;
for (iterm = 0; iterm < nterms; iterm++) {
if (fabs(terms[iterm]) > eps) return 1;
}
return 0;
}

static int StdFace_HasNonZeroComplexTerms(double complex *terms, int nterms, double eps)
{
int iterm;
for (iterm = 0; iterm < nterms; iterm++) {
if (cabs(terms[iterm]) > eps) return 1;
}
return 0;
}

static void CheckMomentumSymmetry(struct StdIntList *StdI)
{
double eps = 1.0e-12;
int idim;
if (!StdFace_UsesMomentumSymmetry(StdI)) return;
if (!StdFace_IsChainLattice(StdI)) {
fprintf(stdout, "\n ERROR ! MomentumIndex currently supports only chain lattice.\n");
StdFace_exit(-1);
}
if (strcmp(StdI->model, "spin") != 0 || StdI->lGC != 0) {
fprintf(stdout, "\n ERROR ! MomentumIndex currently supports only canonical Spin model.\n");
StdFace_exit(-1);
}
if (StdI->S2 != 1) {
fprintf(stdout, "\n ERROR ! MomentumIndex currently supports only Spin-1/2.\n");
StdFace_exit(-1);
}
if (StdI->L <= 0) {
fprintf(stdout, "\n ERROR ! MomentumIndex requires positive L.\n");
StdFace_exit(-1);
}
if (StdI->MomentumIndex < 0 || StdI->MomentumIndex >= StdI->L) {
fprintf(stdout, "\n ERROR ! MomentumIndex must satisfy 0 <= MomentumIndex < L.\n");
StdFace_exit(-1);
}
for (idim = 0; idim < 3; idim++) {
if (fabs(StdI->phase[idim]) > eps || StdI->AntiPeriod[idim] != 0) {
fprintf(stdout, "\n ERROR ! MomentumIndex does not support boundary phase in v1.\n");
StdFace_exit(-1);
}
}
if (StdFace_HasNonZeroComplexTerms(StdI->trans, StdI->ntrans, eps) ||
StdFace_HasNonZeroComplexTerms(StdI->intr, StdI->nintr, eps) ||
StdFace_HasNonZeroRealTerms(StdI->Cinter, StdI->NCinter, eps) ||
StdFace_HasNonZeroRealTerms(StdI->Hund, StdI->NHund, eps) ||
StdFace_HasNonZeroRealTerms(StdI->PairLift, StdI->NPairLift, eps) ||
StdFace_HasNonZeroRealTerms(StdI->PairHopp, StdI->NPairHopp, eps)) {
fprintf(stdout, "\n ERROR ! MomentumIndex currently supports only exchange-only Spin-1/2 chain with Jz = 0 and no field/general/pair terms.\n");
StdFace_exit(-1);
}
}

static void PrintMomentumTransSym(struct StdIntList *StdI)
{
FILE *fp;
int op, site;
if (!StdFace_UsesMomentumSymmetry(StdI)) return;

fp = fopen("qptransidx.def", "w");
fprintf(fp, "=============================================\n");
fprintf(fp, "NQPTrans %10d\n", StdI->L);
fprintf(fp, "=============================================\n");
fprintf(fp, "======== TrIdx_TrWeight_and_TrIdx_i_xi ======\n");
fprintf(fp, "=============================================\n");
for (op = 0; op < StdI->L; op++) {
double angle = -2.0 * StdI->pi * (double)StdI->MomentumIndex * (double)op / (double)StdI->L;
double re = cos(angle);
double im = sin(angle);
if (fabs(re) < 1.0e-14) re = 0.0;
if (fabs(im) < 1.0e-14) im = 0.0;
fprintf(fp, "%d %25.15f %25.15f\n", op, re, im);
}
for (op = 0; op < StdI->L; op++) {
for (site = 0; site < StdI->L; site++) {
fprintf(fp, "%5d %5d %5d %5d\n", op, site, (site + op) % StdI->L, 1);
}
}

fflush(fp);
fclose(fp);
fprintf(stdout, " qptransidx.def is written for MomentumIndex = %d.\n", StdI->MomentumIndex);
}/*static void PrintMomentumTransSym*/
#endif
/**
* @brief Print namelist.def
* @author Mitsuaki Kawamura (The University of Tokyo)
Expand Down Expand Up @@ -1471,6 +1575,8 @@ static void PrintNamelist(struct StdIntList *StdI){
else if (StdI->PumpBody == 2)
fprintf(fp, " TETwoBody tetwo.def\n");
}/*if (strcmp(StdI->method, "timeevolution") == 0)*/
if (StdFace_UsesMomentumSymmetry(StdI))
fprintf( fp, " TransSym qptransidx.def\n");
fprintf( fp, " SpectrumVec %s_eigenvec_0\n",
StdI->CDataFileHead);
if (StdI->lBoost == 1) fprintf( fp, " Boost boost.def\n");
Expand Down Expand Up @@ -2819,6 +2925,7 @@ void StdFace_main(
else if (strcmp(keyword, "lanczos_max") == 0) StoreWithCheckDup_i(keyword, value, &StdI->Lanczos_max);
else if (strcmp(keyword, "largevalue") == 0) StoreWithCheckDup_d(keyword, value, &StdI->LargeValue);
else if (strcmp(keyword, "method") == 0) StoreWithCheckDup_sl(keyword, value, StdI->method);
else if (strcmp(keyword, "momentumindex") == 0) StoreWithCheckDup_i(keyword, value, &StdI->MomentumIndex);
else if (strcmp(keyword, "nomega") == 0) StoreWithCheckDup_i(keyword, value, &StdI->Nomega);
else if (strcmp(keyword, "numave") == 0) StoreWithCheckDup_i(keyword, value, &StdI->NumAve);
else if (strcmp(keyword, "nvec") == 0) StoreWithCheckDup_i(keyword, value, &StdI->nvec);
Expand Down Expand Up @@ -3021,6 +3128,7 @@ void StdFace_main(
/**/
#if defined(_HPhi)
StdFace_LargeValue(StdI);
CheckMomentumSymmetry(StdI);
/*
Generate Hamiltonian for Boost
*/
Expand All @@ -3045,6 +3153,7 @@ void StdFace_main(
PrintLocSpin(StdI);
PrintTrans(StdI);
PrintInteractions(StdI);
PrintMomentumTransSym(StdI);
CheckModPara(StdI);
PrintModPara(StdI);

Expand Down
1 change: 1 addition & 0 deletions src/StdFace_vals.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -289,6 +289,7 @@ struct StdIntList {
int exct;/**<@brief The number of eigenvectors to be computed. input from file*/
int LanczosEps;/**<@brief Convergence threshold for the Lanczos method.*/
int LanczosTarget;/**<@brief Which eigenvector is used for the convergence check.*/
int MomentumIndex;/**<@brief Momentum index for HPhi TransSym generation.*/
int NumAve;/**<@brief Number of trials for TPQ calculation.*/
int ExpecInterval;/**<@brief Interval for the iteration when the expectation
value is computed.*/
Expand Down
Loading