From ecba4f1182b1d9b9f05fddd4867df4e925a664d6 Mon Sep 17 00:00:00 2001
From: Sam Evans <s7evans11@gmail.com>
Date: Wed, 15 Jul 2026 13:45:26 -0400
Subject: [PATCH 1/8] Adds Grid.compute_face_node_angles (WIP)

This provides an initial implementation which seems to be working. See #1566

TODO:
- support assume_convex=False
- type-hint possible UxDataArray return type?
- tests
- examples
---
 uxarray/grid/angles.py | 68 +++++++++++++++++++++++++++++++++++++
 uxarray/grid/grid.py   | 77 ++++++++++++++++++++++++++++++++++++++++++
 2 files changed, 145 insertions(+)
 create mode 100644 uxarray/grid/angles.py

diff --git a/uxarray/grid/angles.py b/uxarray/grid/angles.py
new file mode 100644
index 000000000..a95181dc2
--- /dev/null
+++ b/uxarray/grid/angles.py
@@ -0,0 +1,68 @@
+"""
+Purpose: angle calculations on a grid
+"""
+
+from numba import njit, prange
+import numpy as np
+import xarray as xr
+
+from uxarray.grid.utils import _small_angle_of_2_vectors
+
+
+@njit(cache=True, parallel=True)
+def _compute_face_node_angles_convex(
+    node_x,
+    node_y,
+    node_z,
+    face_node_connectivity,
+    n_nodes_per_face,
+    *,
+    geometry="spherical",
+):
+    """
+    Calculate the angles at each node for each face assuming the faces are convex.
+
+    Parameters
+    ----------
+    node_x : np.ndarray with shape (n_nodes,)
+        X coordinates of the nodes.
+    node_y : np.ndarray with shape (n_nodes,)
+        Y coordinates of the nodes.
+    node_z : np.ndarray with shape (n_nodes,)
+        Z coordinates of the nodes.
+    face_node_connectivity : np.ndarray with shape (n_faces, n_max_face_nodes)
+        Connectivity array defining which nodes form each face.
+    n_nodes_per_face : np.ndarray with shape (n_faces,)
+        Number of nodes for each face.
+    geometry : str, "spherical" or "flat", default "spherical"
+        The geometry to use. "spherical" respects the true underlying geometry,
+        by projecting edges onto the tangent plane at each node.
+        "flat" finds angles between chords, which is less accurate but also less expensive.
+
+    Returns
+    -------
+    np.ndarray with shape (n_faces, n_max_face_nodes)
+        Angles at each node of each face.
+        INT_FILL_VALUE elements from face_node_connectivity correspond with np.nan in the result.
+    """
+    n_faces, n_max_face_nodes = face_node_connectivity.shape
+    result = np.full(face_node_connectivity.shape, np.nan, dtype=np.float64)
+    for i in prange(n_faces):
+        n_nodes = n_nodes_per_face[i]
+        for j in range(n_nodes):
+            ihere = face_node_connectivity[i, j]
+            iprev = face_node_connectivity[i, (j - 1) % n_nodes]
+            inext = face_node_connectivity[i, (j + 1) % n_nodes]
+            xhere = node_x[ihere]
+            yhere = node_y[ihere]
+            zhere = node_z[ihere]
+            v1 = np.array([node_x[iprev] - xhere, node_y[iprev] - yhere, node_z[iprev] - zhere])
+            v2 = np.array([node_x[inext] - xhere, node_y[inext] - yhere, node_z[inext] - zhere])
+            if geometry == "spherical":
+                # Project onto tangent plane at the current node
+                normal = np.array([xhere, yhere, zhere])
+                normal /= np.linalg.norm(normal)
+                v1 -= np.dot(v1, normal) * normal
+                v2 -= np.dot(v2, normal) * normal
+            result[i, j] = _small_angle_of_2_vectors(v1, v2)
+    return result
diff --git a/uxarray/grid/grid.py b/uxarray/grid/grid.py
index f79c6586a..ea1e69b57 100644
--- a/uxarray/grid/grid.py
+++ b/uxarray/grid/grid.py
@@ -17,6 +17,7 @@
 from uxarray.core.utils import _open_dataset_with_fallback
 from uxarray.cross_sections import GridCrossSectionAccessor
 from uxarray.formatting_html import grid_repr
+from uxarray.grid.angles import _compute_face_node_angles_convex
 from uxarray.grid.area import get_all_face_area_from_coords
 from uxarray.grid.bounds import _populate_face_bounds
 from uxarray.grid.connectivity import (
@@ -1929,6 +1930,82 @@ def copy(self):
             source_dims_dict=self._source_dims_dict,
         )
 
+    def compute_face_node_angles(
+        self,
+        geometry: str = "spherical",
+        *,
+        degrees: bool = False,
+        assume_convex: bool = False,
+        cache: bool | None = None,
+        as_uxarray: bool = False,
+    ) -> xr.DataArray:
+        """Compute the angles at each node of each face in the grid.
+
+        Parameters
+        ----------
+        geometry : str, "spherical" or "flat", defaults to "spherical"
+            The geometry to use for angle computation.
+            If "spherical", angles are computed considering the tangent plane at each node.
+            If "flat", angles are computed in 3D Cartesian space,
+            ignoring the true spherical geometry of the grid, but may be faster to compute.
+        degrees : bool, defaults to False
+            Whether to return angles in degrees (if True) or radians (if False).
+        assume_convex : bool, defaults to False
+            Whether to assume that all faces are convex, i.e. all internal angles less than 180 degrees.
+            If True, uses a more efficient algorithm that will produce incorrect results for non-convex faces.
+        cache : None or bool, defaults to None
+            Whether to cache the computed angles in the grid's dataset, in face_node_angles (if "spherical" geometry)
+            or face_node_angles_flat (if "flat" or "euclidean" geometry). Cached angles are always in radians.
+            If None, use cached result if available, else compute but do not cache result.
+            If True, use cached result if available, else compute and cache result.
+            If False, always recompute; do not check or store in cache.
+        as_uxarray : bool, defaults to False
+            Whether to return a uxarray.DataArray (if True) instead of an xarray.DataArray (if False).
+            If True, equivalent to uxarray.DataArray(self.compute_face_node_angles(..., as_uxarray=False), uxgrid=self).
+
+        Returns
+        -------
+        face_node_angles : xr.DataArray or uxarray.UxDataArray (if as_uxarray=True)
+            The internal angles at each node, for each face in the grid.
+            Has 'n_face' and 'n_max_face_nodes' dimensions, with same size as in self.
+        """
+        from uxarray.conventions.ugrid import FACE_DIM, N_MAX_FACE_NODES_DIM
+
+        if geometry not in ("spherical", "flat"):
+            raise ValueError(
+                f"Invalid geometry {geometry!r}; expected 'spherical' or 'flat'."
+            )
+        name = "face_node_angles" if geometry == "spherical" else "face_node_angles_flat"
+        result = None
+        if cache is None or cache:
+            if name in self._ds:
+                result = self._ds[name]
+        if result is None:  # result not in cache; need to compute it
+            if assume_convex:
+                result = _compute_face_node_angles_convex(
+                    self.node_x.values, self.node_y.values, self.node_z.values,
+                    self.face_node_connectivity.values, self.n_nodes_per_face.values,
+                    geometry=geometry,
+                )
+            else:
+                raise NotImplementedError('[TODO] Non-convex face angle computation.')
+            # convert to xr.DataArray
+            result = xr.DataArray(
+                data=result,
+                dims=[FACE_DIM, N_MAX_FACE_NODES_DIM],
+                name=name,
+                attrs={"description": ("Internal angles at each node of each face. "
+                                      f"(geometry={geometry}, assume_convex={assume_convex})")},
+            )
+            if cache:
+                self._ds[name] = result
+        if degrees:
+            result = np.rad2deg(result)
+        if as_uxarray:
+            from uxarray.core.dataarray import UxDataArray  # import at runtime to avoid circular import
+            result = UxDataArray(result, uxgrid=self)
+        return result
+
     def calculate_total_face_area(
         self,
         quadrature_rule: str | None = "triangular",

From f7866b48dea974a4c77c503e5a45a4e39dbf71d2 Mon Sep 17 00:00:00 2001
From: Sam Evans <47793072+Sevans711@users.noreply.github.com>
Date: Tue, 11 Aug 2026 15:21:50 -0400
Subject: [PATCH 2/8] face_node_angles remove unneeded kwargs

removes "geometry" (always use spherical)
removes "assume_convex" (always assume convex)

Also, improves type-hinting for compute_face_node_angles().

Also, formats using pre-commit ruff formatting.
---
 uxarray/grid/angles.py | 31 ++++++++++++++-----------------
 uxarray/grid/grid.py   | 42 ++++++++++++++++--------------------------
 2 files changed, 30 insertions(+), 43 deletions(-)

diff --git a/uxarray/grid/angles.py b/uxarray/grid/angles.py
index a95181dc2..d20bb64e5 100644
--- a/uxarray/grid/angles.py
+++ b/uxarray/grid/angles.py
@@ -2,9 +2,8 @@
 Purpose: angle calculations on a grid
 """
 
-from numba import njit, prange
 import numpy as np
-import xarray as xr
+from numba import njit, prange
 
 from uxarray.grid.utils import _small_angle_of_2_vectors
 
@@ -16,11 +15,10 @@ def _compute_face_node_angles_convex(
     node_z,
     face_node_connectivity,
     n_nodes_per_face,
-    *,
-    geometry="spherical",
 ):
     """
-    Calculate the angles at each node for each face assuming the faces are convex.
+    Calculate the angles at each node for each face, assuming convex faces
+    and a spherical geometry (these assumptions occur throughout uxarray).
 
     Parameters
     ----------
@@ -34,10 +32,6 @@ def _compute_face_node_angles_convex(
         Connectivity array defining which nodes form each face.
     n_nodes_per_face : np.ndarray with shape (n_faces,)
         Number of nodes for each face.
-    geometry : str, "spherical" or "flat", default "spherical"
-        The geometry to use. "spherical" respects the true underlying geometry,
-        by projecting edges onto the tangent plane at each node.
-        "flat" finds angles between chords, which is less accurate but also less expensive.
 
     Returns
     -------
@@ -56,13 +50,16 @@ def _compute_face_node_angles_convex(
             xhere = node_x[ihere]
             yhere = node_y[ihere]
             zhere = node_z[ihere]
-            v1 = np.array([node_x[iprev] - xhere, node_y[iprev] - yhere, node_z[iprev] - zhere])
-            v2 = np.array([node_x[inext] - xhere, node_y[inext] - yhere, node_z[inext] - zhere])
-            if geometry == "spherical":
-                # Project onto tangent plane at the current node
-                normal = np.array([xhere, yhere, zhere])
-                normal /= np.linalg.norm(normal)
-                v1 -= np.dot(v1, normal) * normal
-                v2 -= np.dot(v2, normal) * normal
+            v1 = np.array(
+                [node_x[iprev] - xhere, node_y[iprev] - yhere, node_z[iprev] - zhere]
+            )
+            v2 = np.array(
+                [node_x[inext] - xhere, node_y[inext] - yhere, node_z[inext] - zhere]
+            )
+            # Spherical geometry: project onto tangent plane at the current node
+            normal = np.array([xhere, yhere, zhere])
+            normal /= np.linalg.norm(normal)
+            v1 -= np.dot(v1, normal) * normal
+            v2 -= np.dot(v2, normal) * normal
             result[i, j] = _small_angle_of_2_vectors(v1, v2)
     return result
diff --git a/uxarray/grid/grid.py b/uxarray/grid/grid.py
index cc7ca29f5..5a77f23d6 100644
--- a/uxarray/grid/grid.py
+++ b/uxarray/grid/grid.py
@@ -1952,30 +1952,24 @@ def copy(self):
 
     def compute_face_node_angles(
         self,
-        geometry: str = "spherical",
         *,
         degrees: bool = False,
-        assume_convex: bool = False,
         cache: bool | None = None,
         as_uxarray: bool = False,
-    ) -> xr.DataArray:
+    ) -> xr.DataArray | UxDataArray:
         """Compute the angles at each node of each face in the grid.
+        Assumes convex faces and a spherical geometry (consistent with other uxarray methods).
 
         Parameters
         ----------
-        geometry : str, "spherical" or "flat", defaults to "spherical"
-            The geometry to use for angle computation.
-            If "spherical", angles are computed considering the tangent plane at each node.
-            If "flat", angles are computed in 3D Cartesian space,
-            ignoring the true spherical geometry of the grid, but may be faster to compute.
         degrees : bool, defaults to False
             Whether to return angles in degrees (if True) or radians (if False).
         assume_convex : bool, defaults to False
             Whether to assume that all faces are convex, i.e. all internal angles less than 180 degrees.
             If True, uses a more efficient algorithm that will produce incorrect results for non-convex faces.
         cache : None or bool, defaults to None
-            Whether to cache the computed angles in the grid's dataset, in face_node_angles (if "spherical" geometry)
-            or face_node_angles_flat (if "flat" or "euclidean" geometry). Cached angles are always in radians.
+            Whether to cache the computed angles in the grid's dataset, in face_node_angles data variable.
+            Cached angles are always in radians.
             If None, use cached result if available, else compute but do not cache result.
             If True, use cached result if available, else compute and cache result.
             If False, always recompute; do not check or store in cache.
@@ -1991,38 +1985,34 @@ def compute_face_node_angles(
         """
         from uxarray.conventions.ugrid import FACE_DIM, N_MAX_FACE_NODES_DIM
 
-        if geometry not in ("spherical", "flat"):
-            raise ValueError(
-                f"Invalid geometry {geometry!r}; expected 'spherical' or 'flat'."
-            )
-        name = "face_node_angles" if geometry == "spherical" else "face_node_angles_flat"
+        name = "face_node_angles"
         result = None
         if cache is None or cache:
             if name in self._ds:
                 result = self._ds[name]
         if result is None:  # result not in cache; need to compute it
-            if assume_convex:
-                result = _compute_face_node_angles_convex(
-                    self.node_x.values, self.node_y.values, self.node_z.values,
-                    self.face_node_connectivity.values, self.n_nodes_per_face.values,
-                    geometry=geometry,
-                )
-            else:
-                raise NotImplementedError('[TODO] Non-convex face angle computation.')
+            result = _compute_face_node_angles_convex(
+                self.node_x.values,
+                self.node_y.values,
+                self.node_z.values,
+                self.face_node_connectivity.values,
+                self.n_nodes_per_face.values,
+            )
             # convert to xr.DataArray
             result = xr.DataArray(
                 data=result,
                 dims=[FACE_DIM, N_MAX_FACE_NODES_DIM],
                 name=name,
-                attrs={"description": ("Internal angles at each node of each face. "
-                                      f"(geometry={geometry}, assume_convex={assume_convex})")},
+                attrs={"description": "Internal angles at each node of each face."},
             )
             if cache:
                 self._ds[name] = result
         if degrees:
             result = np.rad2deg(result)
         if as_uxarray:
-            from uxarray.core.dataarray import UxDataArray  # import at runtime to avoid circular import
+            from uxarray.core.dataarray import UxDataArray
+            # ^import at runtime to avoid circular import
+
             result = UxDataArray(result, uxgrid=self)
         return result
 

From f5967ac56522859a40bb7044c1c3751b61c514f6 Mon Sep 17 00:00:00 2001
From: Sam Evans <47793072+Sevans711@users.noreply.github.com>
Date: Tue, 11 Aug 2026 15:58:31 -0400
Subject: [PATCH 3/8] face_node_angles efficiency: tiny arrays in numba

use tuples instead of tiny arrays in numba!

Cleans up docstring in grid.py (forgot to change in previous commit).
---
 uxarray/grid/angles.py | 34 +++++++++++++++++++---------
 uxarray/grid/grid.py   |  4 ----
 uxarray/grid/utils.py  | 51 +++++++++++++++++++++++++++++++++++++-----
 3 files changed, 68 insertions(+), 21 deletions(-)

diff --git a/uxarray/grid/angles.py b/uxarray/grid/angles.py
index d20bb64e5..678f94cb4 100644
--- a/uxarray/grid/angles.py
+++ b/uxarray/grid/angles.py
@@ -5,7 +5,7 @@
 import numpy as np
 from numba import njit, prange
 
-from uxarray.grid.utils import _small_angle_of_2_vectors
+from uxarray.grid.utils import _numba_norm3, _small_angle_of_2_vectors
 
 
 @njit(cache=True, parallel=True)
@@ -39,7 +39,7 @@ def _compute_face_node_angles_convex(
         Angles at each node of each face.
         INT_FILL_VALUE elements from face_node_connectivity correspond with np.nan in the result.
     """
-    n_faces, n_max_face_nodes = face_node_connectivity.shape
+    n_faces, _n_max_face_nodes = face_node_connectivity.shape
     result = np.full(face_node_connectivity.shape, np.nan, dtype=np.float64)
     for i in prange(n_faces):
         n_nodes = n_nodes_per_face[i]
@@ -50,16 +50,28 @@ def _compute_face_node_angles_convex(
             xhere = node_x[ihere]
             yhere = node_y[ihere]
             zhere = node_z[ihere]
-            v1 = np.array(
-                [node_x[iprev] - xhere, node_y[iprev] - yhere, node_z[iprev] - zhere]
+            v1 = (node_x[iprev] - xhere, node_y[iprev] - yhere, node_z[iprev] - zhere)
+            v2 = (node_x[inext] - xhere, node_y[inext] - yhere, node_z[inext] - zhere)
+            # Spherical geometry: project onto tangent plane at the current node
+            normal = (xhere, yhere, zhere)
+            normal_norm = _numba_norm3(normal)  # |normal|
+            normal = (
+                normal[0] / normal_norm,
+                normal[1] / normal_norm,
+                normal[2] / normal_norm,
             )
-            v2 = np.array(
-                [node_x[inext] - xhere, node_y[inext] - yhere, node_z[inext] - zhere]
+            # v1 -= np.dot(v1, normal) * normal
+            v1_dot_normal = v1[0] * normal[0] + v1[1] * normal[1] + v1[2] * normal[2]
+            v2_dot_normal = v2[0] * normal[0] + v2[1] * normal[1] + v2[2] * normal[2]
+            v1 = (
+                v1[0] - v1_dot_normal * normal[0],
+                v1[1] - v1_dot_normal * normal[1],
+                v1[2] - v1_dot_normal * normal[2],
+            )
+            v2 = (
+                v2[0] - v2_dot_normal * normal[0],
+                v2[1] - v2_dot_normal * normal[1],
+                v2[2] - v2_dot_normal * normal[2],
             )
-            # Spherical geometry: project onto tangent plane at the current node
-            normal = np.array([xhere, yhere, zhere])
-            normal /= np.linalg.norm(normal)
-            v1 -= np.dot(v1, normal) * normal
-            v2 -= np.dot(v2, normal) * normal
             result[i, j] = _small_angle_of_2_vectors(v1, v2)
     return result
diff --git a/uxarray/grid/grid.py b/uxarray/grid/grid.py
index 5a77f23d6..ea18f4cf9 100644
--- a/uxarray/grid/grid.py
+++ b/uxarray/grid/grid.py
@@ -1964,9 +1964,6 @@ def compute_face_node_angles(
         ----------
         degrees : bool, defaults to False
             Whether to return angles in degrees (if True) or radians (if False).
-        assume_convex : bool, defaults to False
-            Whether to assume that all faces are convex, i.e. all internal angles less than 180 degrees.
-            If True, uses a more efficient algorithm that will produce incorrect results for non-convex faces.
         cache : None or bool, defaults to None
             Whether to cache the computed angles in the grid's dataset, in face_node_angles data variable.
             Cached angles are always in radians.
@@ -2011,7 +2008,6 @@ def compute_face_node_angles(
             result = np.rad2deg(result)
         if as_uxarray:
             from uxarray.core.dataarray import UxDataArray
-            # ^import at runtime to avoid circular import
 
             result = UxDataArray(result, uxgrid=self)
         return result
diff --git a/uxarray/grid/utils.py b/uxarray/grid/utils.py
index d70962273..2a63f48d5 100644
--- a/uxarray/grid/utils.py
+++ b/uxarray/grid/utils.py
@@ -12,9 +12,9 @@ def _small_angle_of_2_vectors(u, v):
 
     Parameters
     ----------
-    u : numpy.ndarray
+    u : numpy.ndarray or iterable of length 3
         The first 3D vector.
-    v : numpy.ndarray
+    v : numpy.ndarray or iterable of length 3
         The second 3D vector.
 
     Returns
@@ -22,14 +22,53 @@ def _small_angle_of_2_vectors(u, v):
     float
         The smallest angle between `u` and `v` in radians.
     """
-    v_norm_times_u = np.linalg.norm(v) * u
-    u_norm_times_v = np.linalg.norm(u) * v
-    vec_minus = v_norm_times_u - u_norm_times_v
-    vec_sum = v_norm_times_u + u_norm_times_v
+    # don't convert to numpy array if not already numpy array
+    # --> use u[0], u[1], u[2] instead of np.linalg.norm.
+    # TODO: replace with a numba norm function, instead of repeating the logic inline here.
+    # The formula is: angle = 2 * arctan2(| |v|*u - |u|*v |, | |v|*u + |u|*v |)
+    v_norm = _numba_norm3(v)
+    u_norm = _numba_norm3(u)
+    v_norm_times_u = (v_norm * u[0], v_norm * u[1], v_norm * u[2])
+    u_norm_times_v = (u_norm * v[0], u_norm * v[1], u_norm * v[2])
+    vec_minus = (
+        v_norm_times_u[0] - u_norm_times_v[0],
+        v_norm_times_u[1] - u_norm_times_v[1],
+        v_norm_times_u[2] - u_norm_times_v[2],
+    )
+    vec_sum = (
+        v_norm_times_u[0] + u_norm_times_v[0],
+        v_norm_times_u[1] + u_norm_times_v[1],
+        v_norm_times_u[2] + u_norm_times_v[2],
+    )
     angle_u_v_rad = 2 * np.arctan2(np.linalg.norm(vec_minus), np.linalg.norm(vec_sum))
     return angle_u_v_rad
 
 
+# TODO: move _numba_norm3 to a higher-level utils file. For more details, see issue #1648.
+@njit(cache=True)
+def _numba_norm3(u):
+    """
+    Compute the Euclidean norm of a 3D vector.
+    Implementation is currently equivalent to np.linalg.norm:
+        sqrt(u[0]**2 + u[1]**2 + u[2]**2)
+
+    Does NOT internally convert u to a list or numpy array;
+    utilizing tuples in numba instead of many tiny lists/arrays
+    can improve performance significantly.
+
+    Parameters
+    ----------
+    u : iterable of length 3, possibly a numpy array
+        The 3D vector.
+
+    Returns
+    -------
+    float
+        The Euclidean norm of the vector `u`.
+    """
+    return (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 0.5
+
+
 @njit(cache=True)
 def _angle_of_2_vectors(u, v):
     """

From 09bafd7c152223e3423d6b9a72d549f503f01be6 Mon Sep 17 00:00:00 2001
From: Sam Evans <47793072+Sevans711@users.noreply.github.com>
Date: Tue, 11 Aug 2026 16:11:31 -0400
Subject: [PATCH 4/8] face_node_angles fix forgot an optimization

(it was making a type error in numba when trying to run the method)
---
 uxarray/grid/utils.py | 4 +++-
 1 file changed, 3 insertions(+), 1 deletion(-)

diff --git a/uxarray/grid/utils.py b/uxarray/grid/utils.py
index 2a63f48d5..b7c79eec2 100644
--- a/uxarray/grid/utils.py
+++ b/uxarray/grid/utils.py
@@ -40,7 +40,9 @@ def _small_angle_of_2_vectors(u, v):
         v_norm_times_u[1] + u_norm_times_v[1],
         v_norm_times_u[2] + u_norm_times_v[2],
     )
-    angle_u_v_rad = 2 * np.arctan2(np.linalg.norm(vec_minus), np.linalg.norm(vec_sum))
+    norm_vec_minus = _numba_norm3(vec_minus)
+    norm_vec_sum = _numba_norm3(vec_sum)
+    angle_u_v_rad = 2 * np.arctan2(norm_vec_minus, norm_vec_sum)
     return angle_u_v_rad
 
 

From 1e3932e1064ac835e2027a82447ab44c1c4ae546 Mon Sep 17 00:00:00 2001
From: Sam Evans <47793072+Sevans711@users.noreply.github.com>
Date: Wed, 12 Aug 2026 10:36:50 -0400
Subject: [PATCH 5/8] remove face_node_angles cache

compute_face_node_angles() seems to be really fast, so including a cache option in initial implementation isn't worthwhile.
---
 uxarray/grid/grid.py | 42 +++++++++++++-----------------------------
 1 file changed, 13 insertions(+), 29 deletions(-)

diff --git a/uxarray/grid/grid.py b/uxarray/grid/grid.py
index ea18f4cf9..43c9575be 100644
--- a/uxarray/grid/grid.py
+++ b/uxarray/grid/grid.py
@@ -1954,7 +1954,6 @@ def compute_face_node_angles(
         self,
         *,
         degrees: bool = False,
-        cache: bool | None = None,
         as_uxarray: bool = False,
     ) -> xr.DataArray | UxDataArray:
         """Compute the angles at each node of each face in the grid.
@@ -1964,12 +1963,6 @@ def compute_face_node_angles(
         ----------
         degrees : bool, defaults to False
             Whether to return angles in degrees (if True) or radians (if False).
-        cache : None or bool, defaults to None
-            Whether to cache the computed angles in the grid's dataset, in face_node_angles data variable.
-            Cached angles are always in radians.
-            If None, use cached result if available, else compute but do not cache result.
-            If True, use cached result if available, else compute and cache result.
-            If False, always recompute; do not check or store in cache.
         as_uxarray : bool, defaults to False
             Whether to return a uxarray.DataArray (if True) instead of an xarray.DataArray (if False).
             If True, equivalent to uxarray.DataArray(self.compute_face_node_angles(..., as_uxarray=False), uxgrid=self).
@@ -1982,28 +1975,19 @@ def compute_face_node_angles(
         """
         from uxarray.conventions.ugrid import FACE_DIM, N_MAX_FACE_NODES_DIM
 
-        name = "face_node_angles"
-        result = None
-        if cache is None or cache:
-            if name in self._ds:
-                result = self._ds[name]
-        if result is None:  # result not in cache; need to compute it
-            result = _compute_face_node_angles_convex(
-                self.node_x.values,
-                self.node_y.values,
-                self.node_z.values,
-                self.face_node_connectivity.values,
-                self.n_nodes_per_face.values,
-            )
-            # convert to xr.DataArray
-            result = xr.DataArray(
-                data=result,
-                dims=[FACE_DIM, N_MAX_FACE_NODES_DIM],
-                name=name,
-                attrs={"description": "Internal angles at each node of each face."},
-            )
-            if cache:
-                self._ds[name] = result
+        result = _compute_face_node_angles_convex(
+            self.node_x.values,
+            self.node_y.values,
+            self.node_z.values,
+            self.face_node_connectivity.values,
+            self.n_nodes_per_face.values,
+        )
+        result = xr.DataArray(
+            data=result,
+            dims=[FACE_DIM, N_MAX_FACE_NODES_DIM],
+            name="face_node_angles",
+            attrs={"description": "Internal angles at each node of each face."},
+        )
         if degrees:
             result = np.rad2deg(result)
         if as_uxarray:

From 7d610ba0b826dc9dd7e4ff4218527fae9f7052f6 Mon Sep 17 00:00:00 2001
From: Sam Evans <47793072+Sevans711@users.noreply.github.com>
Date: Wed, 12 Aug 2026 10:37:13 -0400
Subject: [PATCH 6/8] add face_node_angles tests

---
 test/grid/geometry/test_angles.py | 67 +++++++++++++++++++++++++++++++
 1 file changed, 67 insertions(+)
 create mode 100644 test/grid/geometry/test_angles.py

diff --git a/test/grid/geometry/test_angles.py b/test/grid/geometry/test_angles.py
new file mode 100644
index 000000000..f818ac576
--- /dev/null
+++ b/test/grid/geometry/test_angles.py
@@ -0,0 +1,67 @@
+"""
+Purpose: tests related to angle calculations on a grid
+"""
+
+import numpy as np
+import xarray as xr
+
+import uxarray as ux
+
+
+def test_face_node_angles_triangle():
+    """ensure Grid.compute_face_node_angles() works as expected for a simple ~30,60,90 triangle."""
+    # make a tiny triangle with known angles (90,60,30 degrees):
+    # (n1)
+    #  |   %%
+    #  |       %%
+    # (n0) ------ (n2)
+    node_lon = [0, 0, np.sqrt(3)]
+    node_lat = [0, 1, 0]
+    face_node_connectivity = [[0, 1, 2]]
+    grid = ux.Grid.from_topology(node_lon, node_lat, face_node_connectivity)
+    angles_rad = grid.compute_face_node_angles()
+    angles_deg = grid.compute_face_node_angles(degrees=True)
+    assert np.allclose(np.rad2deg(angles_rad), angles_deg)
+    angles_uxarr = grid.compute_face_node_angles(as_uxarray=True)
+    assert isinstance(angles_rad, xr.DataArray)
+    assert isinstance(angles_uxarr, ux.UxDataArray)
+    assert np.all(angles_rad == angles_uxarr)
+    angle_at_n0 = angles_deg.isel(n_face=0, n_max_face_nodes=0)
+    angle_at_n1 = angles_deg.isel(n_face=0, n_max_face_nodes=1)
+    angle_at_n2 = angles_deg.isel(n_face=0, n_max_face_nodes=2)
+    assert angle_at_n0 == 90  # this turns out to be exact... from arctan2(any_value, 0)
+    assert np.abs(angle_at_n1 - 60) < 1e-2   # basically 60 degrees
+    assert np.abs(angle_at_n2 - 30) < 1e-2   # basically 30 degrees
+    # on a unit sphere, spherical excess == face area, via Girard's theorem.
+    spherical_excess = angles_rad.sum('n_max_face_nodes') - np.pi
+    face_areas = grid.compute_face_areas()
+    assert np.allclose(spherical_excess, face_areas, atol=0, rtol=1e-12)
+
+def test_face_node_angles_hexagons_and_pentagons():
+    """ensure face_node_angles works as expected on grids with hexagons and pentagons"""
+    grid = ux.tutorial.open_grid('quad-hexagon')  # has multiple faces, all hexagons.
+    angles_deg = grid.compute_face_node_angles(degrees=True)
+    # every hexagon in this grid is close to regular (all 120 degree angles):
+    regular_hex_deviation = angles_deg - 120
+    assert np.max(np.abs(regular_hex_deviation)) < 4.0
+    # generalized spherical excess formula uses (n - 2) * np.pi; n==6 for all of these faces
+    angles = grid.compute_face_node_angles()  # (need to use radians for this formula)
+    spherical_excess = angles.sum('n_max_face_nodes') - (6 - 2) * np.pi
+    face_areas = grid.compute_face_areas()
+    assert np.allclose(spherical_excess, face_areas, atol=0, rtol=1e-10)
+
+    # now test a grid which has pentagons too,
+    # to ensure the implementation works even when the number of nodes per face varies.
+    grid = ux.tutorial.open_grid('mpas-QU-480')
+    # not all close to regular so don't try to check that.
+    # ensure nan values wherever n_max_face_nodes dimension is larger than n_nodes_per_face
+    angles = grid.compute_face_node_angles()
+    assert not np.all(grid.n_nodes_per_face == grid.n_max_face_nodes)
+    should_have_nans = angles.where(grid.n_nodes_per_face < grid.n_max_face_nodes, drop=True)
+    assert should_have_nans.size > 0
+    should_be_nans = should_have_nans.isel(n_max_face_nodes = -1)
+    assert np.all(np.isnan(should_be_nans))
+    # generalized spherical excess formula uses (n_nodes_per_face - 2) * np.pi
+    spherical_excess = angles.sum('n_max_face_nodes') - (grid.n_nodes_per_face - 2) * np.pi
+    face_areas = grid.compute_face_areas()
+    assert np.allclose(spherical_excess, face_areas, atol=0, rtol=1e-9)

From 74b1660036f2a1d9c46a077c5399c1600a9d9b79 Mon Sep 17 00:00:00 2001
From: Sam Evans <47793072+Sevans711@users.noreply.github.com>
Date: Wed, 12 Aug 2026 10:38:28 -0400
Subject: [PATCH 7/8] docs: add compute_face_node_angles to api.rst

---
 docs/api.rst | 1 +
 1 file changed, 1 insertion(+)

diff --git a/docs/api.rst b/docs/api.rst
index 2addaabd8..bc92a7513 100644
--- a/docs/api.rst
+++ b/docs/api.rst
@@ -192,6 +192,7 @@ Methods
    Grid.copy
    Grid.calculate_total_face_area
    Grid.compute_face_areas
+   Grid.compute_face_node_angles
    Grid.construct_face_centers
    Grid.get_ball_tree
    Grid.get_kd_tree

From a6fa5bb5632786e4bf1b6825351e83b556ac4ffc Mon Sep 17 00:00:00 2001
From: Sam Evans <47793072+Sevans711@users.noreply.github.com>
Date: Wed, 12 Aug 2026 11:08:03 -0400
Subject: [PATCH 8/8] clean up comments related to face_node_angles

---
 uxarray/grid/grid.py  | 1 +
 uxarray/grid/utils.py | 4 +---
 2 files changed, 2 insertions(+), 3 deletions(-)

diff --git a/uxarray/grid/grid.py b/uxarray/grid/grid.py
index 43c9575be..31338cf66 100644
--- a/uxarray/grid/grid.py
+++ b/uxarray/grid/grid.py
@@ -1972,6 +1972,7 @@ def compute_face_node_angles(
         face_node_angles : xr.DataArray or uxarray.UxDataArray (if as_uxarray=True)
             The internal angles at each node, for each face in the grid.
             Has 'n_face' and 'n_max_face_nodes' dimensions, with same size as in self.
+            For faces with fewer than n_max_face_nodes, fill value is np.nan.
         """
         from uxarray.conventions.ugrid import FACE_DIM, N_MAX_FACE_NODES_DIM
 
diff --git a/uxarray/grid/utils.py b/uxarray/grid/utils.py
index b7c79eec2..6169d5892 100644
--- a/uxarray/grid/utils.py
+++ b/uxarray/grid/utils.py
@@ -22,9 +22,7 @@ def _small_angle_of_2_vectors(u, v):
     float
         The smallest angle between `u` and `v` in radians.
     """
-    # don't convert to numpy array if not already numpy array
-    # --> use u[0], u[1], u[2] instead of np.linalg.norm.
-    # TODO: replace with a numba norm function, instead of repeating the logic inline here.
+    # don't convert to numpy array if not already numpy array.
     # The formula is: angle = 2 * arctan2(| |v|*u - |u|*v |, | |v|*u + |u|*v |)
     v_norm = _numba_norm3(v)
     u_norm = _numba_norm3(u)