diff --git a/ProcessLib/GroundwaterFlowFEM.h b/ProcessLib/GroundwaterFlowFEM.h
index b2cca22b1a84a39212c615df0f07239d3ef7f88b..748f13f728655a1ea9a445e0e7b34d1c0fe2207b 100644
--- a/ProcessLib/GroundwaterFlowFEM.h
+++ b/ProcessLib/GroundwaterFlowFEM.h
@@ -73,9 +73,10 @@ public:
         IntegrationMethod_ integration_method(_integration_order);
         std::size_t const n_integration_points = integration_method.getNPoints();
 
-        _shape_matrices.resize(n_integration_points);
+        _shape_matrices.reserve(n_integration_points);
         for (std::size_t ip(0); ip < n_integration_points; ip++) {
-            _shape_matrices[ip].resize(ShapeFunction::DIM, ShapeFunction::NPOINTS);
+            _shape_matrices.emplace_back(ShapeFunction::DIM, GlobalDim,
+                                         ShapeFunction::NPOINTS);
             fe.computeShapeFunctions(
                     integration_method.getWeightedPoint(ip).getCoords(),
                     _shape_matrices[ip]);
diff --git a/ProcessLib/NeumannBcAssembler.h b/ProcessLib/NeumannBcAssembler.h
index e4fd38aa35410dc821e3c97d0d000700785cc208..064288059a4f8390b7134f85f8805d8a305b8f06 100644
--- a/ProcessLib/NeumannBcAssembler.h
+++ b/ProcessLib/NeumannBcAssembler.h
@@ -70,9 +70,10 @@ public:
         IntegrationMethod_ integration_method(_integration_order);
         std::size_t const n_integration_points = integration_method.getNPoints();
 
-        _shape_matrices.resize(n_integration_points);
+        _shape_matrices.reserve(n_integration_points);
         for (std::size_t ip(0); ip < n_integration_points; ip++) {
-            _shape_matrices[ip].resize(ShapeFunction::DIM, GlobalDim, ShapeFunction::NPOINTS);
+            _shape_matrices.emplace_back(ShapeFunction::DIM, GlobalDim,
+                                         ShapeFunction::NPOINTS);
             fe.computeShapeFunctions(
                     integration_method.getWeightedPoint(ip).getCoords(),
                     _shape_matrices[ip]);