Fix a readin bug, add pchg test in integrate tests, update some examples

source/source_cell/read_atoms.cpp

@@ -26,7 +26,9 @@ bool unitcell::read_atom_positions(UnitCell& ucell,
     ModuleBase::TITLE("UnitCell","read_atom_positions");
 
     std::string& Coordinate  = ucell.Coordinate;
-    const int    ntype       = ucell.ntype;
+    const int ntype = ucell.ntype;
+    const int nspin = PARAM.inp.nspin; 
+    assert (nspin==1 || nspin==2 || nspin==4);
 
     if( ModuleBase::GlobalFunc::SCAN_LINE_BEGIN(ifpos, "ATOMIC_POSITIONS"))
     {
@@ -95,14 +97,22 @@ bool unitcell::read_atom_positions(UnitCell& ucell,
             if ((PARAM.inp.basis_type == "lcao")||(PARAM.inp.basis_type == "lcao_in_pw"))
             {
                 std::string orbital_file = PARAM.inp.orbital_dir + ucell.orbital_fn[it];
-                elecstate::read_orb_file(it, orbital_file, ofs_running, &(ucell.atoms[it]));
+                bool normal = elecstate::read_orb_file(it, orbital_file, ofs_running, &(ucell.atoms[it]));
+				if(!normal)
+				{
+					return false;
+				}
             }
             else if(PARAM.inp.basis_type == "pw")
             {
                 if ((PARAM.inp.init_wfc.substr(0, 3) == "nao") || PARAM.inp.onsite_radius > 0.0)
                 {
                     std::string orbital_file = PARAM.inp.orbital_dir + ucell.orbital_fn[it];
-                    elecstate::read_orb_file(it, orbital_file, ofs_running, &(ucell.atoms[it]));
+					bool normal = elecstate::read_orb_file(it, orbital_file, ofs_running, &(ucell.atoms[it]));
+					if(!normal)
+					{
+						return false;
+					}
                 }
                 else
                 {
@@ -219,14 +229,14 @@ bool unitcell::read_atom_positions(UnitCell& ucell,
                             ifpos.putback(tmp);
                             ifpos >> tmpid;
                         }
-                        if ( tmpid == "m" )
-                        {
-                                ifpos >> mv.x >> mv.y >> mv.z ;
-                        }
-                        else if ( tmpid == "v" ||tmpid == "vel" || tmpid == "velocity" )
-                        {
-                                ifpos >> ucell.atoms[it].vel[ia].x >> ucell.atoms[it].vel[ia].y >> ucell.atoms[it].vel[ia].z;
-                        }
+						if ( tmpid == "m" )
+						{
+							ifpos >> mv.x >> mv.y >> mv.z ;
+						}
+						else if ( tmpid == "v" ||tmpid == "vel" || tmpid == "velocity" )
+						{
+							ifpos >> ucell.atoms[it].vel[ia].x >> ucell.atoms[it].vel[ia].y >> ucell.atoms[it].vel[ia].z;
+						}
                         else if ( tmpid == "mag" || tmpid == "magmom")
                         {
                             set_element_mag_zero = true;
@@ -352,7 +362,7 @@ bool unitcell::read_atom_positions(UnitCell& ucell,
                         ucell.atoms[it].m_loc_[ia].z = ucell.atoms[it].mag[ia];
                     }
 
-                    if(PARAM.inp.nspin==4)
+                    if(nspin==4)
                     {
                         if(!PARAM.inp.noncolin)
                         {
@@ -382,7 +392,7 @@ bool unitcell::read_atom_positions(UnitCell& ucell,
                         }
                         ModuleBase::GlobalFunc::ZEROS(ucell.magnet.ux_ ,3);
                     }
-                    else if(PARAM.inp.nspin==2)
+                    else if(nspin==2)
                     {// collinear case with nspin = 2, only z component is used
                         ucell.atoms[it].mag[ia] = ucell.atoms[it].m_loc_[ia].z;
                         //print only ia==0 && mag>0 to avoid too much output
@@ -514,7 +524,7 @@ bool unitcell::read_atom_positions(UnitCell& ucell,
         }
         if (autoset_mag)
         {
-            if(PARAM.inp.nspin==4)
+            if(nspin==4)
             {
                 for (int it = 0;it < ntype; it++)
                 {
@@ -530,7 +540,7 @@ bool unitcell::read_atom_positions(UnitCell& ucell,
                     }
                 }
             }
-            else if(PARAM.inp.nspin==2)
+            else if(nspin==2)
             {
                 for (int it = 0;it < ntype; it++)
                 {

source/source_estate/read_orb.cpp

@@ -2,7 +2,7 @@
 #include "source_base/formatter.h"
 
 namespace elecstate {
-    void read_orb_file(int it, std::string &orb_file, std::ofstream &ofs_running, Atom* atom)
+    bool read_orb_file(int it, std::string &orb_file, std::ofstream &ofs_running, Atom* atom)
     {
     // the maximum L is 9 like cc-pV9Z, according to the 
     // basissetexchange https://www.basissetexchange.org/
@@ -14,8 +14,9 @@ namespace elecstate {
     {
         std::cout << " Element index " << it+1 << std::endl;
         std::cout << " orbital file: " << orb_file << std::endl;
-        ModuleBase::WARNING_QUIT("UnitCell::read_orb_file", 
-                                "ABACUS Cannot find the ORBITAL file (basis sets)");
+        ModuleBase::WARNING("elecstate::read_orb_file", 
+                                "cannot open the ORBITAL file (NAO basis sets)");
+        return false;
     }
     std::string word;
     atom->nw = 0;
@@ -54,17 +55,19 @@ namespace elecstate {
             }
             if (!valid)
             {
-                ModuleBase::WARNING_QUIT("UnitCell::read_orb_file", 
-                                         "ABACUS does not support numerical atomic orbital with L > 9, "
+                ModuleBase::WARNING("elecstate::read_orb_file", 
+                                         "ABACUS does not support NAO with L > 9, "
                                          "or an invalid orbital label is found in the ORBITAL file.");
+                return false;
             }
         }
     }
     ifs.close();
     if(!atom->nw)
     {
-        ModuleBase::WARNING_QUIT("UnitCell::read_orb_file","get nw = 0");
+		ModuleBase::WARNING("elecstate::read_orb_file","get nw = 0, check the ORBITAL file");
+		return false;
     }
 }
 
-}
+}

source/source_estate/read_orb.h

@@ -13,11 +13,11 @@ namespace elecstate
      * @param ofs_running ofstream
      * @param atom Atom instance stored in UnitCell
     */
-    void read_orb_file(int it,
+    bool read_orb_file(int it,
                        std::string& orb_file,
                        std::ofstream& ofs_running,
                        Atom* atom);
 
 }
 
-#endif
+#endif

source/source_io/cal_test.cpp

@@ -61,12 +61,6 @@ void Cal_Test::test_memory(const int nat,
 	const int ngmw = Cal_Test::cal_np(GGT,wfcpw->ggecut, rhopw->nx, rhopw->ny, rhopw->nz);
 	const int ngmc = Cal_Test::cal_np(GGT,rhopw->ggecut, rhopw->nx, rhopw->ny, rhopw->nz);
 
-//  const int ecut_wfc = INPUT.ecutwfc;
-//  const int ecut_chg = INPUT.ecutrho;
-
-//	const int ngmw = Cal_Test::cal_np(ecut_wfc, rhopw->nx, rhopw->ny, rhopw->nz);
-//	const int ngmc = Cal_Test::cal_np(ecut_chg, rhopw->nx, rhopw->ny, rhopw->nz);
-
 	std::cout << " number of atoms = " << nat << std::endl;
 	std::cout << " plane wave number for wave functions = " << ngmw << std::endl;
 	std::cout << " plane wave number for chage density  = " << ngmc << std::endl;

source/source_io/read_input.cpp

@@ -95,48 +95,56 @@ bool ReadInput::check_mode = false;
 bool filter_nonascii_and_comment(std::ifstream& ifs,
                        std::stringstream& out_ascii_stream)
 {
-    // 
-    if (!ifs.is_open()) {
-        if (!ifs) return false;
+	if (!ifs.is_open()) 
+	{
+		if (!ifs) return false;
     }
 
     std::streampos old_pos = ifs.tellg();
     ifs.clear();
     ifs.seekg(0, std::ios::beg);
 
-    char c;
-    while (ifs.get(c)) {
-        // If comment start, skip until end of line (but keep the newline)
-        if (c == '#') {
-            char d;
-            bool newline_found = false;
-            while (ifs.get(d)) {
-                if (d == '\n' || d == '\r') {
-                    // preserve line break in output
-                    out_ascii_stream.put('\n');
-                    // If CRLF, consume the LF after CR (already wrote a single '\n')
-                    if (d == '\r' && ifs.peek() == '\n') {
-                        ifs.get(d); // consume '\n'
-                    }
-                    newline_found = true;
-                    break;
-                }
+	char c;
+	while (ifs.get(c)) 
+	{
+		// If comment start, skip until end of line (but keep the newline)
+		if (c == '#') 
+		{
+			char d;
+			bool newline_found = false;
+			while (ifs.get(d)) 
+			{
+				if (d == '\n' || d == '\r') 
+				{
+					// preserve line break in output
+					out_ascii_stream.put('\n');
+					// If CRLF, consume the LF after CR (already wrote a single '\n')
+					if (d == '\r' && ifs.peek() == '\n') 
+					{
+						ifs.get(d); // consume '\n'
+					}
+					newline_found = true;
+					break;
+				}
             }
-            if (!newline_found) {
+			if (!newline_found) 
+			{
                 // reached EOF while skipping comment
                 break;
             }
             continue;
         }
 
         unsigned char uc = static_cast<unsigned char>(c);
-        if (uc <= 0x7F) {
-            // ASCII character
+		if (uc <= 0x7F) 
+		{
+			// ASCII character
             out_ascii_stream.put(c);
         }
-        else {
-            // replace non-ASCII with space character
-            out_ascii_stream.put(' ');
+		else 
+		{
+			// replace non-ASCII with space character
+			out_ascii_stream.put(' ');
         }
     }
 
@@ -151,7 +159,7 @@ bool filter_nonascii_and_comment(std::ifstream& ifs,
 ReadInput::ReadInput(const int& rank)
 {
     this->rank = rank;
-    
+
     // add items
     this->item_system();
     this->item_elec_stru();
@@ -174,6 +182,7 @@ ReadInput::ReadInput(const int& rank)
 void ReadInput::read_parameters(Parameter& param, const std::string& filename_in)
 {
     ModuleBase::TITLE("ReadInput", "read_parameters");
+
     // 1. only rank 0 read the input file
     if (this->rank == 0)
     {
@@ -222,9 +231,6 @@ void ReadInput::read_parameters(Parameter& param, const std::string& filename_in
         param.input.kpar = base_device::information::get_device_kpar(param.inp.kpar, param.inp.bndpar);
     }
 
-
-
-
     if (this->check_mode)
     {
         std::cout << "----------------------------------------------------------" << std::endl;
@@ -287,24 +293,22 @@ void ReadInput::read_txt_input(Parameter& param, const std::string& filename)
 
     std::stringstream ascii_stream;
 
-    {
-        std::ifstream ifs(filename.c_str(), std::ios::in);
+	std::ifstream ifs(filename.c_str(), std::ios::in);
 
-        if (!ifs)
-        {
-            std::cout << " Can't find the INPUT file." << std::endl;
-            ModuleBase::WARNING_QUIT("Input::Init", "Error during readin parameters.", 1);
-        }
+	if (!ifs)
+	{
+		std::cout << " Can't find the INPUT file." << std::endl;
+		ModuleBase::WARNING_QUIT("Input::Init", "Error during readin parameters.", 1);
+	}
 
-        ifs.clear();
-        ifs.seekg(0);
+	ifs.clear();
+	ifs.seekg(0);
 
-        filter_nonascii_and_comment(ifs, ascii_stream);
-        ifs.clear();
+	filter_nonascii_and_comment(ifs, ascii_stream);
+	ifs.clear();
+
+	// file close after reading
 
-        // file close after reading
-    }
-
     int ierr = 0;
     ascii_stream.rdstate();
     while (ascii_stream.good())
@@ -393,8 +397,9 @@ void ReadInput::read_txt_input(Parameter& param, const std::string& filename)
     for (auto& input_item: this->input_lists)
     {
         Input_Item* resetvalue_item = &(input_item.second);
-        if (resetvalue_item->reset_value != nullptr) {
-            resetvalue_item->reset_value(*resetvalue_item, param);
+        if (resetvalue_item->reset_value != nullptr) 
+		{
+			resetvalue_item->reset_value(*resetvalue_item, param);
         }
     }
 }

source/source_io/write_dipole.cpp

@@ -64,7 +64,7 @@ void ModuleIO::write_dipole(const UnitCell& ucell,
     Parallel_Reduce::reduce_pool(dipole_elec_y);
     Parallel_Reduce::reduce_pool(dipole_elec_z);
 
-    ofs << istep << " " << dipole_elec_x << " " << dipole_elec_y << dipole_elec_z;
+    ofs << istep+1 << " " << dipole_elec_x << " " << dipole_elec_y << dipole_elec_z;
 #else
 
     double dipole_elec[3] = {0.0, 0.0, 0.0};
@@ -95,7 +95,7 @@ void ModuleIO::write_dipole(const UnitCell& ucell,
     ModuleBase::GlobalFunc::OUT(GlobalV::ofs_running, "Electronic dipole moment P_elec_y(t)", dipole_elec[1]);
     ModuleBase::GlobalFunc::OUT(GlobalV::ofs_running, "Electronic dipole moment P_elec_z(t)", dipole_elec[2]);
 
-    ofs << std::setprecision(precision) << istep << " " << dipole_elec[0] << " " << dipole_elec[1] << " "
+    ofs << std::setprecision(precision) << istep+1 << " " << dipole_elec[0] << " " << dipole_elec[1] << " "
         << dipole_elec[2] << std::endl;
 
     double dipole_ion[3] = {0.0};

source/source_main/driver.cpp

@@ -115,17 +115,17 @@ void Driver::reading()
     GlobalV::NPROC = PARAM.globalv.nproc;
 
     // (1) read the input file
-    ModuleIO::ReadInput read_input(PARAM.globalv.myrank);
-    read_input.read_parameters(PARAM,  PARAM.globalv.global_in_card);
+    ModuleIO::ReadInput input(PARAM.globalv.myrank);
+    input.read_parameters(PARAM, PARAM.globalv.global_in_card);
 
     // (2) create the output directory, running_*.log and print info
-    read_input.create_directory(PARAM);
+    input.create_directory(PARAM);
     this->print_start_info();
 
     // (3) write the input file
     std::stringstream ss1;
     ss1 << PARAM.globalv.global_out_dir <<  PARAM.globalv.global_in_card;
-    read_input.write_parameters(PARAM, ss1.str());
+    input.write_parameters(PARAM, ss1.str());
 
     // (*temp*) copy the variables from INPUT to each class
     Input_Conv::Convert();

source/source_pw/module_ofdft/evolve_ofdft.cpp

@@ -223,8 +223,9 @@ void Evolve_OFDFT::cal_CD_potential(std::vector<std::complex<double>>& psi_,
 
         for (int ir = 0; ir < pw_rho->nrxx; ++ir)
         {
-            rpot(0, ir) -= mCD_para*2.0*std::real(rCDPotential[ir])*std::pow(ModuleBase::PI,3) / (2.0*std::pow(std::real(kF_r[ir]),2));
-            if (isnan(rpot(0, ir))) 
+            rpot(0, ir) -= mCD_para*2.0*std::real(rCDPotential[ir])*std::pow(ModuleBase::PI,3) 
+                        / (2.0*std::pow(std::real(kF_r[ir]),2));
+            if (std::isnan(rpot(0, ir))) 
             {
                 rpot(0, ir)=0.0;
             }