feat(math): add polynomial vector calculus operations (#1865)

Add polynomial_field domain with 5 operations: gradient, Laplacian,
directional derivative, divergence, and curl. Uses SymPy for exact
symbolic differentiation.

Closes #1865

Author: morluto

src/jacobian/catalog/builtins.py

@@ -218,6 +218,12 @@
 from jacobian.math.optimization._tools import TOOLS as OPTIMIZATION_TOOLS
 from jacobian.math.petri_nets._admission import ADMISSIONS as PETRI_NETS_ADMISSIONS
 from jacobian.math.petri_nets._tools import TOOLS as PETRI_NET_TOOLS
+from jacobian.math.polynomial_vector_calc._admission import (
+    ADMISSIONS as POLYNOMIAL_VECTOR_CALC_ADMISSIONS,
+)
+from jacobian.math.polynomial_vector_calc._tools import (
+    TOOLS as POLYNOMIAL_VECTOR_CALC_TOOLS,
+)
 from jacobian.math.polynomials._admission import ADMISSIONS as POLYNOMIALS_ADMISSIONS
 from jacobian.math.polynomials._tools import TOOLS as POLYNOMIAL_TOOLS
 from jacobian.math.polynomials.maps._admission import (
@@ -309,6 +315,7 @@
     *LATTICES_TOOLS,
     *FORMAL_POWER_SERIES_TOOLS,
     *POLYNOMIAL_TOOLS,
+    *POLYNOMIAL_VECTOR_CALC_TOOLS,
     *MULTIVARIATE_POLYNOMIAL_TOOLS,
     *ANALYSIS_TOOLS,
     *PROBABILITY_TOOLS,
@@ -406,6 +413,7 @@
     *OPTIMIZATION_ADMISSIONS,
     *PETRI_NETS_ADMISSIONS,
     *POLYNOMIALS_ADMISSIONS,
+    *POLYNOMIAL_VECTOR_CALC_ADMISSIONS,
     *POLYNOMIALS_MAPS_ADMISSIONS,
     *POLYNOMIALS_MULTIVARIATE_ADMISSIONS,
     *POLYNOMIALS_REAL_ALGEBRA_ADMISSIONS,

src/jacobian/math/polynomial_vector_calc/__init__.py

@@ -0,0 +1,3 @@
+"""Polynomial vector calculus operations."""
+
+__all__: list[str] = []

src/jacobian/math/polynomial_vector_calc/_admission.py

@@ -0,0 +1,33 @@
+"""Owner-local admission decisions for built-in math operations."""
+
+from __future__ import annotations
+
+from jacobian.catalog.admission import AdmissionDecision, OperationAdmission
+
+ADMISSIONS: tuple[OperationAdmission, ...] = (
+    OperationAdmission(
+        "polynomial_field.scalar.gradient.compute",
+        AdmissionDecision.KEEP,
+        "distinct exact bounded mathematical value or invariant with material computational or reliability leverage",
+    ),
+    OperationAdmission(
+        "polynomial_field.scalar.laplacian.compute",
+        AdmissionDecision.KEEP,
+        "distinct exact bounded mathematical value or invariant with material computational or reliability leverage",
+    ),
+    OperationAdmission(
+        "polynomial_field.scalar.directional_derivative.compute",
+        AdmissionDecision.KEEP,
+        "distinct exact bounded mathematical value or invariant with material computational or reliability leverage",
+    ),
+    OperationAdmission(
+        "polynomial_field.vector.divergence.compute",
+        AdmissionDecision.KEEP,
+        "distinct exact bounded mathematical value or invariant with material computational or reliability leverage",
+    ),
+    OperationAdmission(
+        "polynomial_field.vector.curl.compute",
+        AdmissionDecision.KEEP,
+        "distinct exact bounded mathematical value or invariant with material computational or reliability leverage",
+    ),
+)

src/jacobian/math/polynomial_vector_calc/_models.py

@@ -0,0 +1,75 @@
+"""Typed wire contracts for polynomial vector calculus operations."""
+
+from __future__ import annotations
+
+from typing import Self
+
+from pydantic import Field, model_validator
+
+from jacobian._models import StrictModel
+
+MAX_VARS = 8
+MAX_POLYS = 8
+
+
+class ScalarFieldRequest(StrictModel):
+    """A multivariate polynomial scalar field."""
+
+    variables: tuple[str, ...] = Field(min_length=1, max_length=MAX_VARS)
+    polynomial: str = Field(min_length=1, max_length=4096)
+
+    @model_validator(mode="after")
+    def require_valid(self) -> Self:
+        if any(not v.isidentifier() for v in self.variables):
+            raise ValueError("variable names must be valid identifiers")
+        return self
+
+
+class VectorFieldRequest(StrictModel):
+    """A multivariate polynomial vector field."""
+
+    variables: tuple[str, ...] = Field(min_length=1, max_length=MAX_VARS)
+    components: tuple[str, ...] = Field(
+        min_length=1, max_length=MAX_POLYS
+    )
+
+    @model_validator(mode="after")
+    def require_valid(self) -> Self:
+        if any(not v.isidentifier() for v in self.variables):
+            raise ValueError("variable names must be valid identifiers")
+        return self
+
+
+class DirectionalDerivativeRequest(StrictModel):
+    """Directional derivative of a scalar field along a direction vector."""
+
+    variables: tuple[str, ...] = Field(min_length=1, max_length=MAX_VARS)
+    polynomial: str = Field(min_length=1, max_length=4096)
+    direction: tuple[str, ...] = Field(
+        min_length=1, max_length=MAX_POLYS
+    )
+
+    @model_validator(mode="after")
+    def require_valid(self) -> Self:
+        if len(self.direction) != len(self.variables):
+            raise ValueError("direction vector length must match variables length")
+        return self
+
+
+# Results
+
+
+class ScalarResult(StrictModel):
+    """A scalar polynomial result."""
+
+    result: str
+    variables: tuple[str, ...] = ()
+    method: str = "SYMPY_GRADIENT"
+
+
+class VectorResult(StrictModel):
+    """A vector polynomial result."""
+
+    components: tuple[str, ...]
+    variables: tuple[str, ...] = ()
+    method: str = "SYMPY_GRADIENT"

src/jacobian/math/polynomial_vector_calc/_operations.py

@@ -0,0 +1,99 @@
+"""Domain functions for polynomial vector calculus operations."""
+
+from __future__ import annotations
+
+import sympy
+
+from jacobian.math.polynomial_vector_calc._models import (
+    DirectionalDerivativeRequest,
+    ScalarFieldRequest,
+    ScalarResult,
+    VectorFieldRequest,
+    VectorResult,
+)
+
+
+def _parse_poly(expr_str: str, variables: tuple[str, ...]) -> sympy.Expr:
+    """Parse a polynomial expression string with given variable names."""
+    var_symbols = sympy.symbols(variables)
+    if len(variables) == 1:
+        var_symbols = (var_symbols,)
+    return sympy.sympify(expr_str, locals=dict(zip(variables, var_symbols, strict=True)))
+
+
+def compute_gradient(request: ScalarFieldRequest) -> VectorResult:
+    var_symbols = sympy.symbols(variables := request.variables)
+    if len(variables) == 1:
+        var_symbols = (var_symbols,)
+    poly = _parse_poly(request.polynomial, request.variables)
+    grads = [sympy.diff(poly, v) for v in var_symbols]
+    return VectorResult(
+        components=tuple(str(g) for g in grads),
+        variables=request.variables,
+        method="SYMPY_GRADIENT",
+    )
+
+
+def compute_divergence(request: VectorFieldRequest) -> ScalarResult:
+    var_symbols = sympy.symbols(request.variables)
+    if len(request.variables) == 1:
+        var_symbols = (var_symbols,)
+    polys = [_parse_poly(c, request.variables) for c in request.components]
+    div = sum(
+        sympy.diff(p, v) for p, v in zip(polys, var_symbols, strict=True)
+    )
+    return ScalarResult(
+        result=str(sympy.expand(div)),
+        variables=request.variables,
+        method="SYMPY_DIVERGENCE",
+    )
+
+
+def compute_curl(request: VectorFieldRequest) -> VectorResult:
+    """Curl in 3D: (d_z F_y - d_y F_z, d_x F_z - d_z F_x, d_y F_x - d_x F_y)."""
+    if len(request.variables) != 3:
+        raise ValueError("curl is defined for 3D vector fields")
+    x, y, z = sympy.symbols(request.variables)
+    fx, fy, fz = (sympy.sympify(c) for c in [
+        _parse_poly(request.components[0], request.variables),
+        _parse_poly(request.components[1], request.variables),
+        _parse_poly(request.components[2], request.variables),
+    ])
+    curl_x = sympy.diff(fz, y) - sympy.diff(fy, z)
+    curl_y = sympy.diff(fx, z) - sympy.diff(fz, x)
+    curl_z = sympy.diff(fy, x) - sympy.diff(fx, y)
+    return VectorResult(
+        components=(str(sympy.expand(curl_x)), str(sympy.expand(curl_y)), str(sympy.expand(curl_z))),
+        variables=request.variables,
+        method="SYMPY_CURL",
+    )
+
+
+def compute_laplacian(request: ScalarFieldRequest) -> ScalarResult:
+    var_symbols = sympy.symbols(request.variables)
+    if len(request.variables) == 1:
+        var_symbols = (var_symbols,)
+    poly = _parse_poly(request.polynomial, request.variables)
+    laplacian = sum(sympy.diff(poly, v, 2) for v in var_symbols)
+    return ScalarResult(
+        result=str(sympy.expand(laplacian)),
+        variables=request.variables,
+        method="SYMPY_LAPLACIAN",
+    )
+
+
+def compute_directional_derivative(
+    request: DirectionalDerivativeRequest,
+) -> ScalarResult:
+    var_symbols = sympy.symbols(request.variables)
+    if len(request.variables) == 1:
+        var_symbols = (var_symbols,)
+    poly = _parse_poly(request.polynomial, request.variables)
+    grad = [sympy.diff(poly, v) for v in var_symbols]
+    direction = [sympy.sympify(d) for d in request.direction]
+    result = sum(g * d for g, d in zip(grad, direction, strict=True))
+    return ScalarResult(
+        result=str(sympy.expand(result)),
+        variables=request.variables,
+        method="SYMPY_DIRECTIONAL_DERIVATIVE",
+    )