diff --git a/sumpy/expansion/local.py b/sumpy/expansion/local.py
index b4b435de5ce1aa778b5c176b9a6df18b3d88c8dc..e7e0403e86601d9225fe769e10c8d7eac0be85c5 100644
--- a/sumpy/expansion/local.py
+++ b/sumpy/expansion/local.py
@@ -332,10 +332,10 @@ class _FourierBesselLocalExpansion(LocalExpansionBase):
         source_angle_rel_center = sym.atan2(-avec[1], -avec[0])
         avec_len = sym_real_norm_2(avec)
         return [self.kernel.postprocess_at_source(
-                    hankel_1(l, arg_scale * avec_len)
-                    * rscale ** abs(l)
-                    * sym.exp(sym.I * l * source_angle_rel_center), avec)
-                    for l in self.get_coefficient_identifiers()]
+                    hankel_1(c, arg_scale * avec_len)
+                    * rscale ** abs(c)
+                    * sym.exp(sym.I * c * source_angle_rel_center), avec)
+                    for c in self.get_coefficient_identifiers()]
 
     def evaluate(self, coeffs, bvec, rscale):
         if not self.use_rscale:
@@ -348,12 +348,12 @@ class _FourierBesselLocalExpansion(LocalExpansionBase):
 
         arg_scale = self.get_bessel_arg_scaling()
 
-        return sum(coeffs[self.get_storage_index(l)]
+        return sum(coeffs[self.get_storage_index(c)]
                    * self.kernel.postprocess_at_target(
-                       bessel_j(l, arg_scale * bvec_len)
-                       / rscale ** abs(l)
-                       * sym.exp(sym.I * l * -target_angle_rel_center), bvec)
-                for l in self.get_coefficient_identifiers())
+                       bessel_j(c, arg_scale * bvec_len)
+                       / rscale ** abs(c)
+                       * sym.exp(sym.I * c * -target_angle_rel_center), bvec)
+                for c in self.get_coefficient_identifiers())
 
     def translate_from(self, src_expansion, src_coeff_exprs, src_rscale,
             dvec, tgt_rscale):
@@ -370,13 +370,13 @@ class _FourierBesselLocalExpansion(LocalExpansionBase):
             bessel_j = sym.Function("bessel_j")
             new_center_angle_rel_old_center = sym.atan2(dvec[1], dvec[0])
             translated_coeffs = []
-            for l in self.get_coefficient_identifiers():
+            for j in self.get_coefficient_identifiers():
                 translated_coeffs.append(
                     sum(src_coeff_exprs[src_expansion.get_storage_index(m)]
-                        * bessel_j(m - l, arg_scale * dvec_len)
+                        * bessel_j(m - j, arg_scale * dvec_len)
                         / src_rscale ** abs(m)
-                        * tgt_rscale ** abs(l)
-                        * sym.exp(sym.I * (m - l) * -new_center_angle_rel_old_center)
+                        * tgt_rscale ** abs(j)
+                        * sym.exp(sym.I * (m - j) * -new_center_angle_rel_old_center)
                     for m in src_expansion.get_coefficient_identifiers()))
             return translated_coeffs
 
@@ -385,14 +385,14 @@ class _FourierBesselLocalExpansion(LocalExpansionBase):
             hankel_1 = sym.Function("hankel_1")
             new_center_angle_rel_old_center = sym.atan2(dvec[1], dvec[0])
             translated_coeffs = []
-            for l in self.get_coefficient_identifiers():
+            for j in self.get_coefficient_identifiers():
                 translated_coeffs.append(
                     sum(
-                        (-1) ** l
-                        * hankel_1(m + l, arg_scale * dvec_len)
+                        (-1) ** j
+                        * hankel_1(m + j, arg_scale * dvec_len)
                         * src_rscale ** abs(m)
-                        * tgt_rscale ** abs(l)
-                        * sym.exp(sym.I * (m + l) * new_center_angle_rel_old_center)
+                        * tgt_rscale ** abs(j)
+                        * sym.exp(sym.I * (m + j) * new_center_angle_rel_old_center)
                         * src_coeff_exprs[src_expansion.get_storage_index(m)]
                         for m in src_expansion.get_coefficient_identifiers()))
             return translated_coeffs
diff --git a/sumpy/expansion/multipole.py b/sumpy/expansion/multipole.py
index e3759d4342e42a3d81c5f220ca077681f325001c..25148e2b00d2d7a80c3b8f0ce3885aa94a1b59a7 100644
--- a/sumpy/expansion/multipole.py
+++ b/sumpy/expansion/multipole.py
@@ -290,11 +290,11 @@ class _HankelBased2DMultipoleExpansion(MultipoleExpansionBase):
         source_angle_rel_center = sym.atan2(-avec[1], -avec[0])
         return [
                 self.kernel.postprocess_at_source(
-                    bessel_j(l, arg_scale * avec_len)
-                    / rscale ** abs(l)
-                    * sym.exp(sym.I * l * -source_angle_rel_center),
+                    bessel_j(c, arg_scale * avec_len)
+                    / rscale ** abs(c)
+                    * sym.exp(sym.I * c * -source_angle_rel_center),
                     avec)
-                for l in self.get_coefficient_identifiers()]
+                for c in self.get_coefficient_identifiers()]
 
     def evaluate(self, coeffs, bvec, rscale):
         if not self.use_rscale:
@@ -307,12 +307,12 @@ class _HankelBased2DMultipoleExpansion(MultipoleExpansionBase):
 
         arg_scale = self.get_bessel_arg_scaling()
 
-        return sum(coeffs[self.get_storage_index(l)]
+        return sum(coeffs[self.get_storage_index(c)]
                    * self.kernel.postprocess_at_target(
-                       hankel_1(l, arg_scale * bvec_len)
-                       * rscale ** abs(l)
-                       * sym.exp(sym.I * l * target_angle_rel_center), bvec)
-                for l in self.get_coefficient_identifiers())
+                       hankel_1(c, arg_scale * bvec_len)
+                       * rscale ** abs(c)
+                       * sym.exp(sym.I * c * target_angle_rel_center), bvec)
+                for c in self.get_coefficient_identifiers())
 
     def translate_from(self, src_expansion, src_coeff_exprs, src_rscale,
             dvec, tgt_rscale):
@@ -333,13 +333,13 @@ class _HankelBased2DMultipoleExpansion(MultipoleExpansionBase):
         arg_scale = self.get_bessel_arg_scaling()
 
         translated_coeffs = []
-        for l in self.get_coefficient_identifiers():
+        for j in self.get_coefficient_identifiers():
             translated_coeffs.append(
                 sum(src_coeff_exprs[src_expansion.get_storage_index(m)]
-                    * bessel_j(m - l, arg_scale * dvec_len)
+                    * bessel_j(m - j, arg_scale * dvec_len)
                     * src_rscale ** abs(m)
-                    / tgt_rscale ** abs(l)
-                    * sym.exp(sym.I * (m - l) * new_center_angle_rel_old_center)
+                    / tgt_rscale ** abs(j)
+                    * sym.exp(sym.I * (m - j) * new_center_angle_rel_old_center)
                 for m in src_expansion.get_coefficient_identifiers()))
         return translated_coeffs
 
diff --git a/sumpy/point_calculus.py b/sumpy/point_calculus.py
index f8b4f5bdf38ec7cd9043047f07e5f2c2ce31b69c..fc164d81b226c0ed2e38f2e60c68d4cd466c8b7e 100644
--- a/sumpy/point_calculus.py
+++ b/sumpy/point_calculus.py
@@ -231,9 +231,9 @@ class CalculusPatch(object):
         from pytools.obj_array import make_obj_array
         return make_obj_array([
             sum(
-                levi_civita((l, m, n)) * self.diff(m, arg[n])
+                levi_civita((k, m, n)) * self.diff(m, arg[n])
                 for m in range(3) for n in range(3))
-            for l in range(3)])
+            for k in range(3)])
 
     def eval_at_center(self, f_values):
         """Interpolate *f_values* to the center point.
diff --git a/test/test_cse.py b/test/test_cse.py
index 703e8c9f93ca2069c273ead893a00ed0ba3dc697..2f54f34bab36e9be4047c05f1a2c663523a83509 100644
--- a/test/test_cse.py
+++ b/test/test_cse.py
@@ -253,12 +253,12 @@ def test_issue_4499():
     G = Function('G')  # noqa
     t = Tuple(
         *(a, a + S(1)/2, 2*a, b, 2*a - b + 1, (sqrt(z)/2)**(-2*a + 1)*B(2*a
-        - b, sqrt(z))*B(b - 1, sqrt(z))*G(b)*G(2*a - b + 1),
-        sqrt(z)*(sqrt(z)/2)**(-2*a + 1)*B(b, sqrt(z))*B(2*a - b,
+        - b, sqrt(z))*B(b - 1, sqrt(z))*G(b)*G(2*a - b + 1),  # noqa
+        sqrt(z)*(sqrt(z)/2)**(-2*a + 1)*B(b, sqrt(z))*B(2*a - b,  # noqa
         sqrt(z))*G(b)*G(2*a - b + 1), sqrt(z)*(sqrt(z)/2)**(-2*a + 1)*B(b - 1,
         sqrt(z))*B(2*a - b + 1, sqrt(z))*G(b)*G(2*a - b + 1),
-        (sqrt(z)/2)**(-2*a + 1)*B(b, sqrt(z))*B(2*a - b + 1,
-        sqrt(z))*G(b)*G(2*a - b + 1), 1, 0, S(1)/2, z/2, -b + 1, -2*a + b,
+        (sqrt(z)/2)**(-2*a + 1)*B(b, sqrt(z))*B(2*a - b + 1,  # noqa
+        sqrt(z))*G(b)*G(2*a - b + 1), 1, 0, S(1)/2, z/2, -b + 1, -2*a + b, # noqa
         -2*a))  # noqa
     c = cse(t)
     assert len(c[0]) == 11