diff --git a/src/environment/GPESolver.cpp b/src/environment/GPESolver.cpp
index 5e214b8ba..a1b6151ad 100644
--- a/src/environment/GPESolver.cpp
+++ b/src/environment/GPESolver.cpp
@@ -123,7 +123,7 @@ mrcpp::CompFunction<3> GPESolver::solvePoissonEquation(const mrcpp::CompFunction
     mrcpp::CompFunction<3> first_term;
     mrcpp::CompFunction<3> Vr_np1;
     Vr_np1.func_ptr->isreal = 1;
-    Vr_np1.alloc(0);
+    Vr_np1.alloc(1);
 
     auto eps_inv_func = mrcpp::AnalyticFunction<3>([this](const mrcpp::Coord<3> &r) { return 1.0 / this->epsilon.evalf(r); });
     Density rho_tot(false);
@@ -238,7 +238,7 @@ mrcpp::CompFunction<3> &GPESolver::solveEquation(double prec, const Density &rho
     Timer t_vac;
     mrcpp::CompFunction<3> V_vac;
     V_vac.func_ptr->isreal = 1;
-    V_vac.alloc(0);
+    V_vac.alloc(1);
     mrcpp::apply(this->apply_prec, V_vac.real(), *poisson, rho_tot.real());
     rho_tot.free();
     print_utils::qmfunction(3, "Vacuum potential", V_vac, t_vac);
@@ -250,7 +250,7 @@ mrcpp::CompFunction<3> &GPESolver::solveEquation(double prec, const Density &rho
          mrcpp::CompFunction<3> gamma_0;
          mrcpp::CompFunction<3> V_tot;
          gamma_0.func_ptr->isreal = 1;
-         gamma_0.alloc(0);
+         gamma_0.alloc(1);
 
          computeGamma(V_vac, gamma_0);
          this->Vr_n = solvePoissonEquation(gamma_0, rho_el);
@@ -273,7 +273,7 @@ auto GPESolver::computeEnergies(const Density &rho_el) -> std::tuple<double, dou
 void GPESolver::resetComplexFunction(mrcpp::CompFunction<3> &function) {
     function.func_ptr->isreal = 1;
     function.func_ptr->iscomplex = 0;
-    function.alloc(0);
+    function.alloc(1);
 }
 
 void GPESolver::printParameters() const {
diff --git a/src/environment/PBESolver.cpp b/src/environment/PBESolver.cpp
index 5ea02e2ae..39944dad2 100644
--- a/src/environment/PBESolver.cpp
+++ b/src/environment/PBESolver.cpp
@@ -63,7 +63,7 @@ void PBESolver::computePBTerm(mrcpp::CompFunction<3> &V_tot, const double salt_f
     resetComplexFunction(pb_term);
     mrcpp::CompFunction<3> sinhV;
     sinhV.func_ptr->isreal = 1;
-    sinhV.alloc(0);
+    sinhV.alloc(1);
     mrcpp::map(this->apply_prec / 100, sinhV.real(), V_tot.real(), sinh_f);
 
     mrcpp::multiply(pb_term, sinhV, this->kappa, this->apply_prec);
diff --git a/src/initial_guess/cube.cpp b/src/initial_guess/cube.cpp
index a199df375..ebc5693ce 100644
--- a/src/initial_guess/cube.cpp
+++ b/src/initial_guess/cube.cpp
@@ -120,7 +120,7 @@ bool initial_guess::cube::project_mo(OrbitalVector &Phi, double prec, const std:
         Timer t_i;
         if (mrcpp::mpi::my_func(Phi[i])) {
             CUBEfunction phi_i = CUBEVector[i];
-            Phi[i].alloc(0);
+            Phi[i].alloc(1);
             mrcpp::project(prec, Phi[i].real(), phi_i);
             std::stringstream o_txt;
             o_txt << std::setw(w1 - 1) << i;
diff --git a/src/initial_guess/sad.cpp b/src/initial_guess/sad.cpp
index 779361f9e..489350158 100644
--- a/src/initial_guess/sad.cpp
+++ b/src/initial_guess/sad.cpp
@@ -249,7 +249,7 @@ void initial_guess::sad::project_atomic_densities(double prec, Density &rho_tot,
 
     Timer t_tot;
     Density rho_loc(false);
-    rho_loc.alloc();
+    rho_loc.alloc(1);
     rho_loc.real().setZero();
 
     Timer t_loc;
diff --git a/src/qmfunctions/Density.cpp b/src/qmfunctions/Density.cpp
index f00aa11b2..843e6f50e 100644
--- a/src/qmfunctions/Density.cpp
+++ b/src/qmfunctions/Density.cpp
@@ -111,7 +111,7 @@ void Density::loadDensity(const std::string &file) {
     if (isreal()) {
         std::stringstream fname;
         fname << file << "_re";
-        alloc(0);
+        alloc(1);
         CompD[0]->loadTree(fname.str());
     }
 
@@ -119,7 +119,7 @@ void Density::loadDensity(const std::string &file) {
     if (iscomplex()) {
         std::stringstream fname;
         fname << file << "_cx";
-        alloc(0);
+        alloc(1);
         CompC[0]->loadTree(fname.str());
     }
 }
diff --git a/src/qmfunctions/density_utils.cpp b/src/qmfunctions/density_utils.cpp
index 2613e0017..dd7947dd0 100644
--- a/src/qmfunctions/density_utils.cpp
+++ b/src/qmfunctions/density_utils.cpp
@@ -154,7 +154,7 @@ void density::compute_local_X(double prec, Density &rho, OrbitalVector &Phi, Orb
     double add_prec = prec / N_el; // prec for rho = sum_i rho_i
     if (Phi.size() != X.size()) MSG_ERROR("Size mismatch");
 
-    if (rho.Ncomp() == 0) rho.alloc(0);
+    if (rho.Ncomp() == 0) rho.alloc(1);
 
     // Compute local density from own orbitals
     for (int i = 0; i < Phi.size(); i++) {
@@ -179,7 +179,7 @@ void density::compute_local_XY(double prec, Density &rho, OrbitalVector &Phi, Or
     if (Phi.size() != X.size()) MSG_ERROR("Size mismatch");
     if (Phi.size() != Y.size()) MSG_ERROR("Size mismatch");
 
-    if (rho.Ncomp() == 0) rho.alloc(0);
+    if (rho.Ncomp() == 0) rho.alloc(1);
 
     // Compute local density from own orbitals
     rho.real().setZero();
@@ -205,7 +205,7 @@ void density::compute_local_XY(double prec, Density &rho, OrbitalVector &Phi, Or
 }
 
 void density::compute(double prec, Density &rho, mrcpp::GaussExp<3> &dens_exp) {
-    if (not rho.hasReal()) rho.alloc(NUMBER::Real);
+    if (not rho.hasReal()) rho.alloc(1);
     mrcpp::build_grid(rho.real(), dens_exp);
     mrcpp::project(prec, rho.real(), dens_exp);
 }
@@ -228,7 +228,7 @@ void density::allreduce_density(double prec, Density &rho_tot, Density &rho_loc)
         if (mrcpp::mpi::numerically_exact) rho_loc.crop(prec);
     }
 
-    if (not rho_tot.hasReal()) rho_tot.alloc(NUMBER::Real);
+    if (not rho_tot.hasReal()) rho_tot.alloc(1);
 
     if (rho_tot.isShared()) {
         int tag = 2002;
diff --git a/src/qmfunctions/orbital_utils.cpp b/src/qmfunctions/orbital_utils.cpp
index 4e7f418ac..d4bfcfdce 100644
--- a/src/qmfunctions/orbital_utils.cpp
+++ b/src/qmfunctions/orbital_utils.cpp
@@ -1236,7 +1236,7 @@ void orbital::loadOrbital(const std::string &file, Orbital &orb) {
     if (orb.isreal()) {
         std::stringstream fname;
         fname << file << "_real";
-        orb.alloc(0);
+        orb.alloc(1);
         orb.CompD[0]->loadTree(fname.str());
     }
 
@@ -1244,7 +1244,7 @@ void orbital::loadOrbital(const std::string &file, Orbital &orb) {
     if (orb.iscomplex()) {
         std::stringstream fname;
         fname << file << "_complex";
-        orb.alloc(0);
+        orb.alloc(1);
         orb.CompC[0]->loadTree(fname.str());
     }
     delete mra;
diff --git a/src/qmoperators/QMDerivative.cpp b/src/qmoperators/QMDerivative.cpp
index f37521ee9..3c9046205 100644
--- a/src/qmoperators/QMDerivative.cpp
+++ b/src/qmoperators/QMDerivative.cpp
@@ -97,13 +97,13 @@ QMOperatorVector QMDerivative::apply(QMOperator_p &O) {
         // does not compile?       mrcpp::apply(V_out_base, D, V_inp, d);
         if (this->isReal()) {
             if (V_inp->isreal()) {
-                V_out->alloc(0);
+                V_out->alloc(1);
                 mrcpp::apply(V_out->real(), D, V_inp->real(), d);
             }
             if (V_inp->iscomplex()) {
                 V_out->func_ptr->isreal = 0;
                 V_out->func_ptr->iscomplex = 1;
-                V_out->alloc(0);
+                V_out->alloc(1);
                 mrcpp::apply(V_out->complex(), D, V_inp->complex(), d);
                 //if (V_inp->conjugate()) V_out->imag().rescale(-1.0);
             }
@@ -111,7 +111,7 @@ QMOperatorVector QMDerivative::apply(QMOperator_p &O) {
             if (V_inp->iscomplex()) {
                 V_out->func_ptr->isreal = 0;
                 V_out->func_ptr->iscomplex = 1;
-                V_out->alloc(0);
+                V_out->alloc(1);
                mrcpp::apply(V_out->complex(), D, V_inp->complex(), d);
                if (!V_inp->conjugate()) V_out->CompC[0]->rescale({0.0, 1.0});
             }
diff --git a/src/qmoperators/one_electron/DistanceOperator.h b/src/qmoperators/one_electron/DistanceOperator.h
index 8988a1345..aaa4a0e2a 100644
--- a/src/qmoperators/one_electron/DistanceOperator.h
+++ b/src/qmoperators/one_electron/DistanceOperator.h
@@ -57,7 +57,7 @@ class DistanceOperator final : public RankZeroOperator {
 
         // Project analytic potential, building grid for 1/r
         auto r_pow = std::make_shared<QMPotential>(1);
-        r_pow->alloc(0);
+        r_pow->alloc(1);
         mrcpp::build_grid(r_pow->real(), nuc_func);
         mrcpp::project<3>(proj_prec, r_pow->real(), f);
 
diff --git a/src/qmoperators/one_electron/NuclearOperator.cpp b/src/qmoperators/one_electron/NuclearOperator.cpp
index bbb259eb0..d2a9646f8 100644
--- a/src/qmoperators/one_electron/NuclearOperator.cpp
+++ b/src/qmoperators/one_electron/NuclearOperator.cpp
@@ -170,7 +170,7 @@ void NuclearOperator::allreducePotential(double prec, mrcpp::CompFunction<3> &V_
         if (mrcpp::mpi::numerically_exact) V_loc.crop(prec);
     }
 
-    if (not V_tot.hasReal()) V_tot.alloc(0);
+    if (not V_tot.hasReal()) V_tot.alloc(1);
     if (V_tot.isShared()) {
         int tag = 3141;
         // MPI grand master distributes to shared masters
diff --git a/src/qmoperators/two_electron/CoulombPotential.cpp b/src/qmoperators/two_electron/CoulombPotential.cpp
index 89cc0bfb9..4195634f1 100644
--- a/src/qmoperators/two_electron/CoulombPotential.cpp
+++ b/src/qmoperators/two_electron/CoulombPotential.cpp
@@ -126,7 +126,7 @@ void CoulombPotential::setupGlobalPotential(double prec) {
     bool need_to_apply = not(V.isShared()) or mrcpp::mpi::share_master();
 
     Timer timer;
-    V.alloc(0);
+    V.alloc(1);
     if (need_to_apply) mrcpp::apply(abs_prec, V.real(), P, rho.real());
     mrcpp::mpi::share_function(V, 0, 22445, mrcpp::mpi::comm_share);
     print_utils::qmfunction(3, "Compute global potential", V, timer);
@@ -151,7 +151,7 @@ mrcpp::CompFunction<3> CoulombPotential::setupLocalPotential(double prec) {
 
     Timer timer;
     mrcpp::CompFunction<3> V(false);
-    V.alloc(0);
+    V.alloc(1);
     mrcpp::apply(abs_prec, V.real(), P, rho.real());
     print_utils::qmfunction(3, "Compute local potential", V, timer);
 
@@ -169,7 +169,7 @@ void CoulombPotential::allreducePotential(double prec, mrcpp::CompFunction<3> &V
     // Add up local contributions into the grand master
     mrcpp::mpi::reduce_function(abs_prec, V_loc, mrcpp::mpi::comm_wrk);
 
-    if (not V_tot.hasReal()) V_tot.alloc(0);
+    if (not V_tot.hasReal()) V_tot.alloc(1);
     if (V_tot.isShared()) {
         int tag = 3141;
         // MPI grand master distributes to shared masters
diff --git a/src/qmoperators/two_electron/XCPotentialD1.cpp b/src/qmoperators/two_electron/XCPotentialD1.cpp
index 14c38b7ee..6313e0206 100644
--- a/src/qmoperators/two_electron/XCPotentialD1.cpp
+++ b/src/qmoperators/two_electron/XCPotentialD1.cpp
@@ -54,7 +54,7 @@ mrcpp::FunctionTreeVector<3> XCPotentialD1::setupDensities(double prec, mrcpp::F
             Timer timer;
             Density &rho = getDensity(DensityType::Total, 0);
             if (rho.Ncomp() == 0) {
-                rho.alloc(0);
+                rho.alloc(1);
                 mrcpp::copy_grid(rho.real(), grid);
                 density::compute(prec, rho, *orbitals, DensityType::Total);
             }
@@ -66,7 +66,7 @@ mrcpp::FunctionTreeVector<3> XCPotentialD1::setupDensities(double prec, mrcpp::F
             Timer timer;
             Density &rho = getDensity(DensityType::Alpha, 0);
             if (rho.Ncomp() == 0) {
-                rho.alloc(0);
+                rho.alloc(1);
                 mrcpp::copy_grid(rho.real(), grid);
                 density::compute(prec, rho, *orbitals, DensityType::Alpha);
             }
@@ -77,7 +77,7 @@ mrcpp::FunctionTreeVector<3> XCPotentialD1::setupDensities(double prec, mrcpp::F
             Timer timer;
             Density &rho = getDensity(DensityType::Beta, 0);
             if (rho.Ncomp() == 0) {
-                rho.alloc(0);
+                rho.alloc(1);
                 mrcpp::copy_grid(rho.real(), grid);
                 density::compute(prec, rho, *orbitals, DensityType::Beta);
             }
diff --git a/src/qmoperators/two_electron/XCPotentialD2.cpp b/src/qmoperators/two_electron/XCPotentialD2.cpp
index 94b279c18..f45beb154 100644
--- a/src/qmoperators/two_electron/XCPotentialD2.cpp
+++ b/src/qmoperators/two_electron/XCPotentialD2.cpp
@@ -75,7 +75,7 @@ mrcpp::FunctionTreeVector<3> XCPotentialD2::setupDensities(double prec, mrcpp::F
             Timer timer;
             Density &rho = getDensity(DensityType::Total, 0);
             if (rho.Ncomp() == 0) {
-                rho.alloc(0);
+                rho.alloc(1);
                 mrcpp::copy_grid(rho.real(), grid);
                 density::compute(prec, rho, *orbitals, DensityType::Total);
             }
@@ -86,7 +86,7 @@ mrcpp::FunctionTreeVector<3> XCPotentialD2::setupDensities(double prec, mrcpp::F
             Timer timer;
             Density &rho = getDensity(DensityType::Total, 1);
             if (rho.Ncomp() == 0) {
-                rho.alloc(0);
+                rho.alloc(1);
                 mrcpp::copy_grid(rho.real(), grid);
                 density::compute(prec, rho, *orbitals, *orbitals_x, *orbitals_y, DensityType::Total);
             }
@@ -98,7 +98,7 @@ mrcpp::FunctionTreeVector<3> XCPotentialD2::setupDensities(double prec, mrcpp::F
             Timer timer;
             Density &rho = getDensity(DensityType::Alpha, 0);
             if (rho.Ncomp() == 0) {
-                rho.alloc(0);
+                rho.alloc(1);
                 mrcpp::copy_grid(rho.real(), grid);
                 density::compute(prec, rho, *orbitals, DensityType::Alpha);
             }
@@ -109,7 +109,7 @@ mrcpp::FunctionTreeVector<3> XCPotentialD2::setupDensities(double prec, mrcpp::F
             Timer timer;
             Density &rho = getDensity(DensityType::Beta, 0);
             if (rho.Ncomp() == 0) {
-                rho.alloc(0);
+                rho.alloc(1);
                 mrcpp::copy_grid(rho.real(), grid);
                 density::compute(prec, rho, *orbitals, DensityType::Beta);
             }
@@ -120,7 +120,7 @@ mrcpp::FunctionTreeVector<3> XCPotentialD2::setupDensities(double prec, mrcpp::F
             Timer timer;
             Density &rho = getDensity(DensityType::Alpha, 1);
             if (rho.Ncomp() == 0) {
-                rho.alloc(0);
+                rho.alloc(1);
                 mrcpp::copy_grid(rho.real(), grid);
                 density::compute(prec, rho, *orbitals, *orbitals_x, *orbitals_y, DensityType::Alpha);
             }
@@ -131,7 +131,7 @@ mrcpp::FunctionTreeVector<3> XCPotentialD2::setupDensities(double prec, mrcpp::F
             Timer timer;
             Density &rho = getDensity(DensityType::Beta, 1);
             if (rho.Ncomp() == 0) {
-                rho.alloc(0);
+                rho.alloc(1);
                 mrcpp::copy_grid(rho.real(), grid);
                 density::compute(prec, rho, *orbitals, *orbitals_x, *orbitals_y, DensityType::Beta);
             }
diff --git a/tests/qmoperators/position_operator.cpp b/tests/qmoperators/position_operator.cpp
index 18197306b..47d695684 100644
--- a/tests/qmoperators/position_operator.cpp
+++ b/tests/qmoperators/position_operator.cpp
@@ -63,7 +63,7 @@ TEST_CASE("PositionOperator", "[position_operator]") {
     r.setup(prec);
     SECTION("vector apply") {
         OrbitalVector xPhi = r[0](Phi);
-        ComplexMatrix X = orbital::calc_overlap_matrix(Phi, xPhi);
+        ComplexMatrix X = mrcpp::calc_overlap_matrix(Phi, xPhi);
         for (int i = 0; i < X.rows(); i++) {
             for (int j = 0; j < X.cols(); j++) { REQUIRE(std::abs(X(i, j).real() - ref(i, j)) < thrs); }
         }