Skip to content

Commit ec76661

Browse files
k-yoshimik-yoshimi
authored
Merge pull request #164 from issp-center-dev/ThreeBodyForGeneralSpin
Three body for general spin
2 parents b560cc8 + 007849b commit ec76661

File tree

8 files changed

+599
-14
lines changed

8 files changed

+599
-14
lines changed

‎.github/workflows/main.yml

Lines changed: 76 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -75,4 +75,79 @@ jobs:
7575
- name: ctest
7676
working-directory: ${{runner.workspace}}/build
7777
shell: bash
78-
run: ctest -V
78+
run: ctest -V -E spinone
79+
80+
ctest_spinone:
81+
runs-on: ${{ matrix.os }}
82+
83+
strategy:
84+
matrix:
85+
os: [ubuntu-22.04]
86+
mpisize: [1, 3, 9]
87+
ompsize: [1, 3]
88+
include:
89+
- os: macos-13
90+
mpisize: 1
91+
ompsize: 1
92+
- os: macos-latest
93+
mpisize: 1
94+
ompsize: 1
95+
- os: ubuntu-20.04
96+
mpisize: 1
97+
ompsize: 1
98+
exclude:
99+
- mpisize: 9
100+
ompsize: 3
101+
- mpisize: 3
102+
ompsize: 3
103+
fail-fast: false
104+
105+
env:
106+
MPIRUN: "mpiexec --oversubscribe -np ${{ matrix.mpisize }}"
107+
OMP_NUM_THREADS: ${{ matrix.ompsize }}
108+
109+
steps:
110+
- uses: actions/checkout@v4
111+
112+
- name: apt
113+
if: ${{ runner.os == 'Linux' }}
114+
run: |
115+
sudo apt update
116+
sudo apt install liblapack-dev openmpi-bin libopenmpi-dev libscalapack-openmpi-dev
117+
118+
- name: brew
119+
if: ${{ runner.os == 'macOS' }}
120+
run: |
121+
brew install openmpi scalapack
122+
123+
- name: Setup Python
124+
uses: actions/setup-python@v5
125+
with:
126+
python-version: "3.11"
127+
128+
- name: pip
129+
run: |
130+
python3 -m pip install numpy
131+
132+
- name: make workspace
133+
run: cmake -E make_directory ${{runner.workspace}}/build
134+
135+
- name: cmake
136+
working-directory: ${{runner.workspace}}/build
137+
shell: bash
138+
run: |
139+
if [ ${{ runner.os }} = "macOS" ] ; then
140+
export FC=gfortran-12
141+
fi
142+
cmake -DCMAKE_VERBOSE_MAKEFILE=ON $GITHUB_WORKSPACE
143+
144+
- name: build
145+
working-directory: ${{runner.workspace}}/build
146+
shell: bash
147+
run: cmake --build ./ -j4
148+
149+
- name: ctest
150+
working-directory: ${{runner.workspace}}/build
151+
shell: bash
152+
run: ctest -V -R spinone
153+

‎src/expec_cisajscktaltdc.c

Lines changed: 178 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -771,6 +771,170 @@ int expec_cisajscktalt_SpinHalf(struct BindStruct *X,double complex *vec, FILE *
771771
return 0;
772772
}
773773

774+
/**
775+
* @brief Child function to calculate three-body green functions for general Spin GC model
776+
*
777+
* @param X [in] data list for calculation
778+
* @param vec [in] eigenvectors
779+
* @param _fp [in] output file name
780+
* @retval 0 normally finished
781+
* @retval -1 abnormally finished
782+
*
783+
*/
784+
785+
int expec_Threebody_SpinGeneral(struct BindStruct *X,double complex *vec, FILE **_fp){
786+
long unsigned int i,j;
787+
long unsigned int org_isite1,org_isite2,org_isite3,org_isite4,org_isite5,org_isite6;
788+
long unsigned int org_sigma1,org_sigma2,org_sigma3,org_sigma4,org_sigma5,org_sigma6;
789+
long unsigned int tmp_org_isite1,tmp_org_isite2,tmp_org_isite3,tmp_org_isite4,tmp_org_isite5,tmp_org_isite6;
790+
long unsigned int tmp_org_sigma1,tmp_org_sigma2,tmp_org_sigma3,tmp_org_sigma4,tmp_org_sigma5,tmp_org_sigma6;
791+
long unsigned int tmp_off=0;
792+
long unsigned int tmp_off_2=0;
793+
long unsigned int list1_off=0;
794+
int num1;
795+
double complex tmp_V;
796+
double complex dam_pr;
797+
long int i_max;
798+
int tmp_Sz;
799+
long unsigned int tmp_org=0;
800+
double complex *vec_pr;
801+
vec[0]=0;
802+
i_max=X->Check.idim_max;
803+
X->Large.mode=M_CORR;
804+
805+
vec_pr = cd_1d_allocate(i_max + 1);
806+
for(i=0;i<X->Def.NTBody;i++){
807+
tmp_org_isite1 = X->Def.TBody[i][0]+1;
808+
tmp_org_sigma1 = X->Def.TBody[i][1];
809+
tmp_org_isite2 = X->Def.TBody[i][2]+1;
810+
tmp_org_sigma2 = X->Def.TBody[i][3];
811+
tmp_org_isite3 = X->Def.TBody[i][4]+1;
812+
tmp_org_sigma3 = X->Def.TBody[i][5];
813+
tmp_org_isite4 = X->Def.TBody[i][6]+1;
814+
tmp_org_sigma4 = X->Def.TBody[i][7];
815+
816+
/*[s]For three body*/
817+
tmp_org_isite5 = X->Def.TBody[i][8]+1;
818+
tmp_org_sigma5 = X->Def.TBody[i][9];
819+
tmp_org_isite6 = X->Def.TBody[i][10]+1;
820+
tmp_org_sigma6 = X->Def.TBody[i][11];
821+
/**/
822+
org_isite5 = X->Def.TBody[i][8]+1;
823+
org_sigma5 = X->Def.TBody[i][9];
824+
org_isite6 = X->Def.TBody[i][10]+1;
825+
org_sigma6 = X->Def.TBody[i][11];
826+
dam_pr = 0.0;
827+
828+
/*[s]initialized vec_pr*/
829+
#pragma omp parallel for default(none) private(j) firstprivate(i_max) shared(vec_pr)
830+
for(j=1;j<=i_max;j++){
831+
vec_pr[j] = 0.0+0.0*I;
832+
}
833+
/*[e]initialized vec_pr*/
834+
X->Large.mode = M_MLTPLY;
835+
/* |vec_pr>= c5a6|phi>*/
836+
mltplyGeneralSpinGC_mini(X,tmp_org_isite5-1,tmp_org_sigma5,tmp_org_isite6-1,tmp_org_sigma6,vec_pr,vec);
837+
X->Large.mode = M_CORR;
838+
/*[e]For three body*/
839+
840+
if(Rearray_Interactions(i, &org_isite1, &org_isite2, &org_isite3, &org_isite4, &org_sigma1, &org_sigma2, &org_sigma3, &org_sigma4, &tmp_V, X,3)!=0){
841+
fprintf(*_fp," %4ld %4ld %4ld %4ld %4ld %4ld %4ld %4ld %.10lf %.10lf \n",tmp_org_isite1-1, tmp_org_sigma1, tmp_org_isite2-1, tmp_org_sigma2, tmp_org_isite3-1,tmp_org_sigma3, tmp_org_isite4-1, tmp_org_sigma4,0.0,0.0);
842+
continue;
843+
}
844+
if(org_isite1 > X->Def.Nsite && org_isite3 > X->Def.Nsite){
845+
if(org_sigma1==org_sigma2 && org_sigma3==org_sigma4 ){ //diagonal
846+
dam_pr=child_GC_CisAisCjuAju_GeneralSpin_MPIdouble(org_isite1-1, org_sigma1, org_isite3-1, org_sigma3, tmp_V, X, vec, vec_pr);
847+
}
848+
else if(org_sigma1 == org_sigma2 && org_sigma3 != org_sigma4){
849+
dam_pr=child_GC_CisAisCjuAjv_GeneralSpin_MPIdouble(org_isite1-1, org_sigma1, org_isite3-1, org_sigma3, org_sigma4, tmp_V, X, vec, vec_pr);
850+
}
851+
else if(org_sigma1 != org_sigma2 && org_sigma3 == org_sigma4){
852+
dam_pr=child_GC_CisAitCjuAju_GeneralSpin_MPIdouble(org_isite1-1, org_sigma1, org_sigma2, org_isite3-1, org_sigma3, tmp_V, X, vec, vec_pr);
853+
}
854+
else if(org_sigma1 != org_sigma2 && org_sigma3 != org_sigma4){
855+
dam_pr=child_GC_CisAitCjuAjv_GeneralSpin_MPIdouble(org_isite1-1, org_sigma1, org_sigma2, org_isite3-1, org_sigma3, org_sigma4,tmp_V, X, vec, vec_pr);
856+
}
857+
}
858+
else if(org_isite3 > X->Def.Nsite || org_isite1 > X->Def.Nsite){
859+
if(org_sigma1==org_sigma2 && org_sigma3==org_sigma4 ){ //diagonal
860+
dam_pr=child_GC_CisAisCjuAju_GeneralSpin_MPIsingle(org_isite1-1, org_sigma1, org_isite3-1, org_sigma3, tmp_V, X, vec, vec_pr);
861+
}
862+
else if(org_sigma1 == org_sigma2 && org_sigma3 != org_sigma4){
863+
dam_pr=child_GC_CisAisCjuAjv_GeneralSpin_MPIsingle(org_isite1-1, org_sigma1, org_isite3-1, org_sigma3, org_sigma4, tmp_V, X, vec, vec_pr);
864+
}
865+
else if(org_sigma1 != org_sigma2 && org_sigma3 == org_sigma4){
866+
dam_pr=child_GC_CisAitCjuAju_GeneralSpin_MPIsingle(org_isite1-1, org_sigma1, org_sigma2, org_isite3-1, org_sigma3, tmp_V, X, vec, vec_pr);
867+
}
868+
else if(org_sigma1 != org_sigma2 && org_sigma3 != org_sigma4){
869+
dam_pr=child_GC_CisAitCjuAjv_GeneralSpin_MPIsingle(org_isite1-1, org_sigma1, org_sigma2, org_isite3-1, org_sigma3, org_sigma4,tmp_V, X, vec, vec_pr);
870+
}
871+
}
872+
else{
873+
if(org_sigma1==org_sigma2 && org_sigma3==org_sigma4 ){ //diagonal
874+
#pragma omp parallel for default(none) reduction(+:dam_pr) private(j, num1) firstprivate(i_max,X,org_isite1, org_sigma1,org_isite3, org_sigma3, tmp_V) shared(vec,vec_pr)
875+
for(j=1;j<=i_max;j++){
876+
num1=BitCheckGeneral(j-1, org_isite1, org_sigma1, X->Def.SiteToBit, X->Def.Tpow);
877+
if(num1 != FALSE){
878+
num1=BitCheckGeneral(j-1, org_isite3, org_sigma3, X->Def.SiteToBit, X->Def.Tpow);
879+
if(num1 != FALSE){
880+
dam_pr += tmp_V*conj(vec[j])*vec_pr[j];
881+
}
882+
}
883+
}
884+
}else if(org_sigma1 == org_sigma2 && org_sigma3 != org_sigma4){
885+
#pragma omp parallel for default(none) reduction(+:dam_pr) private(j, num1) firstprivate(i_max,X, org_isite1, org_isite3, org_sigma1,org_sigma3,org_sigma4, tmp_off, tmp_V) shared(vec,vec_pr)
886+
for(j=1;j<=i_max;j++){
887+
num1 = GetOffCompGeneralSpin(j-1, org_isite3, org_sigma4, org_sigma3, &tmp_off, X->Def.SiteToBit, X->Def.Tpow);
888+
if(num1 != FALSE){
889+
num1=BitCheckGeneral(tmp_off, org_isite1, org_sigma1, X->Def.SiteToBit, X->Def.Tpow);
890+
if(num1 != FALSE){
891+
dam_pr += tmp_V*conj(vec[tmp_off+1])*vec_pr[j];
892+
}
893+
}
894+
}
895+
}else if(org_sigma1 != org_sigma2 && org_sigma3 == org_sigma4){
896+
#pragma omp parallel for default(none) reduction(+:dam_pr) private(j, num1) firstprivate(i_max,X, org_isite1, org_isite3, org_sigma1,org_sigma2, org_sigma3, tmp_off, tmp_V) shared(vec,vec_pr)
897+
for(j=1;j<=i_max;j++){
898+
num1 = BitCheckGeneral(j-1, org_isite3, org_sigma3, X->Def.SiteToBit, X->Def.Tpow);
899+
if(num1 != FALSE){
900+
num1 = GetOffCompGeneralSpin(j-1, org_isite1, org_sigma2, org_sigma1, &tmp_off, X->Def.SiteToBit, X->Def.Tpow);
901+
if(num1 != FALSE){
902+
dam_pr += tmp_V*conj(vec[tmp_off+1])*vec_pr[j];
903+
}
904+
}
905+
}
906+
}else if(org_sigma1 != org_sigma2 && org_sigma3 != org_sigma4){
907+
#pragma omp parallel for default(none) reduction(+:dam_pr) private(j, num1) firstprivate(i_max,X, org_isite1, org_isite3, org_sigma1, org_sigma2, org_sigma3, org_sigma4, tmp_off, tmp_off_2, tmp_V) shared(vec,vec_pr)
908+
for(j=1;j<=i_max;j++){
909+
num1 = GetOffCompGeneralSpin(j-1, org_isite3, org_sigma4, org_sigma3, &tmp_off, X->Def.SiteToBit, X->Def.Tpow);
910+
if(num1 != FALSE){
911+
num1 = GetOffCompGeneralSpin(tmp_off, org_isite1, org_sigma2, org_sigma1, &tmp_off_2, X->Def.SiteToBit, X->Def.Tpow);
912+
if(num1 != FALSE){
913+
dam_pr += tmp_V*conj(vec[tmp_off_2+1])*vec_pr[j];
914+
}
915+
}
916+
917+
}
918+
}
919+
}
920+
dam_pr = SumMPI_dc(dam_pr);
921+
fprintf(*_fp, " %4ld %4ld %4ld %4ld "
922+
" %4ld %4ld %4ld %4ld "
923+
" %4ld %4ld %4ld %4ld "
924+
" %.10lf %.10lf \n",
925+
tmp_org_isite1-1, tmp_org_sigma1, tmp_org_isite2-1, tmp_org_sigma2,
926+
tmp_org_isite3-1, tmp_org_sigma3, tmp_org_isite4-1, tmp_org_sigma4,
927+
tmp_org_isite5-1, tmp_org_sigma5, tmp_org_isite6-1, tmp_org_sigma6,
928+
creal(dam_pr),cimag(dam_pr));
929+
}
930+
free_cd_1d_allocate(vec_pr);
931+
return 0;
932+
}
933+
934+
935+
936+
937+
774938
/**
775939
* @brief Child function to calculate two-body green's functions for General Spin model
776940
*
@@ -1104,6 +1268,9 @@ int expec_cisajscktalt_SpinGC(struct BindStruct *X,double complex *vec, FILE **_
11041268
}
11051269
} else {
11061270
info=expec_cisajscktalt_SpinGCGeneral(X,vec, _fp);
1271+
if(X->Def.NTBody>0){
1272+
info = expec_Threebody_SpinGeneral(X,vec,_fp_2);
1273+
}
11071274
}
11081275
return info;
11091276
}
@@ -1546,8 +1713,7 @@ int expec_cisajscktalt_SpinGCGeneral(struct BindStruct *X,double complex *vec, F
15461713
i_max=X->Check.idim_max;
15471714
X->Large.mode=M_CORR;
15481715

1549-
1550-
for(i=0;i<X->Def.NCisAjtCkuAlvDC;i++){
1716+
for(i=0;i<X->Def.NCisAjtCkuAlvDC;i++){
15511717
tmp_org_isite1 = X->Def.CisAjtCkuAlvDC[i][0]+1;
15521718
tmp_org_sigma1 = X->Def.CisAjtCkuAlvDC[i][1];
15531719
tmp_org_isite2 = X->Def.CisAjtCkuAlvDC[i][2]+1;
@@ -1558,13 +1724,13 @@ int expec_cisajscktalt_SpinGCGeneral(struct BindStruct *X,double complex *vec, F
15581724
tmp_org_sigma4 = X->Def.CisAjtCkuAlvDC[i][7];
15591725
dam_pr = 0.0;
15601726

1561-
if(Rearray_Interactions(i, &org_isite1, &org_isite2, &org_isite3, &org_isite4, &org_sigma1, &org_sigma2, &org_sigma3, &org_sigma4, &tmp_V, X,2)!=0){
1562-
//error message will be added
1563-
fprintf(*_fp," %4ld %4ld %4ld %4ld %4ld %4ld %4ld %4ld %.10lf %.10lf \n",tmp_org_isite1-1, tmp_org_sigma1, tmp_org_isite2-1, tmp_org_sigma2, tmp_org_isite3-1,tmp_org_sigma3, tmp_org_isite4-1, tmp_org_sigma4,0.0,0.0);
1564-
continue;
1565-
}
1727+
if(Rearray_Interactions(i, &org_isite1, &org_isite2, &org_isite3, &org_isite4, &org_sigma1, &org_sigma2, &org_sigma3, &org_sigma4, &tmp_V, X,2)!=0){
1728+
//error message will be added
1729+
fprintf(*_fp," %4ld %4ld %4ld %4ld %4ld %4ld %4ld %4ld %.10lf %.10lf \n",tmp_org_isite1-1, tmp_org_sigma1, tmp_org_isite2-1, tmp_org_sigma2, tmp_org_isite3-1,tmp_org_sigma3, tmp_org_isite4-1, tmp_org_sigma4,0.0,0.0);
1730+
continue;
1731+
}
15661732

1567-
if(org_isite1 > X->Def.Nsite && org_isite3 > X->Def.Nsite){
1733+
if(org_isite1 > X->Def.Nsite && org_isite3 > X->Def.Nsite){
15681734
if(org_sigma1==org_sigma2 && org_sigma3==org_sigma4 ){ //diagonal
15691735
dam_pr=child_GC_CisAisCjuAju_GeneralSpin_MPIdouble(org_isite1-1, org_sigma1, org_isite3-1, org_sigma3, tmp_V, X, vec, vec);
15701736
}
@@ -1594,7 +1760,7 @@ int expec_cisajscktalt_SpinGCGeneral(struct BindStruct *X,double complex *vec, F
15941760
}
15951761
else{
15961762
if(org_sigma1==org_sigma2 && org_sigma3==org_sigma4 ){ //diagonal
1597-
#pragma omp parallel for default(none) reduction(+:dam_pr) private(j, num1) firstprivate(i_max,X,org_isite1, org_sigma1,org_isite3, org_sigma3, tmp_V) shared(vec)
1763+
#pragma omp parallel for default(none) reduction(+:dam_pr) private(j, num1) firstprivate(i_max,X,org_isite1, org_sigma1,org_isite3, org_sigma3, tmp_V) shared(vec)
15981764
for(j=1;j<=i_max;j++){
15991765
num1=BitCheckGeneral(j-1, org_isite1, org_sigma1, X->Def.SiteToBit, X->Def.Tpow);
16001766
if(num1 != FALSE){
@@ -1605,7 +1771,7 @@ int expec_cisajscktalt_SpinGCGeneral(struct BindStruct *X,double complex *vec, F
16051771
}
16061772
}
16071773
}else if(org_sigma1 == org_sigma2 && org_sigma3 != org_sigma4){
1608-
#pragma omp parallel for default(none) reduction(+:dam_pr) private(j, num1) firstprivate(i_max,X, org_isite1, org_isite3, org_sigma1,org_sigma3,org_sigma4, tmp_off, tmp_V) shared(vec)
1774+
#pragma omp parallel for default(none) reduction(+:dam_pr) private(j, num1) firstprivate(i_max,X, org_isite1, org_isite3, org_sigma1,org_sigma3,org_sigma4, tmp_off, tmp_V) shared(vec)
16091775
for(j=1;j<=i_max;j++){
16101776
num1 = GetOffCompGeneralSpin(j-1, org_isite3, org_sigma4, org_sigma3, &tmp_off, X->Def.SiteToBit, X->Def.Tpow);
16111777
if(num1 != FALSE){
@@ -1616,7 +1782,7 @@ int expec_cisajscktalt_SpinGCGeneral(struct BindStruct *X,double complex *vec, F
16161782
}
16171783
}
16181784
}else if(org_sigma1 != org_sigma2 && org_sigma3 == org_sigma4){
1619-
#pragma omp parallel for default(none) reduction(+:dam_pr) private(j, num1) firstprivate(i_max,X, org_isite1, org_isite3, org_sigma1,org_sigma2, org_sigma3, tmp_off, tmp_V) shared(vec)
1785+
#pragma omp parallel for default(none) reduction(+:dam_pr) private(j, num1) firstprivate(i_max,X, org_isite1, org_isite3, org_sigma1,org_sigma2, org_sigma3, tmp_off, tmp_V) shared(vec)
16201786
for(j=1;j<=i_max;j++){
16211787
num1 = BitCheckGeneral(j-1, org_isite3, org_sigma3, X->Def.SiteToBit, X->Def.Tpow);
16221788
if(num1 != FALSE){
@@ -1627,7 +1793,7 @@ int expec_cisajscktalt_SpinGCGeneral(struct BindStruct *X,double complex *vec, F
16271793
}
16281794
}
16291795
}else if(org_sigma1 != org_sigma2 && org_sigma3 != org_sigma4){
1630-
#pragma omp parallel for default(none) reduction(+:dam_pr) private(j, num1) firstprivate(i_max,X, org_isite1, org_isite3, org_sigma1, org_sigma2, org_sigma3, org_sigma4, tmp_off, tmp_off_2, tmp_V) shared(vec)
1796+
#pragma omp parallel for default(none) reduction(+:dam_pr) private(j, num1) firstprivate(i_max,X, org_isite1, org_isite3, org_sigma1, org_sigma2, org_sigma3, org_sigma4, tmp_off, tmp_off_2, tmp_V) shared(vec)
16311797
for(j=1;j<=i_max;j++){
16321798
num1 = GetOffCompGeneralSpin(j-1, org_isite3, org_sigma4, org_sigma3, &tmp_off, X->Def.SiteToBit, X->Def.Tpow);
16331799
if(num1 != FALSE){

‎src/include/expec_cisajscktaltdc.h

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -69,5 +69,6 @@ void expec_cisajscktaltdc_alldiag_spin(
6969
double complex *vec
7070
);
7171

72+
int expec_Threebody_SpinGeneral(struct BindStruct *X,double complex *vec, FILE **_fp);
7273
int expec_Threebody_SpinGCHalf(struct BindStruct *X,double complex *vec, FILE **_fp);
7374
int expec_Fourbody_SpinGCHalf(struct BindStruct *X,double complex *vec, FILE **_fp);

‎src/include/mltplySpin.h

Lines changed: 10 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -35,6 +35,16 @@ void mltplyHalfSpinGC_mini(
3535
double complex *tmp_v1//!<[in] Input producted vector
3636
);
3737

38+
void mltplyGeneralSpinGC_mini(
39+
struct BindStruct *X,//!<[inout]
40+
int site_i,
41+
int spin_i,
42+
int site_j,
43+
int spin_j,
44+
double complex *tmp_v0,//!<[inout] Result vector
45+
double complex *tmp_v1//!<[in] Input producted vector
46+
);
47+
3848
int mltplySpinGC(struct BindStruct *X, double complex *tmp_v0,double complex *tmp_v1);
3949

4050
int mltplyHalfSpinGC(struct BindStruct *X, double complex *tmp_v0,double complex *tmp_v1);

0 commit comments

Comments
 (0)