diff --git a/source/Makefile.Objects b/source/Makefile.Objects
index ae6a7bf00b..75c750e654 100644
--- a/source/Makefile.Objects
+++ b/source/Makefile.Objects
@@ -187,6 +187,7 @@ OBJS_CELL=atom_pseudo.o\
     klist.o\
     cell_index.o\
     check_atomic_stru.o\
+    update_cell.o\
 
 OBJS_DEEPKS=LCAO_deepks.o\
         deepks_force.o\
diff --git a/source/module_cell/CMakeLists.txt b/source/module_cell/CMakeLists.txt
index c658450e3f..f4cb5e5991 100644
--- a/source/module_cell/CMakeLists.txt
+++ b/source/module_cell/CMakeLists.txt
@@ -23,6 +23,7 @@ add_library(
     parallel_kpoints.cpp
     cell_index.cpp
     check_atomic_stru.cpp
+    update_cell.cpp
 )
 
 if(ENABLE_COVERAGE)
diff --git a/source/module_cell/test/CMakeLists.txt b/source/module_cell/test/CMakeLists.txt
index c79adcf89d..890b4b776b 100644
--- a/source/module_cell/test/CMakeLists.txt
+++ b/source/module_cell/test/CMakeLists.txt
@@ -103,7 +103,7 @@ add_test(NAME cell_parallel_kpoints_test
 AddTest(
   TARGET cell_unitcell_test
   LIBS parameter  ${math_libs} base device cell_info  symmetry
-  SOURCES unitcell_test.cpp ../../module_io/output.cpp ../../module_elecstate/cal_ux.cpp
+  SOURCES unitcell_test.cpp ../../module_io/output.cpp ../../module_elecstate/cal_ux.cpp ../update_cell.cpp
 
 )
 
@@ -122,7 +122,8 @@ AddTest(
 AddTest(
   TARGET cell_unitcell_test_setupcell
   LIBS parameter  ${math_libs} base device cell_info 
-  SOURCES unitcell_test_setupcell.cpp ../../module_io/output.cpp
+  SOURCES unitcell_test_setupcell.cpp ../../module_io/output.cpp 
+  ../../module_cell/update_cell.cpp
 )
 
 add_test(NAME cell_unitcell_test_parallel
diff --git a/source/module_cell/test/prepare_unitcell.h b/source/module_cell/test/prepare_unitcell.h
index ce6cc65eb9..3890670e18 100644
--- a/source/module_cell/test/prepare_unitcell.h
+++ b/source/module_cell/test/prepare_unitcell.h
@@ -109,7 +109,7 @@ class UcellTestPrepare
 		ucell->latvec.e13 = this->latvec[2];
 		ucell->latvec.e21 = this->latvec[3];
 		ucell->latvec.e22 = this->latvec[4];
-	       	ucell->latvec.e23 = this->latvec[5];
+	    ucell->latvec.e23 = this->latvec[5];
 		ucell->latvec.e31 = this->latvec[6];
 		ucell->latvec.e32 = this->latvec[7];
 		ucell->latvec.e33 = this->latvec[8];
diff --git a/source/module_cell/test/support/mock_unitcell.cpp b/source/module_cell/test/support/mock_unitcell.cpp
index 9a4c28bb65..dfb8e4cedd 100644
--- a/source/module_cell/test/support/mock_unitcell.cpp
+++ b/source/module_cell/test/support/mock_unitcell.cpp
@@ -63,8 +63,6 @@ void UnitCell::print_stru_file(const std::string& fn,
                                const bool& dpks_desc,
                                const int& iproc) const {}
 void UnitCell::check_dtau() {}
-void UnitCell::setup_cell_after_vc(std::ofstream& log) {}
-void UnitCell::remake_cell() {}
 void UnitCell::cal_nwfc(std::ofstream& log) {}
 void UnitCell::cal_meshx() {}
 void UnitCell::cal_natomwfc(std::ofstream& log) {}
diff --git a/source/module_cell/test/unitcell_test.cpp b/source/module_cell/test/unitcell_test.cpp
index 88052a433b..036bbf40ed 100644
--- a/source/module_cell/test/unitcell_test.cpp
+++ b/source/module_cell/test/unitcell_test.cpp
@@ -3,7 +3,6 @@
 #define private public
 #include "module_parameter/parameter.h"
 #undef private
-#include "module_cell/unitcell.h"
 #include "module_elecstate/cal_ux.h"
 #include "module_elecstate/read_pseudo.h"
 
@@ -11,6 +10,7 @@
 #include "module_base/global_variable.h"
 #include "module_base/mathzone.h"
 #include "prepare_unitcell.h"
+#include "module_cell/update_cell.h"
 #include <streambuf>
 #include <valarray>
 #include <vector>
@@ -346,7 +346,7 @@ TEST_F(UcellTest, RemakeCell)
         ucell->latvec.e32 = 0.00;
         ucell->latvec.e33 = 10.0;
         ucell->latName = latname_in[i];
-        ucell->remake_cell();
+        unitcell::remake_cell(ucell->lat);
         if (latname_in[i] == "sc")
         {
             double celldm
@@ -591,7 +591,7 @@ TEST_F(UcellDeathTest, RemakeCellWarnings)
         ucell->latvec.e33 = 10.0;
         ucell->latName = latname_in[i];
         testing::internal::CaptureStdout();
-        EXPECT_EXIT(ucell->remake_cell(), ::testing::ExitedWithCode(1), "");
+        EXPECT_EXIT(unitcell::remake_cell(ucell->lat), ::testing::ExitedWithCode(1), "");
         std::string output = testing::internal::GetCapturedStdout();
         if (latname_in[i] == "none")
         {
diff --git a/source/module_cell/test/unitcell_test_setupcell.cpp b/source/module_cell/test/unitcell_test_setupcell.cpp
index 6b2029adce..7e832c2d99 100644
--- a/source/module_cell/test/unitcell_test_setupcell.cpp
+++ b/source/module_cell/test/unitcell_test_setupcell.cpp
@@ -11,7 +11,7 @@
 #include<valarray>
 #include <streambuf>
 #include "prepare_unitcell.h"
-
+#include "module_cell/update_cell.h"
 #ifdef __LCAO
 #include "module_basis/module_ao/ORB_read.h"
 InfoNonlocal::InfoNonlocal(){}
@@ -145,8 +145,47 @@ TEST_F(UcellTest,SetupCellAfterVC)
 	ucell->ntype = 2;
 	delete[] ucell->magnet.start_magnetization;
 	ucell->magnet.start_magnetization = new double[ucell->ntype];
+
+	
 	ucell->setup_cell(fn,ofs_running);
-	ucell->setup_cell_after_vc(ofs_running);
+	ucell->lat0 = 1.0;
+	ucell->latvec.Zero();
+	ucell->latvec.e11 = 10.0;
+	ucell->latvec.e22 = 10.0;
+	ucell->latvec.e33 = 10.0;
+	for (int i =0;i<ucell->ntype;i++)
+	{
+		ucell->atoms[i].na = 1;
+		ucell->atoms[i].taud.resize(ucell->atoms[i].na);
+		ucell->atoms[i].tau.resize(ucell->atoms[i].na);
+		ucell->atoms[i].taud[0].x = 0.1;
+		ucell->atoms[i].taud[0].y = 0.1;
+		ucell->atoms[i].taud[0].z = 0.1;
+	}
+	
+	unitcell::setup_cell_after_vc(*ucell,ofs_running);
+	EXPECT_EQ(ucell->lat0_angstrom,0.529177);
+	EXPECT_EQ(ucell->tpiba,ModuleBase::TWO_PI);
+	EXPECT_EQ(ucell->tpiba2,ModuleBase::TWO_PI*ModuleBase::TWO_PI);
+	EXPECT_EQ(ucell->a1.x ,10.0);
+	EXPECT_EQ(ucell->a2.y ,10.0);
+	EXPECT_EQ(ucell->a3.z ,10.0);
+	EXPECT_EQ(ucell->omega,1000.0);
+	EXPECT_EQ(ucell->GT.e11,0.1);
+	EXPECT_EQ(ucell->GT.e22,0.1);
+	EXPECT_EQ(ucell->GT.e33,0.1);
+	EXPECT_EQ(ucell->G.e11,0.1);
+	EXPECT_EQ(ucell->G.e22,0.1);
+	EXPECT_EQ(ucell->G.e33,0.1);
+
+	for (int it = 0; it < ucell->ntype; it++) {
+        Atom* atom = &ucell->atoms[it];
+        for (int ia = 0; ia < atom->na; ia++) {
+            EXPECT_EQ(atom->tau[ia].x,1);
+			EXPECT_EQ(atom->tau[ia].y,1);
+			EXPECT_EQ(atom->tau[ia].z,1);
+        }
+    }
 	ofs_running.close();
 	remove("setup_cell.tmp");
 }
diff --git a/source/module_cell/unitcell.cpp b/source/module_cell/unitcell.cpp
index 8cabc76149..635b3624c1 100755
--- a/source/module_cell/unitcell.cpp
+++ b/source/module_cell/unitcell.cpp
@@ -863,72 +863,6 @@ void UnitCell::cal_natomwfc(std::ofstream& log) {
     return;
 }
 
-// LiuXh add a new function here,
-// 20180515
-void UnitCell::setup_cell_after_vc(std::ofstream& log) {
-    ModuleBase::TITLE("UnitCell", "setup_cell_after_vc");
-    assert(lat0 > 0.0);
-    this->omega = std::abs(latvec.Det()) * this->lat0 * lat0 * lat0;
-    if (this->omega <= 0) {
-        ModuleBase::WARNING_QUIT("setup_cell_after_vc", "omega <= 0 .");
-    } else {
-        log << std::endl;
-        ModuleBase::GlobalFunc::OUT(log, "Volume (Bohr^3)", this->omega);
-        ModuleBase::GlobalFunc::OUT(log,
-                                    "Volume (A^3)",
-                                    this->omega
-                                        * pow(ModuleBase::BOHR_TO_A, 3));
-    }
-
-    lat0_angstrom = lat0 * 0.529177;
-    tpiba = ModuleBase::TWO_PI / lat0;
-    tpiba2 = tpiba * tpiba;
-
-    // lattice vectors in another form.
-    a1.x = latvec.e11;
-    a1.y = latvec.e12;
-    a1.z = latvec.e13;
-
-    a2.x = latvec.e21;
-    a2.y = latvec.e22;
-    a2.z = latvec.e23;
-
-    a3.x = latvec.e31;
-    a3.y = latvec.e32;
-    a3.z = latvec.e33;
-
-    //==========================================================
-    // Calculate recip. lattice vectors and dot products
-    // latvec has the unit of lat0, but G has the unit 2Pi/lat0
-    //==========================================================
-    this->GT = latvec.Inverse();
-    this->G = GT.Transpose();
-    this->GGT = G * GT;
-    this->invGGT = GGT.Inverse();
-
-    for (int it = 0; it < ntype; it++) {
-        Atom* atom = &atoms[it];
-        for (int ia = 0; ia < atom->na; ia++) {
-            atom->tau[ia] = atom->taud[ia] * latvec;
-        }
-    }
-
-#ifdef __MPI
-    this->bcast_unitcell();
-#endif
-
-    log << std::endl;
-    output::printM3(log,
-                    "Lattice vectors: (Cartesian coordinate: in unit of a_0)",
-                    latvec);
-    output::printM3(
-        log,
-        "Reciprocal vectors: (Cartesian coordinate: in unit of 2 pi/a_0)",
-        G);
-
-    return;
-}
-
 // check if any atom can be moved
 bool UnitCell::if_atoms_can_move() const {
     for (int it = 0; it < this->ntype; it++) {
@@ -1020,282 +954,6 @@ void UnitCell::setup(const std::string& latname_in,
     return;
 }
 
-void UnitCell::remake_cell() {
-    ModuleBase::TITLE("UnitCell", "rmake_cell");
-
-    // The idea is as follows: for each type of lattice, first calculate
-    // from current latvec the lattice parameters, then use the parameters
-    // to reconstruct latvec
-
-    if (latName == "none") {
-        ModuleBase::WARNING_QUIT(
-            "UnitCell",
-            "to use fixed_ibrav, latname must be provided");
-    } else if (latName == "sc") // ibrav = 1
-    {
-        double celldm = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
-                                  + pow(latvec.e13, 2));
-
-        latvec.Zero();
-        latvec.e11 = latvec.e22 = latvec.e33 = celldm;
-    } else if (latName == "fcc") // ibrav = 2
-    {
-        double celldm = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
-                                  + pow(latvec.e13, 2))
-                        / std::sqrt(2.0);
-
-        latvec.e11 = -celldm;
-        latvec.e12 = 0.0;
-        latvec.e13 = celldm;
-        latvec.e21 = 0.0;
-        latvec.e22 = celldm;
-        latvec.e23 = celldm;
-        latvec.e31 = -celldm;
-        latvec.e32 = celldm;
-        latvec.e33 = 0.0;
-    } else if (latName == "bcc") // ibrav = 3
-    {
-        double celldm = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
-                                  + pow(latvec.e13, 2))
-                        / std::sqrt(3.0);
-
-        latvec.e11 = celldm;
-        latvec.e12 = celldm;
-        latvec.e13 = celldm;
-        latvec.e21 = -celldm;
-        latvec.e22 = celldm;
-        latvec.e23 = celldm;
-        latvec.e31 = -celldm;
-        latvec.e32 = -celldm;
-        latvec.e33 = celldm;
-    } else if (latName == "hexagonal") // ibrav = 4
-    {
-        double celldm1 = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
-                                   + pow(latvec.e13, 2));
-        double celldm3 = std::sqrt(pow(latvec.e31, 2) + pow(latvec.e32, 2)
-                                   + pow(latvec.e33, 2));
-        double e22 = sqrt(3.0) / 2.0;
-
-        latvec.e11 = celldm1;
-        latvec.e12 = 0.0;
-        latvec.e13 = 0.0;
-        latvec.e21 = -0.5 * celldm1;
-        latvec.e22 = celldm1 * e22;
-        latvec.e23 = 0.0;
-        latvec.e31 = 0.0;
-        latvec.e32 = 0.0;
-        latvec.e33 = celldm3;
-    } else if (latName == "trigonal") // ibrav = 5
-    {
-        double celldm1 = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
-                                   + pow(latvec.e13, 2));
-        double celldm2 = std::sqrt(pow(latvec.e21, 2) + pow(latvec.e22, 2)
-                                   + pow(latvec.e23, 2));
-        double celldm12 = (latvec.e11 * latvec.e21 + latvec.e12 * latvec.e22
-                           + latvec.e13 * latvec.e23);
-        double cos12 = celldm12 / celldm1 / celldm2;
-
-        if (cos12 <= -0.5 || cos12 >= 1.0) {
-            ModuleBase::WARNING_QUIT("unitcell", "wrong cos12!");
-        }
-        double t1 = sqrt(1.0 + 2.0 * cos12);
-        double t2 = sqrt(1.0 - cos12);
-
-        double e11 = celldm1 * t2 / sqrt(2.0);
-        double e12 = -celldm1 * t2 / sqrt(6.0);
-        double e13 = celldm1 * t1 / sqrt(3.0);
-        double e22 = celldm1 * sqrt(2.0) * t2 / sqrt(3.0);
-
-        latvec.e11 = e11;
-        latvec.e12 = e12;
-        latvec.e13 = e13;
-        latvec.e21 = 0.0;
-        latvec.e22 = e22;
-        latvec.e23 = e13;
-        latvec.e31 = -e11;
-        latvec.e32 = e12;
-        latvec.e33 = e13;
-    } else if (latName == "st") // ibrav = 6
-    {
-        double celldm1 = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
-                                   + pow(latvec.e13, 2));
-        double celldm3 = std::sqrt(pow(latvec.e31, 2) + pow(latvec.e32, 2)
-                                   + pow(latvec.e33, 2));
-        latvec.e11 = celldm1;
-        latvec.e12 = 0.0;
-        latvec.e13 = 0.0;
-        latvec.e21 = 0.0;
-        latvec.e22 = celldm1;
-        latvec.e23 = 0.0;
-        latvec.e31 = 0.0;
-        latvec.e32 = 0.0;
-        latvec.e33 = celldm3;
-    } else if (latName == "bct") // ibrav = 7
-    {
-        double celldm1 = std::abs(latvec.e11);
-        double celldm2 = std::abs(latvec.e13);
-
-        latvec.e11 = celldm1;
-        latvec.e12 = -celldm1;
-        latvec.e13 = celldm2;
-        latvec.e21 = celldm1;
-        latvec.e22 = celldm1;
-        latvec.e23 = celldm2;
-        latvec.e31 = -celldm1;
-        latvec.e32 = -celldm1;
-        latvec.e33 = celldm2;
-    } else if (latName == "so") // ibrav = 8
-    {
-        double celldm1 = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
-                                   + pow(latvec.e13, 2));
-        double celldm2 = std::sqrt(pow(latvec.e21, 2) + pow(latvec.e22, 2)
-                                   + pow(latvec.e23, 2));
-        double celldm3 = std::sqrt(pow(latvec.e31, 2) + pow(latvec.e32, 2)
-                                   + pow(latvec.e33, 2));
-
-        latvec.e11 = celldm1;
-        latvec.e12 = 0.0;
-        latvec.e13 = 0.0;
-        latvec.e21 = 0.0;
-        latvec.e22 = celldm2;
-        latvec.e23 = 0.0;
-        latvec.e31 = 0.0;
-        latvec.e32 = 0.0;
-        latvec.e33 = celldm3;
-    } else if (latName == "baco") // ibrav = 9
-    {
-        double celldm1 = std::abs(latvec.e11);
-        double celldm2 = std::abs(latvec.e22);
-        double celldm3 = std::abs(latvec.e33);
-
-        latvec.e11 = celldm1;
-        latvec.e12 = celldm2;
-        latvec.e13 = 0.0;
-        latvec.e21 = -celldm1;
-        latvec.e22 = celldm2;
-        latvec.e23 = 0.0;
-        latvec.e31 = 0.0;
-        latvec.e32 = 0.0;
-        latvec.e33 = celldm3;
-    } else if (latName == "fco") // ibrav = 10
-    {
-        double celldm1 = std::abs(latvec.e11);
-        double celldm2 = std::abs(latvec.e22);
-        double celldm3 = std::abs(latvec.e33);
-
-        latvec.e11 = celldm1;
-        latvec.e12 = 0.0;
-        latvec.e13 = celldm3;
-        latvec.e21 = celldm1;
-        latvec.e22 = celldm2;
-        latvec.e23 = 0.0;
-        latvec.e31 = 0.0;
-        latvec.e32 = celldm2;
-        latvec.e33 = celldm3;
-    } else if (latName == "bco") // ibrav = 11
-    {
-        double celldm1 = std::abs(latvec.e11);
-        double celldm2 = std::abs(latvec.e12);
-        double celldm3 = std::abs(latvec.e13);
-
-        latvec.e11 = celldm1;
-        latvec.e12 = celldm2;
-        latvec.e13 = celldm3;
-        latvec.e21 = -celldm1;
-        latvec.e22 = celldm2;
-        latvec.e23 = celldm3;
-        latvec.e31 = -celldm1;
-        latvec.e32 = -celldm2;
-        latvec.e33 = celldm3;
-    } else if (latName == "sm") // ibrav = 12
-    {
-        double celldm1 = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
-                                   + pow(latvec.e13, 2));
-        double celldm2 = std::sqrt(pow(latvec.e21, 2) + pow(latvec.e22, 2)
-                                   + pow(latvec.e23, 2));
-        double celldm3 = std::sqrt(pow(latvec.e31, 2) + pow(latvec.e32, 2)
-                                   + pow(latvec.e33, 2));
-        double celldm12 = (latvec.e11 * latvec.e21 + latvec.e12 * latvec.e22
-                           + latvec.e13 * latvec.e23);
-        double cos12 = celldm12 / celldm1 / celldm2;
-
-        double e21 = celldm2 * cos12;
-        double e22 = celldm2 * std::sqrt(1.0 - cos12 * cos12);
-
-        latvec.e11 = celldm1;
-        latvec.e12 = 0.0;
-        latvec.e13 = 0.0;
-        latvec.e21 = e21;
-        latvec.e22 = e22;
-        latvec.e23 = 0.0;
-        latvec.e31 = 0.0;
-        latvec.e32 = 0.0;
-        latvec.e33 = celldm3;
-    } else if (latName == "bacm") // ibrav = 13
-    {
-        double celldm1 = std::abs(latvec.e11);
-        double celldm2 = std::sqrt(pow(latvec.e21, 2) + pow(latvec.e22, 2)
-                                   + pow(latvec.e23, 2));
-        double celldm3 = std::abs(latvec.e13);
-
-        double cos12 = latvec.e21 / celldm2;
-        if (cos12 >= 1.0) {
-            ModuleBase::WARNING_QUIT("unitcell", "wrong cos12!");
-        }
-
-        double e21 = celldm2 * cos12;
-        double e22 = celldm2 * std::sqrt(1.0 - cos12 * cos12);
-
-        latvec.e11 = celldm1;
-        latvec.e12 = 0.0;
-        latvec.e13 = -celldm3;
-        latvec.e21 = e21;
-        latvec.e22 = e22;
-        latvec.e23 = 0.0;
-        latvec.e31 = celldm1;
-        latvec.e32 = 0.0;
-        latvec.e33 = celldm3;
-    } else if (latName == "triclinic") // ibrav = 14
-    {
-        double celldm1 = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
-                                   + pow(latvec.e13, 2));
-        double celldm2 = std::sqrt(pow(latvec.e21, 2) + pow(latvec.e22, 2)
-                                   + pow(latvec.e23, 2));
-        double celldm3 = std::sqrt(pow(latvec.e31, 2) + pow(latvec.e32, 2)
-                                   + pow(latvec.e33, 2));
-        double celldm12 = (latvec.e11 * latvec.e21 + latvec.e12 * latvec.e22
-                           + latvec.e13 * latvec.e23);
-        double cos12 = celldm12 / celldm1 / celldm2;
-        double celldm13 = (latvec.e11 * latvec.e31 + latvec.e12 * latvec.e32
-                           + latvec.e13 * latvec.e33);
-        double cos13 = celldm13 / celldm1 / celldm3;
-        double celldm23 = (latvec.e21 * latvec.e31 + latvec.e22 * latvec.e32
-                           + latvec.e23 * latvec.e33);
-        double cos23 = celldm23 / celldm2 / celldm3;
-
-        double sin12 = std::sqrt(1.0 - cos12 * cos12);
-        if (cos12 >= 1.0) {
-            ModuleBase::WARNING_QUIT("unitcell", "wrong cos12!");
-        }
-
-        latvec.e11 = celldm1;
-        latvec.e12 = 0.0;
-        latvec.e13 = 0.0;
-        latvec.e21 = celldm2 * cos12;
-        latvec.e22 = celldm2 * sin12;
-        latvec.e23 = 0.0;
-        latvec.e31 = celldm3 * cos13;
-        latvec.e32 = celldm3 * (cos23 - cos13 * cos12) / sin12;
-        double term = 1.0 + 2.0 * cos12 * cos13 * cos23 - cos12 * cos12
-                      - cos13 * cos13 - cos23 * cos23;
-        term = sqrt(term) / sin12;
-        latvec.e33 = celldm3 * term;
-    } else {
-        std::cout << "latname is : " << latName << std::endl;
-        ModuleBase::WARNING_QUIT("UnitCell::read_atom_species",
-                                 "latname not supported!");
-    }
-}
 
 void UnitCell::compare_atom_labels(std::string label1, std::string label2) {
     if (label1
diff --git a/source/module_cell/unitcell.h b/source/module_cell/unitcell.h
index d99eb91b5a..16c5ea2ec7 100644
--- a/source/module_cell/unitcell.h
+++ b/source/module_cell/unitcell.h
@@ -277,11 +277,8 @@ class UnitCell {
                          const bool& dpks_desc = false,
                          const int& iproc = 0) const;
     void check_dtau();
-    void setup_cell_after_vc(std::ofstream& log); // LiuXh add 20180515
-
     // for constrained vc-relaxation where type of lattice
     // is fixed, adjust the lattice vectors
-    void remake_cell();
 
     //================================================================
     // cal_natomwfc : calculate total number of atomic wavefunctions
diff --git a/source/module_cell/unitcell_data.h b/source/module_cell/unitcell_data.h
index 57f81bcd4e..09aab0a0dc 100644
--- a/source/module_cell/unitcell_data.h
+++ b/source/module_cell/unitcell_data.h
@@ -1,3 +1,6 @@
+#ifndef UNITCELL_DATA_H
+#define UNITCELL_DATA_H
+
 #include "module_base/intarray.h"
 #include "module_base/matrix3.h"
 /// @brief info of lattice
@@ -59,4 +62,6 @@ struct Statistics
         delete[] iwt2iat;
         delete[] iwt2iw;
     }
-};
\ No newline at end of file
+};
+
+#endif
\ No newline at end of file
diff --git a/source/module_cell/update_cell.cpp b/source/module_cell/update_cell.cpp
new file mode 100644
index 0000000000..6e989dd82e
--- /dev/null
+++ b/source/module_cell/update_cell.cpp
@@ -0,0 +1,347 @@
+#include "update_cell.h"
+#include "module_base/global_function.h"
+namespace unitcell
+{
+void remake_cell(Lattice& lat)
+{
+    ModuleBase::TITLE("Lattice", "rmake_cell");
+
+    // The idea is as follows: for each type of lattice, first calculate
+    // from current latvec the lattice parameters, then use the parameters
+    // to reconstruct latvec
+    std::string& latName = lat.latName;
+    ModuleBase::Matrix3&  latvec = lat.latvec;
+
+    if (latName == "none") {
+        ModuleBase::WARNING_QUIT("UnitCell",
+            "to use fixed_ibrav, latname must be provided");
+    } else if (latName == "sc") // ibrav = 1
+    {
+        double celldm = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
+                                  + pow(latvec.e13, 2));
+
+        latvec.Zero();
+        latvec.e11 = latvec.e22 = latvec.e33 = celldm;
+    } else if (latName == "fcc") // ibrav = 2
+    {
+        double celldm = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
+                                + pow(latvec.e13, 2)) / std::sqrt(2.0);
+
+        latvec.e11 = -celldm;
+        latvec.e12 = 0.0;
+        latvec.e13 = celldm;
+        latvec.e21 = 0.0;
+        latvec.e22 = celldm;
+        latvec.e23 = celldm;
+        latvec.e31 = -celldm;
+        latvec.e32 = celldm;
+        latvec.e33 = 0.0;
+    } else if (latName == "bcc") // ibrav = 3
+    {
+        double celldm = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
+                                  + pow(latvec.e13, 2))
+                        / std::sqrt(3.0);
+
+        latvec.e11 = celldm;
+        latvec.e12 = celldm;
+        latvec.e13 = celldm;
+        latvec.e21 = -celldm;
+        latvec.e22 = celldm;
+        latvec.e23 = celldm;
+        latvec.e31 = -celldm;
+        latvec.e32 = -celldm;
+        latvec.e33 = celldm;
+    } else if (latName == "hexagonal") // ibrav = 4
+    {
+        double celldm1 = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
+                                   + pow(latvec.e13, 2));
+        double celldm3 = std::sqrt(pow(latvec.e31, 2) + pow(latvec.e32, 2)
+                                   + pow(latvec.e33, 2));
+        double e22 = sqrt(3.0) / 2.0;
+
+        latvec.e11 = celldm1;
+        latvec.e12 = 0.0;
+        latvec.e13 = 0.0;
+        latvec.e21 = -0.5 * celldm1;
+        latvec.e22 = celldm1 * e22;
+        latvec.e23 = 0.0;
+        latvec.e31 = 0.0;
+        latvec.e32 = 0.0;
+        latvec.e33 = celldm3;
+    } else if (latName == "trigonal") // ibrav = 5
+    {
+        double celldm1 = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
+                                   + pow(latvec.e13, 2));
+        double celldm2 = std::sqrt(pow(latvec.e21, 2) + pow(latvec.e22, 2)
+                                   + pow(latvec.e23, 2));
+        double celldm12 = (latvec.e11 * latvec.e21 + latvec.e12 * latvec.e22
+                           + latvec.e13 * latvec.e23);
+        double cos12 = celldm12 / celldm1 / celldm2;
+
+        if (cos12 <= -0.5 || cos12 >= 1.0) {
+            ModuleBase::WARNING_QUIT("unitcell", "wrong cos12!");
+        }
+        double t1 = sqrt(1.0 + 2.0 * cos12);
+        double t2 = sqrt(1.0 - cos12);
+
+        double e11 = celldm1 * t2 / sqrt(2.0);
+        double e12 = -celldm1 * t2 / sqrt(6.0);
+        double e13 = celldm1 * t1 / sqrt(3.0);
+        double e22 = celldm1 * sqrt(2.0) * t2 / sqrt(3.0);
+
+        latvec.e11 = e11;
+        latvec.e12 = e12;
+        latvec.e13 = e13;
+        latvec.e21 = 0.0;
+        latvec.e22 = e22;
+        latvec.e23 = e13;
+        latvec.e31 = -e11;
+        latvec.e32 = e12;
+        latvec.e33 = e13;
+    } else if (latName == "st") // ibrav = 6
+    {
+        double celldm1 = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
+                                   + pow(latvec.e13, 2));
+        double celldm3 = std::sqrt(pow(latvec.e31, 2) + pow(latvec.e32, 2)
+                                   + pow(latvec.e33, 2));
+        latvec.e11 = celldm1;
+        latvec.e12 = 0.0;
+        latvec.e13 = 0.0;
+        latvec.e21 = 0.0;
+        latvec.e22 = celldm1;
+        latvec.e23 = 0.0;
+        latvec.e31 = 0.0;
+        latvec.e32 = 0.0;
+        latvec.e33 = celldm3;
+    } else if (latName == "bct") // ibrav = 7
+    {
+        double celldm1 = std::abs(latvec.e11);
+        double celldm2 = std::abs(latvec.e13);
+
+        latvec.e11 = celldm1;
+        latvec.e12 = -celldm1;
+        latvec.e13 = celldm2;
+        latvec.e21 = celldm1;
+        latvec.e22 = celldm1;
+        latvec.e23 = celldm2;
+        latvec.e31 = -celldm1;
+        latvec.e32 = -celldm1;
+        latvec.e33 = celldm2;
+    } else if (latName == "so") // ibrav = 8
+    {
+        double celldm1 = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
+                                   + pow(latvec.e13, 2));
+        double celldm2 = std::sqrt(pow(latvec.e21, 2) + pow(latvec.e22, 2)
+                                   + pow(latvec.e23, 2));
+        double celldm3 = std::sqrt(pow(latvec.e31, 2) + pow(latvec.e32, 2)
+                                   + pow(latvec.e33, 2));
+
+        latvec.e11 = celldm1;
+        latvec.e12 = 0.0;
+        latvec.e13 = 0.0;
+        latvec.e21 = 0.0;
+        latvec.e22 = celldm2;
+        latvec.e23 = 0.0;
+        latvec.e31 = 0.0;
+        latvec.e32 = 0.0;
+        latvec.e33 = celldm3;
+    } else if (latName == "baco") // ibrav = 9
+    {
+        double celldm1 = std::abs(latvec.e11);
+        double celldm2 = std::abs(latvec.e22);
+        double celldm3 = std::abs(latvec.e33);
+
+        latvec.e11 = celldm1;
+        latvec.e12 = celldm2;
+        latvec.e13 = 0.0;
+        latvec.e21 = -celldm1;
+        latvec.e22 = celldm2;
+        latvec.e23 = 0.0;
+        latvec.e31 = 0.0;
+        latvec.e32 = 0.0;
+        latvec.e33 = celldm3;
+    } else if (latName == "fco") // ibrav = 10
+    {
+        double celldm1 = std::abs(latvec.e11);
+        double celldm2 = std::abs(latvec.e22);
+        double celldm3 = std::abs(latvec.e33);
+
+        latvec.e11 = celldm1;
+        latvec.e12 = 0.0;
+        latvec.e13 = celldm3;
+        latvec.e21 = celldm1;
+        latvec.e22 = celldm2;
+        latvec.e23 = 0.0;
+        latvec.e31 = 0.0;
+        latvec.e32 = celldm2;
+        latvec.e33 = celldm3;
+    } else if (latName == "bco") // ibrav = 11
+    {
+        double celldm1 = std::abs(latvec.e11);
+        double celldm2 = std::abs(latvec.e12);
+        double celldm3 = std::abs(latvec.e13);
+
+        latvec.e11 = celldm1;
+        latvec.e12 = celldm2;
+        latvec.e13 = celldm3;
+        latvec.e21 = -celldm1;
+        latvec.e22 = celldm2;
+        latvec.e23 = celldm3;
+        latvec.e31 = -celldm1;
+        latvec.e32 = -celldm2;
+        latvec.e33 = celldm3;
+    } else if (latName == "sm") // ibrav = 12
+    {
+        double celldm1 = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
+                                   + pow(latvec.e13, 2));
+        double celldm2 = std::sqrt(pow(latvec.e21, 2) + pow(latvec.e22, 2)
+                                   + pow(latvec.e23, 2));
+        double celldm3 = std::sqrt(pow(latvec.e31, 2) + pow(latvec.e32, 2)
+                                   + pow(latvec.e33, 2));
+        double celldm12 = (latvec.e11 * latvec.e21 + latvec.e12 * latvec.e22
+                           + latvec.e13 * latvec.e23);
+        double cos12 = celldm12 / celldm1 / celldm2;
+
+        double e21 = celldm2 * cos12;
+        double e22 = celldm2 * std::sqrt(1.0 - cos12 * cos12);
+
+        latvec.e11 = celldm1;
+        latvec.e12 = 0.0;
+        latvec.e13 = 0.0;
+        latvec.e21 = e21;
+        latvec.e22 = e22;
+        latvec.e23 = 0.0;
+        latvec.e31 = 0.0;
+        latvec.e32 = 0.0;
+        latvec.e33 = celldm3;
+    } else if (latName == "bacm") // ibrav = 13
+    {
+        double celldm1 = std::abs(latvec.e11);
+        double celldm2 = std::sqrt(pow(latvec.e21, 2) + pow(latvec.e22, 2)
+                                   + pow(latvec.e23, 2));
+        double celldm3 = std::abs(latvec.e13);
+
+        double cos12 = latvec.e21 / celldm2;
+        if (cos12 >= 1.0) {
+            ModuleBase::WARNING_QUIT("unitcell", "wrong cos12!");
+        }
+
+        double e21 = celldm2 * cos12;
+        double e22 = celldm2 * std::sqrt(1.0 - cos12 * cos12);
+
+        latvec.e11 = celldm1;
+        latvec.e12 = 0.0;
+        latvec.e13 = -celldm3;
+        latvec.e21 = e21;
+        latvec.e22 = e22;
+        latvec.e23 = 0.0;
+        latvec.e31 = celldm1;
+        latvec.e32 = 0.0;
+        latvec.e33 = celldm3;
+    } else if (latName == "triclinic") // ibrav = 14
+    {
+        double celldm1 = std::sqrt(pow(latvec.e11, 2) + pow(latvec.e12, 2)
+                                   + pow(latvec.e13, 2));
+        double celldm2 = std::sqrt(pow(latvec.e21, 2) + pow(latvec.e22, 2)
+                                   + pow(latvec.e23, 2));
+        double celldm3 = std::sqrt(pow(latvec.e31, 2) + pow(latvec.e32, 2)
+                                   + pow(latvec.e33, 2));
+        double celldm12 = (latvec.e11 * latvec.e21 + latvec.e12 * latvec.e22
+                           + latvec.e13 * latvec.e23);
+        double cos12 = celldm12 / celldm1 / celldm2;
+        double celldm13 = (latvec.e11 * latvec.e31 + latvec.e12 * latvec.e32
+                           + latvec.e13 * latvec.e33);
+        double cos13 = celldm13 / celldm1 / celldm3;
+        double celldm23 = (latvec.e21 * latvec.e31 + latvec.e22 * latvec.e32
+                           + latvec.e23 * latvec.e33);
+        double cos23 = celldm23 / celldm2 / celldm3;
+
+        double sin12 = std::sqrt(1.0 - cos12 * cos12);
+        if (cos12 >= 1.0) {
+            ModuleBase::WARNING_QUIT("unitcell", "wrong cos12!");
+        }
+
+        latvec.e11 = celldm1;
+        latvec.e12 = 0.0;
+        latvec.e13 = 0.0;
+        latvec.e21 = celldm2 * cos12;
+        latvec.e22 = celldm2 * sin12;
+        latvec.e23 = 0.0;
+        latvec.e31 = celldm3 * cos13;
+        latvec.e32 = celldm3 * (cos23 - cos13 * cos12) / sin12;
+        double term = 1.0 + 2.0 * cos12 * cos13 * cos23 - cos12 * cos12
+                      - cos13 * cos13 - cos23 * cos23;
+        term = sqrt(term) / sin12;
+        latvec.e33 = celldm3 * term;
+    } else {
+        std::cout << "latname is : " << latName << std::endl;
+        ModuleBase::WARNING_QUIT("UnitCell::read_atom_species",
+                                 "latname not supported!");
+    }
+}
+
+// LiuXh add a new function here,
+// 20180515
+void setup_cell_after_vc(UnitCell& ucell, std::ofstream& log) {
+    ModuleBase::TITLE("UnitCell", "setup_cell_after_vc");
+    assert(ucell.lat0 > 0.0);
+    ucell.omega = std::abs(ucell.latvec.Det()) * 
+                           pow(ucell.lat0, 3);
+    if (ucell.omega <= 0) 
+    {
+        ModuleBase::WARNING_QUIT("setup_cell_after_vc", "omega <= 0 .");
+    } else {
+        log << std::endl;
+        ModuleBase::GlobalFunc::OUT(log, "Volume (Bohr^3)", ucell.omega);
+        ModuleBase::GlobalFunc::OUT(log, "Volume (A^3)",
+                                    ucell.omega * pow(ModuleBase::BOHR_TO_A, 3));
+    }
+
+    ucell.lat0_angstrom = ucell.lat0 * 0.529177;
+    ucell.tpiba = ModuleBase::TWO_PI / ucell.lat0;
+    ucell.tpiba2 = ucell.tpiba * ucell.tpiba;
+
+    // lattice vectors in another form.
+    ucell.a1.x = ucell.latvec.e11;
+    ucell.a1.y = ucell.latvec.e12;
+    ucell.a1.z = ucell.latvec.e13;
+
+    ucell.a2.x = ucell.latvec.e21;
+    ucell.a2.y = ucell.latvec.e22;
+    ucell.a2.z = ucell.latvec.e23;
+
+    ucell.a3.x = ucell.latvec.e31;
+    ucell.a3.y = ucell.latvec.e32;
+    ucell.a3.z = ucell.latvec.e33;
+
+    //==========================================================
+    // Calculate recip. lattice vectors and dot products
+    // latvec has the unit of lat0, but G has the unit 2Pi/lat0
+    //==========================================================
+    ucell.GT = ucell.latvec.Inverse();
+    ucell.G = ucell.GT.Transpose();
+    ucell.GGT = ucell.G * ucell.GT;
+    ucell.invGGT = ucell.GGT.Inverse();
+
+    for (int it = 0; it < ucell.ntype; it++) {
+        Atom* atom = &ucell.atoms[it];
+        for (int ia = 0; ia < atom->na; ia++) {
+            atom->tau[ia] = atom->taud[ia] * ucell.latvec;
+        }
+    }
+
+#ifdef __MPI
+    ucell.bcast_unitcell();
+#endif
+
+    log << std::endl;
+    output::printM3(log,
+                    "Lattice vectors: (Cartesian coordinate: in unit of a_0)",
+                    ucell.latvec);
+    output::printM3(log,
+                   "Reciprocal vectors: (Cartesian coordinate: in unit of 2 pi/a_0)",
+                    ucell.G);
+
+    return;
+}
+}
diff --git a/source/module_cell/update_cell.h b/source/module_cell/update_cell.h
new file mode 100644
index 0000000000..3f5a02c941
--- /dev/null
+++ b/source/module_cell/update_cell.h
@@ -0,0 +1,25 @@
+#ifndef UPDATE_CELL_H
+#define UPDATE_CELL_H
+
+#include "unitcell_data.h"
+#include "unitcell.h"  
+
+/*
+this file is used to update the cell,contains the following functions:
+1. remake_cell: for constrained vc-relaxation where type of lattice 
+is fixed, adjust the lattice vectors
+2. setup_cell_after_vc: setup cell after vc-relaxation
+the functions are defined in the namespace UnitCell,
+Accually, the functions are focused on the cell-relax part functions
+of the UnitCell class.
+*/
+namespace unitcell
+{
+    // for constrained vc-relaxation where type of lattice
+    // is fixed, adjust the lattice vectors
+    void remake_cell(Lattice& lat);
+
+    void setup_cell_after_vc(UnitCell& ucell, std::ofstream& log);
+}
+//
+#endif // UPDATE_CELL_H
\ No newline at end of file
diff --git a/source/module_md/msst.cpp b/source/module_md/msst.cpp
index 67484f0f88..faa9eada19 100644
--- a/source/module_md/msst.cpp
+++ b/source/module_md/msst.cpp
@@ -1,5 +1,6 @@
 #include "msst.h"
 
+#include "module_cell/update_cell.h"
 #include "md_func.h"
 #ifdef __MPI
 #include "mpi.h"
@@ -128,7 +129,7 @@ void MSST::first_half(std::ofstream& ofs)
     return;
 }
 
-void MSST::second_half(void)
+void MSST::second_half()
 {
     ModuleBase::TITLE("MSST", "second_half");
     ModuleBase::timer::tick("MSST", "second_half");
@@ -233,7 +234,7 @@ void MSST::restart(const std::string& global_readin_dir)
     return;
 }
 
-double MSST::vel_sum(void)
+double MSST::vel_sum()
 {
     double vsum = 0;
 
@@ -256,7 +257,7 @@ void MSST::rescale(std::ofstream& ofs, const double& volume)
     ucell.latvec.e22 *= dilation[1];
     ucell.latvec.e33 *= dilation[2];
 
-    ucell.setup_cell_after_vc(ofs);
+    unitcell::setup_cell_after_vc(ucell,ofs);
 
     /// rescale velocity
     for (int i = 0; i < ucell.nat; ++i)
@@ -266,7 +267,7 @@ void MSST::rescale(std::ofstream& ofs, const double& volume)
 }
 
 
-void MSST::propagate_vel(void)
+void MSST::propagate_vel()
 {
     if (my_rank == 0)
     {
@@ -306,7 +307,7 @@ void MSST::propagate_vel(void)
 }
 
 
-void MSST::propagate_voldot(void)
+void MSST::propagate_voldot()
 {
     const int sd = mdp.msst_direction;
     const double dthalf = 0.5 * md_dt;
diff --git a/source/module_md/nhchain.cpp b/source/module_md/nhchain.cpp
index 215e94327a..db4662f44c 100644
--- a/source/module_md/nhchain.cpp
+++ b/source/module_md/nhchain.cpp
@@ -5,7 +5,7 @@
 #include "mpi.h"
 #endif
 #include "module_base/timer.h"
-
+#include "module_cell/update_cell.h"
 Nose_Hoover::Nose_Hoover(const Parameter& param_in, UnitCell& unit_in) : MD_base(param_in, unit_in)
 {
     const double unit_transform = ModuleBase::HARTREE_SI / pow(ModuleBase::BOHR_RADIUS_SI, 3) * 1.0e-8;
@@ -272,7 +272,7 @@ void Nose_Hoover::first_half(std::ofstream& ofs)
 }
 
 
-void Nose_Hoover::second_half(void)
+void Nose_Hoover::second_half()
 {
     ModuleBase::TITLE("Nose_Hoover", "second_half");
     ModuleBase::timer::tick("Nose_Hoover", "second_half");
@@ -809,7 +809,7 @@ void Nose_Hoover::update_volume(std::ofstream& ofs)
     }
 
     /// reset ucell and pos due to change of lattice
-    ucell.setup_cell_after_vc(ofs);
+    unitcell::setup_cell_after_vc(ucell,ofs);
 }
 
 void Nose_Hoover::target_stress()
diff --git a/source/module_md/test/CMakeLists.txt b/source/module_md/test/CMakeLists.txt
index 1d197a5350..6510f0033d 100644
--- a/source/module_md/test/CMakeLists.txt
+++ b/source/module_md/test/CMakeLists.txt
@@ -87,6 +87,8 @@ AddTest(
   SOURCES nhchain_test.cpp
   ../md_base.cpp
   ../nhchain.cpp
+  ../../module_cell/update_cell.cpp
+  ../../module_io/output.cpp
   ${depend_files}
 )
 
@@ -96,6 +98,8 @@ AddTest(
   SOURCES msst_test.cpp
   ../md_base.cpp
   ../msst.cpp
+  ../../module_cell/update_cell.cpp
+  ../../module_io/output.cpp
   ${depend_files}
 )
 
diff --git a/source/module_relax/relax_new/relax.cpp b/source/module_relax/relax_new/relax.cpp
index e873a7aa50..88800c88d8 100644
--- a/source/module_relax/relax_new/relax.cpp
+++ b/source/module_relax/relax_new/relax.cpp
@@ -1,8 +1,10 @@
 #include "relax.h"
 
+
 #include "module_base/matrix3.h"
 #include "module_base/parallel_common.h"
 #include "module_base/tool_title.h"
+#include "module_cell/update_cell.h"
 #include "module_hamilt_pw/hamilt_pwdft/global.h"
 #include "module_parameter/parameter.h"
 #include "module_relax/relax_old/ions_move_basic.h"
@@ -583,7 +585,7 @@ void Relax::move_cell_ions(UnitCell& ucell, const bool is_new_dir)
         }
         if (PARAM.inp.fixed_ibrav)
         {
-            ucell.remake_cell();
+            unitcell::remake_cell(ucell.lat);
         }
     }
 
@@ -693,7 +695,7 @@ void Relax::move_cell_ions(UnitCell& ucell, const bool is_new_dir)
     // I do not want to change it
     if (if_cell_moves)
     {
-        ucell.setup_cell_after_vc(GlobalV::ofs_running);
+        unitcell::setup_cell_after_vc(ucell,GlobalV::ofs_running);
         ModuleBase::GlobalFunc::DONE(GlobalV::ofs_running, "SETUP UNITCELL");
     }
 }
\ No newline at end of file
diff --git a/source/module_relax/relax_new/test/CMakeLists.txt b/source/module_relax/relax_new/test/CMakeLists.txt
index 92fb60770c..8a2f5f7ff5 100644
--- a/source/module_relax/relax_new/test/CMakeLists.txt
+++ b/source/module_relax/relax_new/test/CMakeLists.txt
@@ -18,6 +18,7 @@ AddTest(
     ../../../module_base/matrix3.cpp ../../../module_base/intarray.cpp ../../../module_base/tool_title.cpp
     ../../../module_base/global_function.cpp ../../../module_base/complexmatrix.cpp ../../../module_base/matrix.cpp
     ../../../module_base/complexarray.cpp ../../../module_base/tool_quit.cpp ../../../module_base/realarray.cpp ../../../module_base/blas_connector.cpp
-  LIBS parameter ${math_libs} 
+    ../../../module_cell/update_cell.cpp  ../../../module_io/output.cpp
+    LIBS parameter ${math_libs} 
 )
 
diff --git a/source/module_relax/relax_new/test/relax_test.cpp b/source/module_relax/relax_new/test/relax_test.cpp
index a95516b210..74328026a4 100644
--- a/source/module_relax/relax_new/test/relax_test.cpp
+++ b/source/module_relax/relax_new/test/relax_test.cpp
@@ -72,6 +72,8 @@ class Test_SETGRAD : public testing::Test
             ucell.atoms[0].taud[1] = 0.0;
             ucell.atoms[0].taud[2] = 0.0;
 
+            ucell.atoms[0].tau.resize(nat);
+
             ucell.lc[0] = 1;
             ucell.lc[1] = 1;
             ucell.lc[2] = 1;
@@ -105,11 +107,13 @@ class Test_SETGRAD : public testing::Test
             ucell.latvec.Identity();
             input.fixed_axes = "a"; //anything other than "None"
             input.fixed_ibrav = true;
+            ucell.latName = "sc";
             ucell.lc[0] = 0;
             ucell.lc[1] = 0;
             ucell.lc[2] = 0;
             rl.init_relax(nat);
             rl.relax_step(ucell,force_in,stress_in,0.0);
+            
             push_result();
         }
 
@@ -167,7 +171,6 @@ class Test_RELAX : public testing::Test
             this->setup_cell();
 
             ModuleBase::matrix force_in, stress_in;
-
             force_in.create(nat,3);
             stress_in.create(3,3);
 
@@ -192,6 +195,7 @@ class Test_RELAX : public testing::Test
 
                 energy_file >> energy;
 
+                PARAM.input.fixed_ibrav = false;
                 rl.relax_step(ucell,force_in,stress_in,energy);
 
                 result.push_back(ucell.atoms[0].taud[0].x);
@@ -253,6 +257,7 @@ class Test_RELAX : public testing::Test
                 int na = ucell.atoms[i].na;
                 ucell.atoms[i].mbl.resize(na);
                 ucell.atoms[i].taud.resize(na);
+                ucell.atoms[i].tau.resize(na);
                 for (int j=0;j<na;j++)
                 {
                     ucell.atoms[i].mbl[j] = {1,1,1};
diff --git a/source/module_relax/relax_new/test/relax_test.h b/source/module_relax/relax_new/test/relax_test.h
index 9516cd2edc..d60d37adf9 100644
--- a/source/module_relax/relax_new/test/relax_test.h
+++ b/source/module_relax/relax_new/test/relax_test.h
@@ -10,8 +10,6 @@ namespace GlobalC
 UnitCell::UnitCell(){};
 UnitCell::~UnitCell(){};
 
-void UnitCell::remake_cell(){};
-
 void UnitCell::update_pos_taud(double* posd_in)
 {
     int iat = 0;
@@ -38,7 +36,6 @@ void UnitCell::print_stru_file(const std::string& fn,
                                const bool& dpks_desc,
                                const int& iproc) const {};
 void UnitCell::print_tau()const{};
-void UnitCell::setup_cell_after_vc(std::ofstream &log){};
 
 Magnetism::Magnetism(){};
 Magnetism::~Magnetism(){};
diff --git a/source/module_relax/relax_old/relax_old.cpp b/source/module_relax/relax_old/relax_old.cpp
index b18219c6a7..82c77b58e6 100644
--- a/source/module_relax/relax_old/relax_old.cpp
+++ b/source/module_relax/relax_old/relax_old.cpp
@@ -3,7 +3,7 @@
 #include "module_base/global_function.h"
 #include "module_base/global_variable.h"
 #include "module_parameter/parameter.h"
-
+#include "module_cell/update_cell.h"
 void Relax_old::init_relax(const int& natom)
 {
     // Geometry optimization algorithm setup.
@@ -70,7 +70,7 @@ bool Relax_old::relax_step(const int& istep,
             force_step = 1;
             stress_step++;
             ucell.cell_parameter_updated = true;
-            ucell.setup_cell_after_vc(GlobalV::ofs_running);
+            unitcell::setup_cell_after_vc(ucell,GlobalV::ofs_running);
             ModuleBase::GlobalFunc::DONE(GlobalV::ofs_running, "SETUP UNITCELL");
         }
         return converged;