diff --git a/meanas/eigensolvers.py b/meanas/eigensolvers.py
index 4b96f60..6b3fba2 100644
--- a/meanas/eigensolvers.py
+++ b/meanas/eigensolvers.py
@@ -25,7 +25,7 @@ def power_iteration(
     Returns:
         (Largest-magnitude eigenvalue, Corresponding eigenvector estimate)
     """
-    if numpy.any(numpy.equal(guess_vector, None)):
+    if guess_vector is None:
         v = numpy.random.rand(operator.shape[0]) + 1j * numpy.random.rand(operator.shape[0])
     else:
         v = guess_vector
diff --git a/meanas/fdfd/functional.py b/meanas/fdfd/functional.py
index 1a1506b..cfc6187 100644
--- a/meanas/fdfd/functional.py
+++ b/meanas/fdfd/functional.py
@@ -45,7 +45,7 @@ def e_full(
         curls = ch(mu * ce(e))          # type: ignore   # mu = None ok because we don't return the function
         return curls - omega ** 2 * epsilon * e
 
-    if numpy.any(numpy.equal(mu, None)):
+    if mu is None:
         return op_1
     else:
         return op_mu
@@ -82,7 +82,7 @@ def eh_full(
         return (ch(h) - 1j * omega * epsilon * e,
                 ce(e) + 1j * omega * mu * h)            # type: ignore   # mu=None ok
 
-    if numpy.any(numpy.equal(mu, None)):
+    if mu is None:
         return op_1
     else:
         return op_mu
@@ -114,7 +114,7 @@ def e2h(
     def e2h_mu(e: cfdfield_t) -> cfdfield_t:
         return ce(e) / (-1j * omega * mu)       # type: ignore   # mu=None ok
 
-    if numpy.any(numpy.equal(mu, None)):
+    if mu is None:
         return e2h_1_1
     else:
         return e2h_mu
@@ -149,7 +149,7 @@ def m2j(
         J = ch(m) / (-1j * omega)
         return J
 
-    if numpy.any(numpy.equal(mu, None)):
+    if mu is None:
         return m2j_1
     else:
         return m2j_mu
diff --git a/meanas/fdfd/operators.py b/meanas/fdfd/operators.py
index ea45cba..86d9fc8 100644
--- a/meanas/fdfd/operators.py
+++ b/meanas/fdfd/operators.py
@@ -77,18 +77,18 @@ def e_full(
     ch = curl_back(dxes[1])
     ce = curl_forward(dxes[0])
 
-    if numpy.any(numpy.equal(pec, None)):
+    if pec is None:
         pe = sparse.eye(epsilon.size)
     else:
         pe = sparse.diags(numpy.where(pec, 0, 1))     # Set pe to (not PEC)
 
-    if numpy.any(numpy.equal(pmc, None)):
+    if pmc is None:
         pm = sparse.eye(epsilon.size)
     else:
         pm = sparse.diags(numpy.where(pmc, 0, 1))     # set pm to (not PMC)
 
     e = sparse.diags(epsilon)
-    if numpy.any(numpy.equal(mu, None)):
+    if mu is None:
         m_div = sparse.eye(epsilon.size)
     else:
         m_div = sparse.diags(1 / mu)
@@ -161,12 +161,12 @@ def h_full(
     ch = curl_back(dxes[1])
     ce = curl_forward(dxes[0])
 
-    if numpy.any(numpy.equal(pec, None)):
+    if pec is None:
         pe = sparse.eye(epsilon.size)
     else:
         pe = sparse.diags(numpy.where(pec, 0, 1))    # set pe to (not PEC)
 
-    if numpy.any(numpy.equal(pmc, None)):
+    if pmc is None:
         pm = sparse.eye(epsilon.size)
     else:
         pm = sparse.diags(numpy.where(pmc, 0, 1))    # Set pe to (not PMC)
@@ -227,19 +227,19 @@ def eh_full(
     Returns:
         Sparse matrix containing the wave operator.
     """
-    if numpy.any(numpy.equal(pec, None)):
+    if pec is None:
         pe = sparse.eye(epsilon.size)
     else:
         pe = sparse.diags(numpy.where(pec, 0, 1))    # set pe to (not PEC)
 
-    if numpy.any(numpy.equal(pmc, None)):
+    if pmc is None:
         pm = sparse.eye(epsilon.size)
     else:
         pm = sparse.diags(numpy.where(pmc, 0, 1))    # set pm to (not PMC)
 
     iwe = pe @ (1j * omega * sparse.diags(epsilon)) @ pe
     iwm = 1j * omega
-    if not numpy.any(numpy.equal(mu, None)):
+    if mu is not None:
         iwm *= sparse.diags(mu)
     iwm = pm @ iwm @ pm
 
@@ -274,10 +274,10 @@ def e2h(
     """
     op = curl_forward(dxes[0]) / (-1j * omega)
 
-    if not numpy.any(numpy.equal(mu, None)):
+    if mu is not None:
         op = sparse.diags(1 / mu) @ op
 
-    if not numpy.any(numpy.equal(pmc, None)):
+    if pmc is not None:
         op = sparse.diags(numpy.where(pmc, 0, 1)) @ op
 
     return op
@@ -302,7 +302,7 @@ def m2j(
     """
     op = curl_back(dxes[1]) / (1j * omega)
 
-    if not numpy.any(numpy.equal(mu, None)):
+    if mu is not None:
         op = op @ sparse.diags(1 / mu)
 
     return op
diff --git a/meanas/fdfd/waveguide_2d.py b/meanas/fdfd/waveguide_2d.py
index 69d0b19..3ac95d5 100644
--- a/meanas/fdfd/waveguide_2d.py
+++ b/meanas/fdfd/waveguide_2d.py
@@ -239,7 +239,7 @@ def operator_e(
     Returns:
         Sparse matrix representation of the operator.
     """
-    if numpy.any(numpy.equal(mu, None)):
+    if mu is None:
         mu = numpy.ones_like(epsilon)
 
     Dfx, Dfy = deriv_forward(dxes[0])
@@ -306,7 +306,7 @@ def operator_h(
     Returns:
         Sparse matrix representation of the operator.
     """
-    if numpy.any(numpy.equal(mu, None)):
+    if mu is None:
         mu = numpy.ones_like(epsilon)
 
     Dfx, Dfy = deriv_forward(dxes[0])
@@ -513,7 +513,7 @@ def hxy2h(
     Dfx, Dfy = deriv_forward(dxes[0])
     hxy2hz = sparse.hstack((Dfx, Dfy)) / (1j * wavenumber)
 
-    if not numpy.any(numpy.equal(mu, None)):
+    if mu is not None:
         mu_parts = numpy.split(mu, 3)
         mu_xy = sparse.diags(numpy.hstack((mu_parts[0], mu_parts[1])))
         mu_z_inv = sparse.diags(1 / mu_parts[2])
@@ -547,7 +547,7 @@ def exy2e(
     Dbx, Dby = deriv_back(dxes[1])
     exy2ez = sparse.hstack((Dbx, Dby)) / (1j * wavenumber)
 
-    if not numpy.any(numpy.equal(epsilon, None)):
+    if epsilon is not None:
         epsilon_parts = numpy.split(epsilon, 3)
         epsilon_xy = sparse.diags(numpy.hstack((epsilon_parts[0], epsilon_parts[1])))
         epsilon_z_inv = sparse.diags(1 / epsilon_parts[2])
@@ -580,7 +580,7 @@ def e2h(
         Sparse matrix representation of the operator.
     """
     op = curl_e(wavenumber, dxes) / (-1j * omega)
-    if not numpy.any(numpy.equal(mu, None)):
+    if mu is not None:
         op = sparse.diags(1 / mu) @ op
     return op
 
@@ -675,7 +675,7 @@ def h_err(
 
     eps_inv = sparse.diags(1 / epsilon)
 
-    if numpy.any(numpy.equal(mu, None)):
+    if mu is None:
         op = ce @ eps_inv @ ch @ h - omega ** 2 * h
     else:
         op = ce @ eps_inv @ ch @ h - omega ** 2 * (mu * h)
@@ -708,7 +708,7 @@ def e_err(
     ce = curl_e(wavenumber, dxes)
     ch = curl_h(wavenumber, dxes)
 
-    if numpy.any(numpy.equal(mu, None)):
+    if mu is None:
         op = ch @ ce @ e - omega ** 2 * (epsilon * e)
     else:
         mu_inv = sparse.diags(1 / mu)
diff --git a/meanas/fdmath/vectorization.py b/meanas/fdmath/vectorization.py
index ef97f7c..5d9e932 100644
--- a/meanas/fdmath/vectorization.py
+++ b/meanas/fdmath/vectorization.py
@@ -40,7 +40,7 @@ def vec(f: Union[fdfield_t, cfdfield_t, ArrayLike, None]) -> Union[vfdfield_t, v
     Returns:
         1D ndarray containing the linearized field (or `None`)
     """
-    if numpy.any(numpy.equal(f, None)):
+    if f is None:
         return None
     return numpy.ravel(f, order='C')
 
@@ -72,7 +72,7 @@ def unvec(v: Union[vfdfield_t, vcfdfield_t, None], shape: Sequence[int]) -> Unio
     Returns:
         `[f_x, f_y, f_z]` where each `f_` is a `len(shape)` dimensional ndarray (or `None`)
     """
-    if numpy.any(numpy.equal(v, None)):
+    if v is None:
         return None
     return v.reshape((3, *shape), order='C')