diff --git a/cashocs/_forms/control_form_handler.py b/cashocs/_forms/control_form_handler.py
index bce5ce60..b8aeff47 100644
--- a/cashocs/_forms/control_form_handler.py
+++ b/cashocs/_forms/control_form_handler.py
@@ -173,7 +173,7 @@ def scalar_product(
             x = fenics.as_backend_type(a[i].vector()).vec()
             y = fenics.as_backend_type(b[i].vector()).vec()
 
-            temp, _ = self.riesz_projection_matrices[i].getVecs()
+            temp = self.riesz_projection_matrices[i].createVecRight()
             self.riesz_projection_matrices[i].mult(x, temp)
             result += temp.dot(y)
 
diff --git a/cashocs/_forms/shape_form_handler.py b/cashocs/_forms/shape_form_handler.py
index 34657c30..1f4b8b3e 100644
--- a/cashocs/_forms/shape_form_handler.py
+++ b/cashocs/_forms/shape_form_handler.py
@@ -686,7 +686,7 @@ def scalar_product(
             x = fenics.as_backend_type(a[0].vector()).vec()
             y = fenics.as_backend_type(b[0].vector()).vec()
 
-            temp, _ = self.scalar_product_matrix.getVecs()
+            temp = self.scalar_product_matrix.createVecRight()
             self.scalar_product_matrix.mult(x, temp)
             result = temp.dot(y)
 
diff --git a/cashocs/_optimization/line_search/line_search.py b/cashocs/_optimization/line_search/line_search.py
index 213376ed..c9f3fabe 100644
--- a/cashocs/_optimization/line_search/line_search.py
+++ b/cashocs/_optimization/line_search/line_search.py
@@ -140,9 +140,9 @@ def perform(
 
             transfer_matrix = cast(PETSc.Mat, transfer_matrix)
 
-            _, temp = transfer_matrix.getVecs()
-            transfer_matrix.mult(x, temp)
-            self.global_deformation_vector.axpy(1.0, temp)
+            transfer_matrix.multAdd(
+                x, self.global_deformation_vector, self.global_deformation_vector
+            )
             self.deformation_function.vector().apply("")
 
         self.post_line_search()
diff --git a/cashocs/_utils/linalg.py b/cashocs/_utils/linalg.py
index 5f4f37fe..a37304dd 100644
--- a/cashocs/_utils/linalg.py
+++ b/cashocs/_utils/linalg.py
@@ -518,10 +518,10 @@ def solve(
         A = setup_matrix_and_preconditioner(ksp, A, P)
 
         if b is None:
-            return A.getVecs()[0]
+            return A.createVecRight()
 
         if fun is None:
-            x, _ = A.getVecs()
+            x = A.createVecRight()
         else:
             x = fun.vector().vec()
 
@@ -635,9 +635,7 @@ def interpolate(self, u: fenics.Function) -> fenics.Function:
         """
         v = fenics.Function(self.target_space)
         x = fenics.as_backend_type(u.vector()).vec()
-        _, temp = self.transfer_matrix.getVecs()
-        self.transfer_matrix.mult(x, temp)
-        v.vector().vec().aypx(0.0, temp)
+        self.transfer_matrix.mult(x, v.vector().vec())
         v.vector().apply("")
 
         return v
diff --git a/docs/source/_static/version_switcher.json b/docs/source/_static/version_switcher.json
index 2765dc44..de86e2a2 100755
--- a/docs/source/_static/version_switcher.json
+++ b/docs/source/_static/version_switcher.json
@@ -6,8 +6,8 @@
     },
     {
         "name": "2.1 (stable)",
-        "version": "2.1.3",
-        "url": "https://cashocs.readthedocs.io/en/v2.1.3/",
+        "version": "2.1.4",
+        "url": "https://cashocs.readthedocs.io/en/v2.1.4/",
         "preferred": true
     },
     {
diff --git a/pyproject.toml b/pyproject.toml
index cdf2df8e..c349916a 100755
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -45,7 +45,7 @@ classifiers = [
 dependencies = [
 	"meshio >= 4.1.0",
 	"numpy",
-    "scipy",
+        "scipy",
 	"typing_extensions",
         "matplotlib",
 ]