diff --git a/src/openfermion/contrib/representability/_dualbasis.py b/src/openfermion/contrib/representability/_dualbasis.py
index 1369388e3..80f602ce4 100644
--- a/src/openfermion/contrib/representability/_dualbasis.py
+++ b/src/openfermion/contrib/representability/_dualbasis.py
@@ -1,4 +1,4 @@
-from typing import List, Optional, Union, Tuple
+from collections.abc import Sequence
 import copy
 
 
@@ -18,14 +18,20 @@ class DualBasisElement:
     for all i in [1, dim(`M')].
     """
 
+    primal_tensors_names: list[str]
+    primal_elements: list[tuple[int, ...]]
+    primal_coeffs: list[float]
+    constant_bias: float | int
+    dual_scalar: float | int
+
     def __init__(
         self,
         *,
-        tensor_names: Optional[Union[None, List[str]]] = None,
-        tensor_elements: Optional[Union[None, List[Tuple[int, ...]]]] = None,
-        tensor_coeffs: Optional[Union[None, List[float]]] = None,
-        bias: Optional[int] = 0,
-        scalar: Optional[int] = 0,
+        tensor_names: list[str] | None = None,
+        tensor_elements: Sequence[tuple[int, ...]] | None = None,
+        tensor_coeffs: list[float] | None = None,
+        bias: float | int = 0,
+        scalar: float | int = 0,
     ):
         """
         Define a linear operator on a tensor `A', a bias `b', and a result `c'
@@ -52,7 +58,7 @@ def __init__(
         if tensor_elements is None:
             self.primal_elements = []
         else:
-            self.primal_elements = tensor_elements
+            self.primal_elements = list(tensor_elements)
 
         if tensor_coeffs is None:
             self.primal_coeffs = []
@@ -155,7 +161,9 @@ def __add__(self, other):
 
 
 class DualBasis:
-    def __init__(self, elements: Optional[Union[None, List[DualBasisElement]]] = None):
+    elements: list[DualBasisElement]
+
+    def __init__(self, elements: list[DualBasisElement] | None = None):
         """
         A collection of DualBasisElements
 
diff --git a/src/openfermion/contrib/representability/constraints/spin_orbital_2pos_constraints.py b/src/openfermion/contrib/representability/constraints/spin_orbital_2pos_constraints.py
index add6fd4b0..3219c5d1c 100644
--- a/src/openfermion/contrib/representability/constraints/spin_orbital_2pos_constraints.py
+++ b/src/openfermion/contrib/representability/constraints/spin_orbital_2pos_constraints.py
@@ -236,9 +236,7 @@ def tpdm_to_phdm_mapping(dim: int) -> DualBasis:
     """
     dbe_list = []
 
-    def g2d2map(
-        p: int, q: int, r: int, s: int, factor: Optional[Union[float, int]] = 1
-    ) -> DualBasisElement:
+    def g2d2map(p: int, q: int, r: int, s: int, factor: float | int = 1) -> DualBasisElement:
         """
         Build the dual basis element for a symmetric 2-marginal
 
diff --git a/src/openfermion/hamiltonians/special_operators.py b/src/openfermion/hamiltonians/special_operators.py
index 92f2f0cb2..1030f6b9d 100644
--- a/src/openfermion/hamiltonians/special_operators.py
+++ b/src/openfermion/hamiltonians/special_operators.py
@@ -11,8 +11,6 @@
 #   limitations under the License.
 """Commonly used operators (mainly instances of SymbolicOperator)."""
 
-from typing import Optional, Union, Tuple
-
 import openfermion.config as config
 from openfermion.ops.operators import BosonOperator, FermionOperator
 from openfermion.utils.indexing import down_index, up_index
@@ -213,7 +211,7 @@ def s_squared_operator(n_spatial_orbitals: int) -> FermionOperator:
 
 
 def majorana_operator(
-    term: Optional[Union[Tuple[int, int], str]] = None, coefficient=1.0
+    term: tuple[int, int] | str | None = None, coefficient: int | float | complex = 1.0
 ) -> FermionOperator:
     r"""Initialize a Majorana operator.
 
@@ -244,15 +242,15 @@ def majorana_operator(
 
     # If term is a string, convert it to a tuple
     if isinstance(term, str):
-        operator_type = term[0]
+        operator_type_char = term[0]
         mode = int(term[1:])
-        if operator_type == 'c':
-            operator_type = 0
-        elif operator_type == 'd':
-            operator_type = 1
+        if operator_type_char == 'c':
+            operator_type_int = 0
+        elif operator_type_char == 'd':
+            operator_type_int = 1
         else:
-            raise ValueError('Invalid operator type: {}'.format(operator_type))
-        term = (mode, operator_type)
+            raise ValueError('Invalid operator type: {}'.format(operator_type_char))
+        term = (mode, operator_type_int)
 
     # Process term
 
@@ -268,8 +266,9 @@ def majorana_operator(
             majorana_op = FermionOperator(((mode, 1),), coefficient)
             majorana_op += FermionOperator(((mode, 0),), coefficient)
         elif operator_type == 1:
-            majorana_op = FermionOperator(((mode, 1),), 1.0j * coefficient)
-            majorana_op -= FermionOperator(((mode, 0),), 1.0j * coefficient)
+            imag_coeff = 1.0j * coefficient
+            majorana_op = FermionOperator(((mode, 1),), imag_coeff)
+            majorana_op -= FermionOperator(((mode, 0),), imag_coeff)
         else:
             raise ValueError('Invalid operator type: {}'.format(str(operator_type)))
 
@@ -281,8 +280,8 @@ def majorana_operator(
 
 
 def number_operator(
-    n_modes: int, mode: Optional[int] = None, coefficient=1.0, parity: int = -1
-) -> Union[BosonOperator, FermionOperator]:
+    n_modes: int, mode: int | None = None, coefficient=1.0, parity: int = -1
+) -> BosonOperator | FermionOperator:
     """Return a fermionic or bosonic number operator.
 
     Args:
@@ -299,7 +298,7 @@ def number_operator(
     """
 
     if parity == -1:
-        Op = FermionOperator
+        Op: type[FermionOperator] | type[BosonOperator] = FermionOperator
     elif parity == 1:
         Op = BosonOperator
     else: