Merge branch 'develop' into tensor-lapack-extend

Author: Cstandardlib

source/module_basis/module_pw/pw_basis.cpp

@@ -38,7 +38,8 @@ PW_Basis:: ~PW_Basis()
     delete[] ig2igg;
     delete[] gg_uniq;
 #if defined(__CUDA) || defined(__ROCM)
-    if (this->device == "gpu") {
+    if (this->device == "gpu")
+    {
         delmem_int_op()(this->d_is2fftixy);
     }
 #endif
@@ -70,9 +71,14 @@ void PW_Basis::setuptransform()
 
 void PW_Basis::getstartgr()
 {
-    if(this->gamma_only)    this->nmaxgr = ( this->npw > (this->nrxx+1)/2 ) ? this->npw : (this->nrxx+1)/2;
-    else                    this->nmaxgr = ( this->npw > this->nrxx ) ? this->npw : this->nrxx;
-    this->nmaxgr = (this->nz * this->nst > this->nxy * nplane) ? this->nz * this->nst : this->nxy * nplane;
+    if(this->gamma_only)    
+    {
+        this->nmaxgr = ( this->npw > (this->nrxx+1)/2 ) ? this->npw : (this->nrxx+1)/2;
+    }
+    else
+    {
+        this->nmaxgr = ( this->npw > this->nrxx ) ? this->npw : this->nrxx;
+    }
 
     //---------------------------------------------
 	// sum : starting plane of FFT box.
@@ -85,23 +91,35 @@ void PW_Basis::getstartgr()
 	// Each processor has a set of full sticks,
 	// 'rank_use' processor send a piece(npps[ip]) of these sticks(nst_per[rank_use])
 	// to all the other processors in this pool
-	for (int ip = 0;ip < poolnproc; ++ip) this->numg[ip] = this->nst_per[poolrank] * this->numz[ip];
+	for (int ip = 0;ip < poolnproc; ++ip)
+    {
+        this->numg[ip] = this->nst_per[poolrank] * this->numz[ip];
+    }
 
 
 	// Each processor in a pool send a piece of each stick(nst_per[ip]) to
 	// other processors in this pool
 	// rank_use processor receive datas in npps[rank_p] planes.
-	for (int ip = 0;ip < poolnproc; ++ip) this->numr[ip] = this->nst_per[ip] * this->numz[poolrank];
+	for (int ip = 0;ip < poolnproc; ++ip)
+    {
+        this->numr[ip] = this->nst_per[ip] * this->numz[poolrank];
+    }
 
 
 	// startg record the starting 'numg' position in each processor.
 	this->startg[0] = 0;
-	for (int ip = 1;ip < poolnproc; ++ip) this->startg[ip] = this->startg[ip-1] + this->numg[ip-1];
+	for (int ip = 1;ip < poolnproc; ++ip)
+    {
+        this->startg[ip] = this->startg[ip-1] + this->numg[ip-1];
+    }
 
 
 	// startr record the starting 'numr' position
 	this->startr[0] = 0;
-	for (int ip = 1;ip < poolnproc; ++ip) this->startr[ip] = this->startr[ip-1] + this->numr[ip-1];
+	for (int ip = 1;ip < poolnproc; ++ip)
+    {
+        this->startr[ip] = this->startr[ip-1] + this->numr[ip-1];
+    }
     return;
 }
 
@@ -112,7 +130,10 @@ void PW_Basis::getstartgr()
 /// 
 void PW_Basis::collect_local_pw()
 {
-    if(this->npw <= 0) return;
+    if(this->npw <= 0)
+    {
+        return;
+    }
     this->ig_gge0 = -1;
     delete[] this->gg; this->gg = new double[this->npw];
     delete[] this->gdirect; this->gdirect = new ModuleBase::Vector3<double>[this->npw];
@@ -127,16 +148,28 @@ void PW_Basis::collect_local_pw()
         int ixy = this->is2fftixy[is];
         int ix = ixy / this->fftny;
         int iy = ixy % this->fftny;
-        if (ix >= int(this->nx/2) + 1) ix -= this->nx;
-        if (iy >= int(this->ny/2) + 1) iy -= this->ny;
-        if (iz >= int(this->nz/2) + 1) iz -= this->nz;
+        if (ix >= int(this->nx/2) + 1)
+        {
+            ix -= this->nx;
+        }
+        if (iy >= int(this->ny/2) + 1)
+        {
+            iy -= this->ny;
+        }
+        if (iz >= int(this->nz/2) + 1)
+        {
+            iz -= this->nz;
+        }
         f.x = ix;
         f.y = iy;
         f.z = iz;
         this->gg[ig] = f * (this->GGT * f);
         this->gdirect[ig] = f;
         this->gcar[ig] = f * this->G;
-        if(this->gg[ig] < 1e-8) this->ig_gge0 = ig;
+        if(this->gg[ig] < 1e-8)
+        {
+            this->ig_gge0 = ig;
+        }
     }
     return;
 }
@@ -148,10 +181,13 @@ void PW_Basis::collect_local_pw()
 /// 
 void PW_Basis::collect_uniqgg()
 {
-    if(this->npw <= 0) return;
+    if(this->npw <= 0)
+    {
+        return;
+    }
     this->ig_gge0 = -1;
     delete[] this->ig2igg; this->ig2igg = new int [this->npw];
-//add by A.s 202406
+    
     int *sortindex = new int [this->npw];//Reconstruct the mapping of the plane wave index ig according to the energy size of the plane waves
     double *tmpgg = new double [this->npw];//Ranking the plane waves by energy size while ensuring that the same energy is preserved for each wave to correspond
     double *tmpgg2 = new double [this->npw];//ranking the plane waves by energy size and removing the duplicates
@@ -164,14 +200,26 @@ void PW_Basis::collect_uniqgg()
         int ixy = this->is2fftixy[is];
         int ix = ixy / this->fftny;
         int iy = ixy % this->fftny;
-        if (ix >= int(this->nx/2) + 1) ix -= this->nx;
-        if (iy >= int(this->ny/2) + 1) iy -= this->ny;
-        if (iz >= int(this->nz/2) + 1) iz -= this->nz;
+        if (ix >= int(this->nx/2) + 1)
+        {
+            ix -= this->nx;
+        }
+        if (iy >= int(this->ny/2) + 1)
+        {
+            iy -= this->ny;
+        }
+        if (iz >= int(this->nz/2) + 1)
+        {
+            iz -= this->nz;
+        }
         f.x = ix;
         f.y = iy;
         f.z = iz;
         tmpgg[ig] = f * (this->GGT * f);
-        if(tmpgg[ig] < 1e-8) this->ig_gge0 = ig;
+        if(tmpgg[ig] < 1e-8)
+        {
+            this->ig_gge0 = ig;
+        }
     }
 
     ModuleBase::GlobalFunc::ZEROS(sortindex, this->npw);
@@ -215,7 +263,10 @@ void PW_Basis::collect_uniqgg()
 void PW_Basis::getfftixy2is(int * fftixy2is) const
 {
 //Note: please assert when is1 >= is2, fftixy2is[is1] >= fftixy2is[is2]!
-    for(int ixy = 0 ; ixy < this->fftnxy ; ++ixy)   fftixy2is[ixy] = -1;
+    for(int ixy = 0 ; ixy < this->fftnxy ; ++ixy)
+    {
+        fftixy2is[ixy] = -1;
+    }
     int ixy = 0;
     for(int is = 0; is < this->nst; ++is)
     {

source/module_basis/module_pw/pw_basis.h

@@ -237,11 +237,11 @@ class PW_Basis
     int nx=0, ny=0, nz=0, nxyz=0, nxy=0; // Gamma_only: fftny = int(ny/2)-1 , others: fftny = ny
     int liy=0, riy=0;// liy: the left edge of the pw ball; riy: the right edge of the pw ball in the y direction
     int lix=0, rix=0;// lix: the left edge of the pw ball; rix: the right edge of the pw ball in the x direction
-    bool xprime = true; // true: when do recip2real, x-fft will be done last and when doing real2recip, x-fft will be done first; false: y-fft
-                         // For gamma_only, true: we use half x; false: we use half y
+    bool xprime = true; // true: when do recip2real, x-fft will be done last and when doing real2recip, x-fft will be
+                        // done first; false: y-fft For gamma_only, true: we use half x; false: we use half y
     int ng_xeq0 = 0; //only used when xprime = true, number of g whose gx = 0
-    int nmaxgr=0; // Gamma_only: max between npw and (nrxx+1)/2, others: max between npw and nrxx
-                // Thus complex<double>[nmaxgr] is able to contain either reciprocal or real data
+    int nmaxgr = 0;  // Gamma_only: max between npw and (nrxx+1)/2, others: max between npw and nrxx
+                     // Thus complex<double>[nmaxgr] is able to contain either reciprocal or real data
     // FFT ft;
     FFT_Bundle fft_bundle;
     //The position of pointer in and out can be equal(in-place transform) or different(out-of-place transform).

source/module_basis/module_pw/test_serial/pw_basis_test.cpp

@@ -278,7 +278,7 @@ TEST_F(PWBasisTEST,GetStartGR)
 	EXPECT_EQ(pwb.nrxx,8000);
 	EXPECT_EQ(pwb.nxy,400);
 	EXPECT_EQ(pwb.nplane,20);
-	EXPECT_EQ(pwb.nmaxgr,8000);
+	EXPECT_EQ(pwb.nmaxgr,4000);
 	EXPECT_EQ(pwb.numg[0],3120);
 	EXPECT_EQ(pwb.numr[0],3120);
 	EXPECT_EQ(pwb.startg[0],0);

source/module_io/unk_overlap_pw.cpp

@@ -75,7 +75,7 @@ std::complex<double> unkOverlap_pw::unkdotp_G0(const ModulePW::PW_Basis* rhopw,
     std::complex<double>* phase = new std::complex<double>[rhopw->nmaxgr];
 
     // get the phase value in realspace
-    for (int ig = 0; ig < rhopw->nmaxgr; ig++)
+    for (int ig = 0; ig < rhopw->npw; ig++)
     {
         ModuleBase::Vector3<double> delta_G = rhopw->gdirect[ig] - G;
         if (delta_G.norm2() < 1e-10) // rhopw->gdirect[ig] == G
@@ -89,7 +89,7 @@ std::complex<double> unkOverlap_pw::unkdotp_G0(const ModulePW::PW_Basis* rhopw,
     rhopw->recip2real(phase, phase);
     wfcpw->recip2real(&evc[0](ik_L, iband_L, 0), psi_r, ik_L);
 
-    for (int ir = 0; ir < rhopw->nmaxgr; ir++)
+    for (int ir = 0; ir < rhopw->nrxx; ir++)
     {
         psi_r[ir] = psi_r[ir] * phase[ir];
     }

tests/integrate/tools/catch_properties.sh

@@ -401,7 +401,7 @@ if ! test -z "$out_dm"  && [ $out_dm == 1 ]; then
 fi
 
 if ! test -z "$out_mul"  && [ $out_mul == 1 ]; then
-    python3 ../tools/CompareFile.py mulliken.txt.ref OUT.autotest/mulliken.txt 6
+    python3 ../tools/CompareFile.py mulliken.txt.ref OUT.autotest/mulliken.txt 4
 	echo "Compare_mulliken_pass $?" >>$1
 fi