diff --git a/src/colvaratoms.cpp b/src/colvaratoms.cpp
index 3818818e5..6e09845a2 100644
--- a/src/colvaratoms.cpp
+++ b/src/colvaratoms.cpp
@@ -1562,32 +1562,30 @@ cvm::ag_vector_real_t cvm::atom_group::positions() const
   return atoms_pos;
 }
 
-cvm::ag_vector_real_t cvm::atom_group::positions_shifted(cvm::rvector const &shift) const
+int cvm::atom_group::positions_shifted(cvm::rvector const &shift, cvm::ag_vector_real_t& out_soa) const
 {
   if (b_dummy) {
-    cvmodule->error("Error: positions are not available "
+    return cvmodule->error("Error: positions are not available "
                "from a dummy atom group.\n", COLVARS_INPUT_ERROR);
   }
 
   if (is_enabled(f_ag_scalable)) {
-    cvmodule->error("Error: atomic positions are not available "
+    return cvmodule->error("Error: atomic positions are not available "
                "from a scalable atom group.\n", COLVARS_INPUT_ERROR);
   }
 
-  // std::vector<cvm::atom_pos> x(this->size(), 0.0);
-  // cvmodule->atom_const_iter ai = this->begin();
-  // std::vector<cvm::atom_pos>::iterator xi = x.begin();
-  // for ( ; ai != this->end(); ++xi, ++ai) {
-  //   *xi = (ai->pos + shift);
-  // }
-  // return x;
-  cvm::ag_vector_real_t shifted = atoms_pos;
+  if (out_soa.size() != 3 * num_atoms) {
+    out_soa.resize(3 * num_atoms);
+  }
+  auto* out_x = out_soa.data();
+  auto* out_y = out_x + num_atoms;
+  auto* out_z = out_y + num_atoms;
   for (size_t i = 0; i < num_atoms; ++i) {
-    shifted[i]             += shift.x;
-    shifted[i+num_atoms]   += shift.y;
-    shifted[i+2*num_atoms] += shift.z;
+    out_x[i] = pos_x(i) + shift.x;
+    out_y[i] = pos_y(i) + shift.y;
+    out_z[i] = pos_z(i) + shift.z;
   }
-  return shifted;
+  return COLVARS_OK;
 }
 
 std::vector<cvm::real> cvm::atom_group::velocities() const {
diff --git a/src/colvaratoms.h b/src/colvaratoms.h
index 717d9a281..a95469408 100644
--- a/src/colvaratoms.h
+++ b/src/colvaratoms.h
@@ -404,7 +404,7 @@ class cvm::atom_group: public colvardeps {
   /**
    * @brief Return a copy of the current atom positions, shifted by a constant vector
    */
-  cvm::ag_vector_real_t positions_shifted(cvm::rvector const &shift) const;
+  int positions_shifted(cvm::rvector const &shift, cvm::ag_vector_real_t& out) const;
   /**
    * @brief Return a copy of the current atom velocities
    */
diff --git a/src/colvarcomp_rotations.cpp b/src/colvarcomp_rotations.cpp
index 7c2b5b2b6..dc42fd4c1 100644
--- a/src/colvarcomp_rotations.cpp
+++ b/src/colvarcomp_rotations.cpp
@@ -117,7 +117,7 @@ int colvar::orientation::init(std::string const &conf)
 void colvar::orientation::calc_value()
 {
   atoms_cog = atoms->center_of_geometry();
-  shifted_pos_soa = atoms->positions_shifted(-1.0 * atoms_cog);
+  (void)atoms->positions_shifted(-1.0 * atoms_cog, shifted_pos_soa);
   rot.calc_optimal_rotation_soa(ref_pos_soa, shifted_pos_soa, num_ref_pos, atoms->size());
 
   if ((rot.q).inner(ref_quat) >= 0.0) {
@@ -190,7 +190,7 @@ colvar::orientation_angle::orientation_angle()
 void colvar::orientation_angle::calc_value()
 {
   atoms_cog = atoms->center_of_geometry();
-  shifted_pos_soa = atoms->positions_shifted(-1.0 * atoms_cog);
+  (void)atoms->positions_shifted(-1.0 * atoms_cog, shifted_pos_soa);
   rot.calc_optimal_rotation_soa(ref_pos_soa, shifted_pos_soa, num_ref_pos, atoms->size());
 
   if ((rot.q).q0 >= 0.0) {
@@ -262,7 +262,7 @@ colvar::orientation_proj::orientation_proj()
 void colvar::orientation_proj::calc_value()
 {
   atoms_cog = atoms->center_of_geometry();
-  shifted_pos_soa = atoms->positions_shifted(-1.0 * atoms_cog);
+  (void)atoms->positions_shifted(-1.0 * atoms_cog, shifted_pos_soa);
   rot.calc_optimal_rotation_soa(ref_pos_soa, shifted_pos_soa, num_ref_pos, atoms->size());
   x.real_value = 2.0 * (rot.q).q0 * (rot.q).q0 - 1.0;
 }
@@ -310,7 +310,7 @@ void colvar::tilt::calc_value()
 {
   atoms_cog = atoms->center_of_geometry();
 
-  shifted_pos_soa = atoms->positions_shifted(-1.0 * atoms_cog);
+  (void)atoms->positions_shifted(-1.0 * atoms_cog, shifted_pos_soa);
   rot.calc_optimal_rotation_soa(ref_pos_soa, shifted_pos_soa, num_ref_pos, atoms->size());
 
   x.real_value = rot.cos_theta(axis);
@@ -345,7 +345,7 @@ void colvar::spin_angle::calc_value()
 {
   atoms_cog = atoms->center_of_geometry();
 
-  shifted_pos_soa = atoms->positions_shifted(-1.0 * atoms_cog);
+  (void)atoms->positions_shifted(-1.0 * atoms_cog, shifted_pos_soa);
   rot.calc_optimal_rotation_soa(ref_pos_soa, shifted_pos_soa, num_ref_pos, atoms->size());
 
   x.real_value = rot.spin_angle(axis);
@@ -380,7 +380,7 @@ colvar::euler_phi::euler_phi()
 void colvar::euler_phi::calc_value()
 {
   atoms_cog = atoms->center_of_geometry();
-  shifted_pos_soa = atoms->positions_shifted(-1.0 * atoms_cog);
+  (void)atoms->positions_shifted(-1.0 * atoms_cog, shifted_pos_soa);
   rot.calc_optimal_rotation_soa(ref_pos_soa, shifted_pos_soa, num_ref_pos, atoms->size());
 
   const cvm::real& q0 = rot.q.q0;
@@ -427,7 +427,7 @@ colvar::euler_psi::euler_psi()
 void colvar::euler_psi::calc_value()
 {
   atoms_cog = atoms->center_of_geometry();
-  shifted_pos_soa = atoms->positions_shifted(-1.0 * atoms_cog);
+  (void)atoms->positions_shifted(-1.0 * atoms_cog, shifted_pos_soa);
   rot.calc_optimal_rotation_soa(ref_pos_soa, shifted_pos_soa, num_ref_pos, atoms->size());
 
   const cvm::real& q0 = rot.q.q0;
@@ -474,7 +474,7 @@ colvar::euler_theta::euler_theta()
 void colvar::euler_theta::calc_value()
 {
   atoms_cog = atoms->center_of_geometry();
-  shifted_pos_soa = atoms->positions_shifted(-1.0 * atoms_cog);
+  (void)atoms->positions_shifted(-1.0 * atoms_cog, shifted_pos_soa);
   rot.calc_optimal_rotation_soa(ref_pos_soa, shifted_pos_soa, num_ref_pos, atoms->size());
 
   const cvm::real& q0 = rot.q.q0;