3232
3333// #define COMPARE_ORIGINAL 1
3434
35+ #ifdef COMPARE_ORIGINAL
36+ #define TOLERANCE 1e-10
37+
38+ // if value is greater than 1, use relative value, else
39+ // use absolute value for comparison to avoid artificially big differences
40+ inline bool comp (Real val, Real test) {
41+ bool comp = true ;
42+ if (std::abs (test) > 1.0 ) {
43+ if (std ::abs ((val - test) / test) > TOLERANCE ) {
44+ comp = false ;
45+ }
46+ }
47+ else {
48+ if (std ::abs (val - test) > TOLERANCE ) {
49+ comp = false ;
50+ }
51+ }
52+ return comp;
53+ }
54+ #endif
55+
3556// Fortran layout for matrix
3657// #define IDX(i, j, lDim) (((j)*(lDim))+(i))
3758#define IDX3 (i, j, k, lDim, mDim ) \
@@ -50,9 +71,8 @@ __device__ inline void calculateHankelHip(deviceDoubleComplex prel, double r,
5071 const deviceDoubleComplex sqrtm1 = make_hipDoubleComplex (0.0 , 1.0 );
5172 deviceDoubleComplex z = hipCmul (prel , make_hipDoubleComplex (r, 0.0 ) );
5273 // deviceDoubleComplex z = prel * (deviceDoubleComplex) make_hipDoubleComplex(r, 0.0);
53- hfn[0 ] = (deviceDoubleComplex) make_hipDoubleComplex (0.0 , -1.0 ); // -sqrtm1;
54- hfn[1 ] = hipCdiv (make_hipDoubleComplex (-1.0 , 0.0 ) - sqrtm1, z);
55- // hfn[1] = (deviceDoubleComplex) make_hipDoubleComplex(-1.0, 0.0) - sqrtm1 / z;
74+ hfn[0 ] = make_hipDoubleComplex (0.0 , -1.0 ); // -sqrtm1;
75+ hfn[1 ] = make_hipDoubleComplex (-1.0 , 0.0 ) - hipCdiv (sqrtm1, z);
5676 for (int l = 1 ; l < lend; l++) {
5777 hfn[l + 1 ] = hipCdiv ((2.0 * l + 1.0 ) * hfn[l], z) - hfn[l - 1 ];
5878 // hfn[l + 1] = ((2.0 * l + 1.0) * hfn[l] / z) - hfn[l - 1];
@@ -639,11 +659,9 @@ void buildKKRMatrixLMaxIdenticalHip(LSMSSystemParameters &lsms,
639659 if (atom.LIZPos (0 , i) != testLIZPos (0 , i) ||
640660 atom.LIZPos (1 , i) != testLIZPos (1 , i) ||
641661 atom.LIZPos (2 , i) != testLIZPos (2 , i)) {
642- printf (
643- " atom.LIZPos(*,%d) [%lf,%lf,%lf] != devAtom.LIZPos(*,%d) "
644- " [%lf,%lf,%lf]\n " ,
645- i, atom.LIZPos (0 , i), atom.LIZPos (1 , i), atom.LIZPos (2 , i), i,
646- testLIZPos (0 , i), testLIZPos (1 , i), testLIZPos (2 , i));
662+ printf (" atom.LIZPos(*,%d) [%.15f,%.15f,%.15f] != devAtom.LIZPos(*,%d) [%.15f,%.15f,%.15f]\n " ,
663+ i,atom.LIZPos (0 ,i),atom.LIZPos (1 ,i),atom.LIZPos (2 ,i),
664+ i,testLIZPos (0 ,i),testLIZPos (1 ,i),testLIZPos (2 ,i));
647665 }
648666 }
649667 // loop over the LIZ blocks
@@ -688,73 +706,51 @@ void buildKKRMatrixLMaxIdenticalHip(LSMSSystemParameters &lsms,
688706 &IFactors::illp (0 , 0 ), &hfn[0 ], &dlm[0 ], &gijTest (0 , 0 ), &pi4,
689707 &lsms.global .iprint , lsms.global .istop , 32 );
690708
691- if (ir1 == 0 && ir2 == 10 ) {
709+ if (ir1 == 0 && ir2 == 1 ) {
692710 for (int l = 0 ; l <= atom.LIZlmax [ir1] + atom.LIZlmax [ir2]; l++) {
693- if (sinmp[l] != testSinmp[l])
694- printf (" sinmp[%d] (%g) != testSinmp[%d] (%g)\n " , l, sinmp[l], l,
695- testSinmp[l]);
696- if (cosmp[l] != testCosmp[l])
697- printf (" cosmp[%d] (%g) != testCosmp[%d] (%g)\n " , l, cosmp[l], l,
698- testCosmp[l]);
699- if (hfn[l] != testHfn[l])
700- printf (" hfn[%d] (%g + %gi) != testHfn[%d] (%g + %gi)\n " , l,
701- hfn[l].real (), hfn[l].imag (), l, testHfn[l].real (),
702- testHfn[l].imag ());
711+ if (!comp (sinmp[l], testSinmp[l]))
712+ printf (" sinmp[%d] (%.15f) != testSinmp[%d] (%.15f)\n " , l, sinmp[l], l, testSinmp[l]);
713+ if (!comp (cosmp[l], testCosmp[l]))
714+ printf (" cosmp[%d] (%.15f) != testCosmp[%d] (%.15f)\n " , l, cosmp[l], l, testCosmp[l]);
715+ if (!comp (hfn[l].real (), testHfn[l].real ()) ||
716+ !comp (hfn[l].imag (), testHfn[l].imag ()))
717+ printf (" hfn[%d] (%.15f + %.15fi) != testHfn[%d] (%.15f + %.15fi)\n " , l, hfn[l].real (), hfn[l].imag (), l, testHfn[l].real (), testHfn[l].imag ());
703718 }
704719 }
705720
706721 int idx = 0 ;
707722 for (int i = 0 ; i < kkri; i++)
708723 for (int j = 0 ; j < kkrj; j++) {
709- if (bgij (iOffset + i, jOffset + j) != gijTest (i, j))
710- // if (bgij[idx] != gijTest[idx] )
724+ if (! comp ( bgij (iOffset + i, jOffset + j). real (), gijTest (i, j). real ()) ||
725+ ! comp (bgij (iOffset + i, jOffset + j). imag (), gijTest (i, j). imag ()) )
711726 {
712- printf (
713- " buildBGijCPU [idx=%d]: bgij(%d + %d, %d + %d) [%g + %gi] != "
714- " gijTest(%d, %d) [%g + %gi]\n " ,
715- idx, iOffset, i, jOffset, j,
716- bgij (iOffset + i, jOffset + j).real (),
717- bgij (iOffset + i, jOffset + j).imag (), i, j,
718- gijTest (i, j).real (), gijTest (i, j).imag ());
727+ printf (" buildBGijCPU [idx=%d]: bgij(%d + %d, %d + %d) [%.15f + %.15fi] != gijTest(%d, %d) [%.15f + %.15fi]\n " , idx,
728+ iOffset, i, jOffset, j, bgij (iOffset + i, jOffset + j).real (), bgij (iOffset + i, jOffset + j).imag (),
729+ i, j, gijTest (i, j).real (), gijTest (i, j).imag ());
719730 exitCompare = true ;
720731 }
721- if (bgij (iOffset + kkri + i, jOffset + kkrj + j) != gijTest (i, j))
722- // if (bgij[idx] != gijTest[idx] )
732+ if (! comp ( bgij (iOffset + kkri + i, jOffset + kkrj + j). real (), gijTest (i, j). real ()) ||
733+ ! comp (bgij (iOffset + kkri + i, jOffset + kkrj + j). imag (), gijTest (i, j). imag ()) )
723734 {
724- printf (
725- " buildBGijCPU : bgij(%d + %d, %d + %d) [%g + %gi] != "
726- " gijTest(%d, %d) [%g + %gi]\n " ,
727- iOffset, i + kkri, jOffset, j + kkrj,
728- bgij (iOffset + kkri + i, jOffset + kkrj + j).real (),
729- bgij (iOffset + kkri + i, jOffset + kkrj + j).imag (), i, j,
730- gijTest (i, j).real (), gijTest (i, j).imag ());
735+ printf (" buildBGijCPU : bgij(%d + %d, %d + %d) [%.15f + %.15fi] != gijTest(%d, %d) [%.15f + %.15fi]\n " ,
736+ iOffset, i + kkri, jOffset, j + kkrj, bgij (iOffset + kkri + i, jOffset + kkrj + j).real (), bgij (iOffset + kkri + i, jOffset + kkrj + j).imag (),
737+ i, j, gijTest (i, j).real (), gijTest (i, j).imag ());
731738 exitCompare = true ;
732739 }
733- if (bgij (iOffset + kkri + i, jOffset + j) !=
734- 0.0 ) // gijTest(i+kkri,j))
735- // if(bgij[idx] != gijTest[idx])
740+ if (bgij (iOffset + kkri + i, jOffset + j) != 0.0 )
736741 {
737- printf (
738- " buildBGijCPU : bgij(%d + %d, %d + %d) [%g + %gi] != 0.0\n " ,
739- iOffset, i + kkri, jOffset, j,
740- bgij (iOffset + kkri + i, jOffset + j).real (),
741- bgij (iOffset + kkri + i, jOffset + j).imag ());
742+ printf (" buildBGijCPU : bgij(%d + %d, %d + %d) [%.15f + %.15fi] != 0.0\n " ,
743+ iOffset, i + kkri, jOffset, j, bgij (iOffset + kkri + i, jOffset + j).real (), bgij (iOffset + kkri + i, jOffset + j).imag ());
742744 exitCompare = true ;
743745 }
744- if (bgij (iOffset + i, jOffset + kkrj + j) !=
745- 0.0 ) // gijTest(i,j+kkrj))
746- // if(bgij[idx] != gijTest[idx])
746+ if (bgij (iOffset + i, jOffset + kkrj + j) != 0.0 )
747747 {
748- printf (
749- " buildBGijCPU : bgij(%d + %d, %d + %d) [%g + %gi] != 0.0\n " ,
750- iOffset, i, jOffset, j + kkrj,
751- bgij (iOffset + i, jOffset + kkrj + j).real (),
752- bgij (iOffset + i, jOffset + kkrj + j).imag ());
748+ printf (" buildBGijCPU : bgij(%d + %d, %d + %d) [%.15f + %.15fi] != 0.0\n " ,
749+ iOffset, i, jOffset, j + kkrj, bgij (iOffset + i, jOffset + kkrj + j).real (), bgij (iOffset + i, jOffset + kkrj + j).imag ());
753750 exitCompare = true ;
754751 }
755752 idx++;
756753 }
757- // if(exitCompare) exit(1);
758754 }
759755 }
760756 if (exitCompare) exit (1 );
@@ -866,22 +862,22 @@ local.tmatStore(iie*local.blkSizeTmatStore + , atom.LIZStoreIdx[ir1]) *
866862 buildKKRMatrixCPU (lsms, local, atom, ispin, iie, energy, prel, mCPU );
867863
868864 for (int i = 0 ; i < nrmat_ns; i++)
869- for (int j = 0 ; j < nrmat_ns; j++) {
870- if (mCPU (i, j) != mGPU (i, j))
871- // if(bgij[idx] != gijTest[idx])
865+ for (int j = 0 ; j < nrmat_ns; j++)
866+ {
867+ if (!comp (mCPU (i, j).real (), mGPU (i, j).real ()) ||
868+ !comp (mCPU (i, j).imag (), mGPU (i, j).imag ()))
872869 {
873- printf (
874- " buildBGijCPU : mCPU(%d, %d) [%g + %gi] != mGPU(%d, %d) [%g + "
875- " %gi]\n " ,
876- i, j, mCPU (i, j).real (), mCPU (i, j).imag (), i, j, mGPU (i, j).real (),
877- mGPU (i, j).imag ());
870+ printf (" buildBGijCPU : mCPU(%d, %d) [%.15f + %.15fi] != mGPU(%d, %d) [%.15f + %.15fi]\n " ,
871+ i, j, mCPU (i, j).real (), mCPU (i, j).imag (),
872+ i, j, mGPU (i, j).real (), mGPU (i, j).imag ());
878873 exitCompare = true ;
879874 }
880875 }
881876
882877 if (exitCompare) exit (1 );
883878#endif
884879}
880+
885881void buildKKRMatrixLMaxDifferentHip (LSMSSystemParameters &lsms,
886882 LocalTypeInfo &local, AtomData &atom,
887883 DeviceStorage &d, DeviceAtom &devAtom,
0 commit comments