update the update_pos_taud

A-006 · A-006 · commit c4676b08b05b · 2025-01-02T18:41:19.000+08:00
diff --git a/source/module_cell/test/support/mock_unitcell.cpp b/source/module_cell/test/support/mock_unitcell.cpp
@@ -33,7 +33,6 @@ bool UnitCell::read_atom_positions(std::ifstream& ifpos,
                                    std::ofstream& ofs_warning) {
     return true;
 }
-void UnitCell::update_pos_taud(double* posd_in) {}
 void UnitCell::update_pos_taud(const ModuleBase::Vector3<double>* posd_in) {}
 void UnitCell::update_vel(const ModuleBase::Vector3<double>* vel_in) {}
 void UnitCell::bcast_atoms_tau() {}
diff --git a/source/module_cell/test/unitcell_test_para.cpp b/source/module_cell/test/unitcell_test_para.cpp
@@ -167,7 +167,8 @@ TEST_F(UcellTest, UpdatePosTaud)
         ucell->iat2iait(iat, &ia, &it);
         tmp[iat] = ucell->atoms[it].taud[ia];
     }
-    ucell->update_pos_taud(pos_in);
+    unitcell::update_pos_taud(ucell->lat,pos_in,ucell->ntype,
+                              ucell->nat,ucell->atoms);
     for (int iat = 0; iat < ucell->nat; ++iat)
     {
         int it, ia;
diff --git a/source/module_cell/unitcell.cpp b/source/module_cell/unitcell.cpp
@@ -314,25 +314,6 @@ std::vector<ModuleBase::Vector3<int>> UnitCell::get_constrain() const
 	return constrain;
 }
 
-
-
-void UnitCell::update_pos_taud(double* posd_in) {
-    int iat = 0;
-    for (int it = 0; it < this->ntype; it++) {
-        Atom* atom = &this->atoms[it];
-        for (int ia = 0; ia < atom->na; ia++) {
-            for (int ik = 0; ik < 3; ++ik) {
-                atom->taud[ia][ik] += posd_in[3 * iat + ik];
-                atom->dis[ia][ik] = posd_in[3 * iat + ik];
-            }
-            iat++;
-        }
-    }
-    assert(iat == this->nat);
-    unitcell::periodic_boundary_adjustment(this->atoms,this->latvec, this->ntype);
-    this->bcast_atoms_tau();
-}
-
 // posd_in is atomic displacements here  liuyu 2023-03-22
 void UnitCell::update_pos_taud(const ModuleBase::Vector3<double>* posd_in) {
     int iat = 0;
diff --git a/source/module_cell/unitcell.h b/source/module_cell/unitcell.h
@@ -201,7 +201,6 @@ class UnitCell {
     void print_cell_xyz(const std::string& fn) const;
 
     void update_pos_taud(const ModuleBase::Vector3<double>* posd_in);
-    void update_pos_taud(double* posd_in);
     void update_vel(const ModuleBase::Vector3<double>* vel_in);
     void bcast_atoms_tau();
     bool judge_big_cell() const;
diff --git a/source/module_cell/update_cell.cpp b/source/module_cell/update_cell.cpp
@@ -374,6 +374,32 @@ void update_pos_tau(const Lattice& lat,
     bcast_atoms_tau(atoms, ntype);
 }
 
+void update_pos_taud(const Lattice& lat,
+                     const double* posd_in,
+                     const int ntype,
+                     const int nat,
+                     Atom* atoms)
+{
+    int iat = 0;
+    for (int it = 0; it < ntype; it++) 
+    {
+        Atom* atom = &atoms[it];
+        for (int ia = 0; ia < atom->na; ia++) 
+        {
+            for (int ik = 0; ik < 3; ++ik) 
+            {
+                atom->taud[ia][ik] += posd_in[3 * iat + ik];
+                atom->dis[ia][ik] = posd_in[3 * iat + ik];
+            }
+            iat++;
+        }
+    }
+    assert(iat == nat);
+    periodic_boundary_adjustment(atoms,lat.latvec,ntype);
+    bcast_atoms_tau(atoms, ntype);
+}
+
+
 void periodic_boundary_adjustment(Atom* atoms,
                                   const ModuleBase::Matrix3& latvec,
                                   const int ntype) 
@@ -384,12 +410,6 @@ void periodic_boundary_adjustment(Atom* atoms,
     // first adjust direct coordinates,
     // then update them into cartesian coordinates,
     //----------------------------------------------
-    for (int i=0;i<ntype;i++) {
-        Atom* atom = &atoms[i];
-        for (int j=0;j<atom->na;j++) {
-            printf("the taud is %f %f %f\n",atom->taud[j].x,atom->taud[j].y,atom->taud[j].z);
-        }
-    }
     for (int it = 0; it < ntype; it++) {
         Atom* atom = &atoms[it];
         for (int ia = 0; ia < atom->na; ia++) {
diff --git a/source/module_cell/update_cell.h b/source/module_cell/update_cell.h
@@ -48,6 +48,31 @@ namespace unitcell
                         const int ntype,
                         const int nat,
                         Atom* atoms);
+    
+    /**
+     * @brief update the position and tau of the atoms
+     * 
+     * @param lat: the lattice of the atoms [in]
+     * @param pos_in: the position of the atoms in direct coordinate system [in]
+     * @param ntype: the number of types of the atoms [in]
+     * @param nat: the number of atoms [in]
+     * @param atoms: the atoms to be updated [out]
+     */
+    void update_pos_taud(const Lattice& lat,
+                         const double* posd_in,
+                         const int ntype,
+                         const int nat,
+                         Atom* atoms);
+    /**
+     * @brief update the velocity of the atoms
+     * 
+     * @param lat: the lattice of the atoms [in]
+     * @param pos_in: the position of the atoms in direct coordinate system 
+     * in ModuleBase::Vector3 version [in]
+     * @param ntype: the number of types of the atoms [in]
+     * @param nat: the number of atoms [in]
+     * @param atoms: the atoms to be updated [out]
+     */
 }
 //
 #endif // UPDATE_CELL_H
diff --git a/source/module_relax/relax_new/relax.cpp b/source/module_relax/relax_new/relax.cpp
@@ -633,7 +633,7 @@ void Relax::move_cell_ions(UnitCell& ucell, const bool is_new_dir)
         ucell.symm.symmetrize_vec3_nat(move_ion);
     }
 
-    ucell.update_pos_taud(move_ion);
+    unitcell::update_pos_taud(ucell.lat,move_ion,ucell.ntype,ucell.nat,ucell.atoms);
 
     // Print the structure file.
     ucell.print_tau();
diff --git a/source/module_relax/relax_new/test/relax_test.cpp b/source/module_relax/relax_new/test/relax_test.cpp
@@ -1,4 +1,5 @@
 #include "gtest/gtest.h"
+#include <iomanip>
 #define private public
 #include "module_parameter/parameter.h"
 #undef private
@@ -27,13 +28,13 @@ class Test_SETGRAD : public testing::Test
             force_in.create(nat,3);
             stress_in.create(3,3);
 
-            force_in(0,0) = 1; force_in(0,1) = 2; force_in(0,2)= 3;
-            force_in(1,0) = 4; force_in(1,1) = 5; force_in(1,2)= 6;
-            force_in(2,0) = 7; force_in(2,1) = 8; force_in(2,2)= 9;
+            force_in(0,0) = 0.1; force_in(0,1) = 0.1; force_in(0,2)= 0.1;
+            force_in(1,0) = 0; force_in(1,1) = 0.1; force_in(1,2)= 0.1;
+            force_in(2,0) = 0; force_in(2,1) = 0; force_in(2,2)= 0.1;
 
-            stress_in(0,0) = 1; stress_in(0,1) = 2; stress_in(0,2)= 3;
-            stress_in(1,0) = 4; stress_in(1,1) = 5; stress_in(1,2)= 6;
-            stress_in(2,0) = 7; stress_in(2,1) = 8; stress_in(2,2)= 9;
+            stress_in(0,0) = 1; stress_in(0,1) = 1; stress_in(0,2)= 1;
+            stress_in(1,0) = 0; stress_in(1,1) = 1; stress_in(1,2)= 1;
+            stress_in(2,0) = 0; stress_in(2,1) = 0; stress_in(2,2)= 1;
 
             ucell.ntype = 1;
             ucell.nat = nat;
@@ -46,6 +47,8 @@ class Test_SETGRAD : public testing::Test
             ucell.iat2ia = new int[nat];
             ucell.atoms[0].mbl.resize(nat);
             ucell.atoms[0].taud.resize(nat);
+            ucell.atoms[0].tau.resize(nat);
+            ucell.atoms[0].dis.resize(nat);
             ucell.lc = new int[3];
 
             ucell.iat2it[0] = 0;
@@ -136,14 +139,25 @@ TEST_F(Test_SETGRAD, relax_new)
 {
     std::vector<double> result_ref = 
     {
-        0,0,0.1709672056,0,0.2849453427,0,0.3989234797,0,0,1.005319517,
-        0.01063903455,0.01595855183,0.0212780691,1.026597586,0.03191710366,
-        0.03723662093,0.04255613821,1.047875655,1.059181731,0,0,0,1.059181731,
-        0,0,0,1.059181731,1.034363264,0.01301504537,0.01952256806,
-        0.02603009074,1.060393355,0.03904513611,0.0455526588,0.05206018148,
-        1.086423445,1,0,0,0,1,0,0,0,1
+        0,           0,            0.24293434145, 
+        0,           0.242934341453,           0,
+        0,           0,           0,
+        //paramter for taud
+        1.2267616333,0.2267616333,0.22676163333, 
+        0,            1.2267616333 ,0.2267616333,
+        0,            0,           1.22676163333,
+        // paramter for fisrt time after relaxation
+        1.3677603495, 0,            0,
+        0,            1.36776034956, 0,
+        0,            0,            1.36776034956,
+        // paramter for second time after relaxation
+        1.3677603495  ,0.3633367476,0.36333674766,
+        0,            1.3677603495 ,0.36333674766,
+        0,            0,            1.3677603495 ,
+        // paramter for third time after relaxation
+        1,0,0,0,1,0,0,0,1
+        // paramter for fourth time after relaxation
     };
-
     for(int i=0;i<result.size();i++)
     {
         EXPECT_NEAR(result[i],result_ref[i],1e-8);
@@ -258,6 +272,7 @@ class Test_RELAX : public testing::Test
                 ucell.atoms[i].mbl.resize(na);
                 ucell.atoms[i].taud.resize(na);
                 ucell.atoms[i].tau.resize(na);
+                ucell.atoms[i].dis.resize(na);
                 for (int j=0;j<na;j++)
                 {
                     ucell.atoms[i].mbl[j] = {1,1,1};
@@ -291,7 +306,6 @@ TEST_F(Test_RELAX, relax_new)
     int size = 1584;
     double tmp;
     std::ifstream result_ref("./support/result_ref.txt");
-
     for(int i=0;i<size;i++)
     {
         result_ref >> tmp;
diff --git a/source/module_relax/relax_new/test/relax_test.h b/source/module_relax/relax_new/test/relax_test.h
@@ -2,30 +2,13 @@
 
 namespace GlobalC
 {
-	UnitCell ucell;
     Structure_Factor sf;
     ModulePW::PW_Basis* rhopw;
 }
 
 UnitCell::UnitCell(){};
 UnitCell::~UnitCell(){};
 
-void UnitCell::update_pos_taud(double* posd_in)
-{
-    int iat = 0;
-    for (int it = 0; it < this->ntype; it++)
-    {
-        Atom* atom = &this->atoms[it];
-        for (int ia = 0; ia < atom->na; ia++)
-        {
-            this->atoms[it].taud[ia].x += posd_in[iat*3];
-            this->atoms[it].taud[ia].y += posd_in[iat*3 + 1];
-            this->atoms[it].taud[ia].z += posd_in[iat*3 + 2];
-            iat++;
-        }
-    }
-    assert(iat == this->nat);
-}
 
 void UnitCell::print_stru_file(const std::string& fn, 
                                const int& nspin,
diff --git a/source/module_relax/relax_new/test/support/result_ref.txt b/source/module_relax/relax_new/test/support/result_ref.txt

Original file line number	Diff line number	Diff line change
`@@ -33,7 +33,6 @@ bool UnitCell::read_atom_positions(std::ifstream& ifpos,`
`33`	`33`	`std::ofstream& ofs_warning) {`
`34`	`34`	`return true;`
`35`	`35`	`}`
`36`		`-void UnitCell::update_pos_taud(double* posd_in) {}`
`37`	`36`	`void UnitCell::update_pos_taud(const ModuleBase::Vector3<double>* posd_in) {}`
`38`	`37`	`void UnitCell::update_vel(const ModuleBase::Vector3<double>* vel_in) {}`
`39`	`38`	`void UnitCell::bcast_atoms_tau() {}`
Original file line number	Diff line number	Diff line change
`@@ -167,7 +167,8 @@ TEST_F(UcellTest, UpdatePosTaud)`
`167`	`167`	`ucell->iat2iait(iat, &ia, &it);`
`168`	`168`	`tmp[iat] = ucell->atoms[it].taud[ia];`
`169`	`169`	`}`
`170`		`- ucell->update_pos_taud(pos_in);`
	`170`	`+ unitcell::update_pos_taud(ucell->lat,pos_in,ucell->ntype,`
	`171`	`+ ucell->nat,ucell->atoms);`
`171`	`172`	`for (int iat = 0; iat < ucell->nat; ++iat)`
`172`	`173`	`{`
`173`	`174`	`int it, ia;`
Original file line number	Diff line number	Diff line change
`@@ -633,7 +633,7 @@ void Relax::move_cell_ions(UnitCell& ucell, const bool is_new_dir)`
`633`	`633`	`ucell.symm.symmetrize_vec3_nat(move_ion);`
`634`	`634`	`}`
`635`	`635`
`636`		`- ucell.update_pos_taud(move_ion);`
	`636`	`+ unitcell::update_pos_taud(ucell.lat,move_ion,ucell.ntype,ucell.nat,ucell.atoms);`
`637`	`637`
`638`	`638`	`// Print the structure file.`
`639`	`639`	`ucell.print_tau();`