diff --git a/src/pylife/mesh/gradient.py b/src/pylife/mesh/gradient.py
index ad36b600..763f4973 100644
--- a/src/pylife/mesh/gradient.py
+++ b/src/pylife/mesh/gradient.py
@@ -311,27 +311,27 @@ def _compute_gradient_simplex(self, df):
         
         nodal_values = df.iloc[:4,3]
         self._initialize_ansatz_function_derivative_simplex()
+    
+        # compute jacobian matrix of the mapping from reference space [0,1]^3 to world space
+        J11 = -x11 + x21
+        J21 = -x12 + x22
+        J31 = -x13 + x23
+        J12 = -x11 + x31
+        J22 = -x12 + x32
+        J32 = -x13 + x33
+        J13 = -x11 + x41
+        J23 = -x12 + x42
+        J33 = -x13 + x43
+        J = np.array([[J11,J12,J13],[J21,J22,J23],[J31,J32,J33]])
+        
+        # invert jacobian to map from world space to reference domain
+        try:
+            Jinv = np.linalg.inv(J)
+        except np.linalg.LinAlgError:
+            return df
         
         for node_index in range(4):
             
-            # compute jacobian matrix of the mapping from reference space [0,1]^3 to world space
-            J11 = -x11 + x21
-            J21 = -x12 + x22
-            J31 = -x13 + x23
-            J12 = -x11 + x31
-            J22 = -x12 + x32
-            J32 = -x13 + x33
-            J13 = -x11 + x41
-            J23 = -x12 + x42
-            J33 = -x13 + x43
-            J = np.array([[J11,J12,J13],[J21,J22,J23],[J31,J32,J33]])
-            
-            # invert jacobian to map from world space to reference domain
-            try:
-                Jinv = np.linalg.inv(J)
-            except np.linalg.LinAlgError:
-                continue
-            
             # loop over derivative index d/dx_k
             df.iloc[node_index,4:7] \
                 = self._compute_gradient_simplex_single_node(node_index, nodal_values, Jinv)