diff --git a/source/module_cell/atom_spec.cpp b/source/module_cell/atom_spec.cpp
index f00a062c1d..0bf5045f2d 100644
--- a/source/module_cell/atom_spec.cpp
+++ b/source/module_cell/atom_spec.cpp
@@ -102,51 +102,54 @@ void Atom::bcast_atom(void)
     Parallel_Common::bcast_bool(this->flag_empty_element);
     Parallel_Common::bcast_double(mass);
 
-    if (GlobalV::MY_RANK != 0)
+    if (na > 0)
     {
-        assert(na != 0);
-        this->tau.resize(na, ModuleBase::Vector3<double>(0, 0, 0));
-        this->dis.resize(na, ModuleBase::Vector3<double>(0, 0, 0));
-        this->taud.resize(na, ModuleBase::Vector3<double>(0, 0, 0));
-        this->vel.resize(na, ModuleBase::Vector3<double>(0, 0, 0));
-        this->mag.resize(na, 0);
-        this->angle1.resize(na, 0);
-        this->angle2.resize(na, 0);
-        this->m_loc_.resize(na, ModuleBase::Vector3<double>(0, 0, 0));
-        this->mbl.resize(na, ModuleBase::Vector3<int>(0, 0, 0));
-        this->lambda.resize(na, ModuleBase::Vector3<double>(0, 0, 0));
-        this->constrain.resize(na, ModuleBase::Vector3<int>(0, 0, 0));
-    }
+        if (GlobalV::MY_RANK != 0)
+        {
+            assert(na != 0);
+            this->tau.resize(na, ModuleBase::Vector3<double>(0, 0, 0));
+            this->dis.resize(na, ModuleBase::Vector3<double>(0, 0, 0));
+            this->taud.resize(na, ModuleBase::Vector3<double>(0, 0, 0));
+            this->vel.resize(na, ModuleBase::Vector3<double>(0, 0, 0));
+            this->mag.resize(na, 0);
+            this->angle1.resize(na, 0);
+            this->angle2.resize(na, 0);
+            this->m_loc_.resize(na, ModuleBase::Vector3<double>(0, 0, 0));
+            this->mbl.resize(na, ModuleBase::Vector3<int>(0, 0, 0));
+            this->lambda.resize(na, ModuleBase::Vector3<double>(0, 0, 0));
+            this->constrain.resize(na, ModuleBase::Vector3<int>(0, 0, 0));
+        }
 
-    for (int i = 0; i < na; i++)
-    {
-        Parallel_Common::bcast_double(tau[i].x);
-        Parallel_Common::bcast_double(tau[i].y);
-        Parallel_Common::bcast_double(tau[i].z);
-        Parallel_Common::bcast_double(taud[i].x);
-        Parallel_Common::bcast_double(taud[i].y);
-        Parallel_Common::bcast_double(taud[i].z);
-        Parallel_Common::bcast_double(dis[i].x);
-        Parallel_Common::bcast_double(dis[i].y);
-        Parallel_Common::bcast_double(dis[i].z);
-        Parallel_Common::bcast_double(vel[i].x);
-        Parallel_Common::bcast_double(vel[i].y);
-        Parallel_Common::bcast_double(vel[i].z);
-        Parallel_Common::bcast_double(mag[i]);
-        Parallel_Common::bcast_double(angle1[i]);
-        Parallel_Common::bcast_double(angle2[i]);
-        Parallel_Common::bcast_double(m_loc_[i].x);
-        Parallel_Common::bcast_double(m_loc_[i].y);
-        Parallel_Common::bcast_double(m_loc_[i].z);
-        Parallel_Common::bcast_int(mbl[i].x);
-        Parallel_Common::bcast_int(mbl[i].y);
-        Parallel_Common::bcast_int(mbl[i].z);
-        Parallel_Common::bcast_double(lambda[i].x);
-        Parallel_Common::bcast_double(lambda[i].y);
-        Parallel_Common::bcast_double(lambda[i].z);
-        Parallel_Common::bcast_int(constrain[i].x);
-        Parallel_Common::bcast_int(constrain[i].y);
-        Parallel_Common::bcast_int(constrain[i].z);
+        for (int i = 0; i < na; i++)
+        {
+            Parallel_Common::bcast_double(tau[i].x);
+            Parallel_Common::bcast_double(tau[i].y);
+            Parallel_Common::bcast_double(tau[i].z);
+            Parallel_Common::bcast_double(taud[i].x);
+            Parallel_Common::bcast_double(taud[i].y);
+            Parallel_Common::bcast_double(taud[i].z);
+            Parallel_Common::bcast_double(dis[i].x);
+            Parallel_Common::bcast_double(dis[i].y);
+            Parallel_Common::bcast_double(dis[i].z);
+            Parallel_Common::bcast_double(vel[i].x);
+            Parallel_Common::bcast_double(vel[i].y);
+            Parallel_Common::bcast_double(vel[i].z);
+            Parallel_Common::bcast_double(mag[i]);
+            Parallel_Common::bcast_double(angle1[i]);
+            Parallel_Common::bcast_double(angle2[i]);
+            Parallel_Common::bcast_double(m_loc_[i].x);
+            Parallel_Common::bcast_double(m_loc_[i].y);
+            Parallel_Common::bcast_double(m_loc_[i].z);
+            Parallel_Common::bcast_int(mbl[i].x);
+            Parallel_Common::bcast_int(mbl[i].y);
+            Parallel_Common::bcast_int(mbl[i].z);
+            Parallel_Common::bcast_double(lambda[i].x);
+            Parallel_Common::bcast_double(lambda[i].y);
+            Parallel_Common::bcast_double(lambda[i].z);
+            Parallel_Common::bcast_int(constrain[i].x);
+            Parallel_Common::bcast_int(constrain[i].y);
+            Parallel_Common::bcast_int(constrain[i].z);
+        }
     }
 
     return;
diff --git a/source/module_cell/read_atoms.cpp b/source/module_cell/read_atoms.cpp
index 81608ce609..30ae5c504b 100644
--- a/source/module_cell/read_atoms.cpp
+++ b/source/module_cell/read_atoms.cpp
@@ -504,7 +504,18 @@ bool UnitCell::read_atom_positions(std::ifstream &ifpos, std::ofstream &ofs_runn
                 ModuleBase::WARNING("read_atom_positions", " atom number < 0.");
                 return false;
             }
-            if (na > 0)
+            else if (na == 0)
+            {
+                std::cout << "%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%" << std::endl;
+                std::cout << " Warning: atom number is 0 for atom type: " << atoms[it].label << std::endl;
+                std::cout << " If you are confident that this is not a mistake, please ignore this warning." << std::endl;
+                std::cout << "%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%" << std::endl;
+                ofs_running << "%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%" << std::endl;
+                ofs_running << " Warning: atom number is 0 for atom type: " << atoms[it].label << std::endl;
+                ofs_running << " If you are confident that this is not a mistake, please ignore this warning." << std::endl;
+                ofs_running << "%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%" << std::endl;
+            }
+            else if (na > 0)
             {
                 atoms[it].tau.resize(na, ModuleBase::Vector3<double>(0,0,0));
                 atoms[it].dis.resize(na, ModuleBase::Vector3<double>(0,0,0));
@@ -891,6 +902,12 @@ bool UnitCell::read_atom_positions(std::ifstream &ifpos, std::ofstream &ofs_runn
     ofs_running << std::endl;
     ModuleBase::GlobalFunc::OUT(ofs_running,"TOTAL ATOM NUMBER",nat);
 
+    if (nat == 0)
+    {
+        ModuleBase::WARNING("read_atom_positions","no atom in the system!");
+        return false;
+    }
+
     // mohan add 2010-06-30    
     //xiaohui modify 2015-03-15, cancel outputfile "STRU_READIN.xyz"
     //this->print_cell_xyz("STRU_READIN.xyz");
diff --git a/source/module_hamilt_pw/hamilt_pwdft/forces.cpp b/source/module_hamilt_pw/hamilt_pwdft/forces.cpp
index c1fcd2299c..dd24ce82c7 100644
--- a/source/module_hamilt_pw/hamilt_pwdft/forces.cpp
+++ b/source/module_hamilt_pw/hamilt_pwdft/forces.cpp
@@ -595,21 +595,24 @@ void Forces<FPTYPE, Device>::cal_force_ew(const UnitCell& ucell,
         }
         for (int it = 0; it < ucell.ntype; it++)
         {
-            double dzv;
-            if (PARAM.inp.use_paw)
+            if (ucell.atoms[it].na != 0)
             {
-#ifdef USE_PAW
-                dzv = GlobalC::paw_cell.get_val(it);
-#endif
-            }
-            else
-            {
-                dzv = ucell.atoms[it].ncpp.zv;
-            }
+                double dzv;
+                if (PARAM.inp.use_paw)
+                {
+    #ifdef USE_PAW
+                    dzv = GlobalC::paw_cell.get_val(it);
+    #endif
+                }
+                else
+                {
+                    dzv = ucell.atoms[it].ncpp.zv;
+                }
 
-            for (int ig = igb; ig < ig_end; ++ig)
-            { // accumulate aux
-                aux[ig] += dzv * conj(p_sf->strucFac(it, ig));
+                for (int ig = igb; ig < ig_end; ++ig)
+                { // accumulate aux
+                    aux[ig] += dzv * conj(p_sf->strucFac(it, ig));
+                }
             }
         }
     }
@@ -714,30 +717,33 @@ void Forces<FPTYPE, Device>::cal_force_ew(const UnitCell& ucell,
                 last_it = it;
             }
 
-            const auto ig_loop = [&](int ig_beg, int ig_end) {
-                for (int ig = ig_beg; ig < ig_end; ig++)
-                {
-                    const ModuleBase::Vector3<double> gcar = rho_basis->gcar[ig];
-                    const double arg = ModuleBase::TWO_PI * (gcar * ucell.atoms[it].tau[ia]);
-                    double sinp, cosp;
-                    ModuleBase::libm::sincos(arg, &sinp, &cosp);
-                    double sumnb = -cosp * aux[ig].imag() + sinp * aux[ig].real();
-                    forceion(iat, 0) += gcar[0] * sumnb;
-                    forceion(iat, 1) += gcar[1] * sumnb;
-                    forceion(iat, 2) += gcar[2] * sumnb;
-                }
-            };
+            if (ucell.atoms[it].na != 0)
+            {
+                const auto ig_loop = [&](int ig_beg, int ig_end) {
+                    for (int ig = ig_beg; ig < ig_end; ig++)
+                    {
+                        const ModuleBase::Vector3<double> gcar = rho_basis->gcar[ig];
+                        const double arg = ModuleBase::TWO_PI * (gcar * ucell.atoms[it].tau[ia]);
+                        double sinp, cosp;
+                        ModuleBase::libm::sincos(arg, &sinp, &cosp);
+                        double sumnb = -cosp * aux[ig].imag() + sinp * aux[ig].real();
+                        forceion(iat, 0) += gcar[0] * sumnb;
+                        forceion(iat, 1) += gcar[1] * sumnb;
+                        forceion(iat, 2) += gcar[2] * sumnb;
+                    }
+                };
 
-            // skip ig_gge0 point by separating ig loop into two part
-            ig_loop(0, ig_gap);
-            ig_loop(ig_gap + 1, rho_basis->npw);
+                // skip ig_gge0 point by separating ig loop into two part
+                ig_loop(0, ig_gap);
+                ig_loop(ig_gap + 1, rho_basis->npw);
 
-            forceion(iat, 0) *= it_fact;
-            forceion(iat, 1) *= it_fact;
-            forceion(iat, 2) *= it_fact;
+                forceion(iat, 0) *= it_fact;
+                forceion(iat, 1) *= it_fact;
+                forceion(iat, 2) *= it_fact;
 
-            ++iat;
-            ucell.step_iait(&ia, &it);
+                ++iat;
+                ucell.step_iait(&ia, &it);
+            }
         }
 
         // means that the processor contains G=0 term.
@@ -771,7 +777,7 @@ void Forces<FPTYPE, Device>::cal_force_ew(const UnitCell& ucell,
                 int T2 = 0;
                 while (iat2 < this->nat)
                 {
-                    if (iat1 != iat2)
+                    if (iat1 != iat2 && ucell.atoms[T2].na != 0 && ucell.atoms[T1].na != 0)
                     {
                         ModuleBase::Vector3<double> d_tau
                             = ucell.atoms[T1].tau[I1] - ucell.atoms[T2].tau[I2];
diff --git a/source/module_hamilt_pw/hamilt_pwdft/stress_func_ewa.cpp b/source/module_hamilt_pw/hamilt_pwdft/stress_func_ewa.cpp
index fd6bf8dc0a..5dbe2cd6b4 100644
--- a/source/module_hamilt_pw/hamilt_pwdft/stress_func_ewa.cpp
+++ b/source/module_hamilt_pw/hamilt_pwdft/stress_func_ewa.cpp
@@ -137,25 +137,28 @@ void Stress_Func<FPTYPE, Device>::stress_ewa(const UnitCell& ucell,
 
 		while (ijat < ijat_end)
 		{
-			//calculate tau[na]-tau[nb]
-			d_tau = ucell.atoms[it].tau[i] - ucell.atoms[jt].tau[j];
-			//generates nearest-neighbors shells 
-			H_Ewald_pw::rgen(d_tau, rmax, irr, ucell.latvec, ucell.G, r, r2, nrm);
-			for(int nr=0 ; nr<nrm ; nr++)
+			if (ucell.atoms[it].na != 0 && ucell.atoms[jt].na != 0)
 			{
-				rr=sqrt(r2[nr]) * ucell.lat0;
-				fac = -ModuleBase::e2/2.0/ucell.omega*pow(ucell.lat0,2)*ucell.atoms[it].ncpp.zv * ucell.atoms[jt].ncpp.zv / pow(rr,3) * (erfc(sqa*rr)+rr * sq8a_2pi *  ModuleBase::libm::exp(-alpha * pow(rr,2)));
-				for(int l=0; l<3; l++)
+				//calculate tau[na]-tau[nb]
+				d_tau = ucell.atoms[it].tau[i] - ucell.atoms[jt].tau[j];
+				//generates nearest-neighbors shells 
+				H_Ewald_pw::rgen(d_tau, rmax, irr, ucell.latvec, ucell.G, r, r2, nrm);
+				for(int nr=0 ; nr<nrm ; nr++)
 				{
-					for(int m=0; m<l+1; m++)
+					rr=sqrt(r2[nr]) * ucell.lat0;
+					fac = -ModuleBase::e2/2.0/ucell.omega*pow(ucell.lat0,2)*ucell.atoms[it].ncpp.zv * ucell.atoms[jt].ncpp.zv / pow(rr,3) * (erfc(sqa*rr)+rr * sq8a_2pi *  ModuleBase::libm::exp(-alpha * pow(rr,2)));
+					for(int l=0; l<3; l++)
 					{
-						r0[0] = r[nr].x;
-						r0[1] = r[nr].y;
-						r0[2] = r[nr].z;
-						local_sigma(l,m) += fac * r0[l] * r0[m];
-					}//end m
-				}//end l
-			}//end nr
+						for(int m=0; m<l+1; m++)
+						{
+							r0[0] = r[nr].x;
+							r0[1] = r[nr].y;
+							r0[2] = r[nr].z;
+							local_sigma(l,m) += fac * r0[l] * r0[m];
+						}//end m
+					}//end l
+				}//end nr
+			}
 
 			++ijat;
 			ucell.step_jajtiait(&j, &jt, &i, &it);
diff --git a/tests/integrate/123_PW_zero_atom/INPUT b/tests/integrate/123_PW_zero_atom/INPUT
new file mode 100644
index 0000000000..8488f61585
--- /dev/null
+++ b/tests/integrate/123_PW_zero_atom/INPUT
@@ -0,0 +1,28 @@
+INPUT_PARAMETERS
+#Parameters (1.General)
+suffix			autotest
+calculation     	scf
+
+nbands			6
+symmetry		0
+pseudo_dir      	../../PP_ORB/
+
+#Parameters (2.Iteration)
+ecutwfc			20
+
+#Parameters (3.Basis)
+basis_type		pw
+
+#Parameters (4.Smearing)
+smearing_method		gauss
+smearing_sigma		0.0002
+
+#Parameters (5.Mixing)
+mixing_type		broyden
+mixing_beta		0.7
+cal_force		1
+test_force		1
+cal_stress		1
+test_stress		1
+
+pw_seed         1
\ No newline at end of file
diff --git a/tests/integrate/123_PW_zero_atom/KPT b/tests/integrate/123_PW_zero_atom/KPT
new file mode 100644
index 0000000000..f5f7f4ec34
--- /dev/null
+++ b/tests/integrate/123_PW_zero_atom/KPT
@@ -0,0 +1,4 @@
+K_POINTS
+0
+Gamma
+2 2 2 0 0 0
diff --git a/tests/integrate/123_PW_zero_atom/STRU b/tests/integrate/123_PW_zero_atom/STRU
new file mode 100644
index 0000000000..49a2ef9cb7
--- /dev/null
+++ b/tests/integrate/123_PW_zero_atom/STRU
@@ -0,0 +1,23 @@
+ATOMIC_SPECIES
+Al 13 Al_ONCV_PBE-1.0.upf upf201
+Si 14 Si_ONCV_PBE-1.0.upf upf201
+
+LATTICE_CONSTANT
+10.2  // add lattice constant
+
+LATTICE_VECTORS
+0.5 0.5 0.0
+0.5 0.0 0.5
+0.0 0.5 0.5
+
+ATOMIC_POSITIONS
+Direct 
+
+Al
+0.0
+0
+Si      // Element type
+0.0     // magnetism
+2
+0.00 0.00 0.00 1 1 1
+0.25 0.25 0.25 1 1 1
diff --git a/tests/integrate/123_PW_zero_atom/jd b/tests/integrate/123_PW_zero_atom/jd
new file mode 100644
index 0000000000..ebecbd81ec
--- /dev/null
+++ b/tests/integrate/123_PW_zero_atom/jd
@@ -0,0 +1 @@
+test calculation when atom number is zero.
diff --git a/tests/integrate/123_PW_zero_atom/result.ref b/tests/integrate/123_PW_zero_atom/result.ref
new file mode 100644
index 0000000000..c93cf82f0f
--- /dev/null
+++ b/tests/integrate/123_PW_zero_atom/result.ref
@@ -0,0 +1,5 @@
+etotref -211.8003327929410489
+etotperatomref -105.9001663965
+totalforceref 0.000014
+totalstressref 368.726447
+totaltimeref 0.61
diff --git a/tests/integrate/CASES_CPU.txt b/tests/integrate/CASES_CPU.txt
index 8277d51eab..1c76104ae0 100644
--- a/tests/integrate/CASES_CPU.txt
+++ b/tests/integrate/CASES_CPU.txt
@@ -98,6 +98,7 @@
 121_PW_KPAR
 121_PW_kspacing
 127_PW_15_PK_AF
+123_PW_zero_atom
 128_PW_zero_ntype
 133_PW_DJ_PK
 135_PW_15_PK