From 470262c7e6087557df9279eb4242603b3ab06ac8 Mon Sep 17 00:00:00 2001 From: Pablo Brubeck Date: Wed, 4 Mar 2026 10:40:19 +0000 Subject: [PATCH 1/2] Hypergeometric2F1 --- ufl/__init__.py | 2 ++ ufl/algorithms/apply_derivatives.py | 15 ++++++++++++ ufl/mathfunctions.py | 36 +++++++++++++++++++++++++++++ ufl/operators.py | 10 ++++++++ 4 files changed, 63 insertions(+) diff --git a/ufl/__init__.py b/ufl/__init__.py index 415914d77..c6fc6ec8b 100644 --- a/ufl/__init__.py +++ b/ufl/__init__.py @@ -388,6 +388,7 @@ ge, grad, gt, + hyp2f1, imag, inner, inv, @@ -608,6 +609,7 @@ "grad", "gt", "hexahedron", + "hyp2f1", "i", "identity_pullback", "imag", diff --git a/ufl/algorithms/apply_derivatives.py b/ufl/algorithms/apply_derivatives.py index 1194ace3a..4747e66e0 100644 --- a/ufl/algorithms/apply_derivatives.py +++ b/ufl/algorithms/apply_derivatives.py @@ -90,6 +90,7 @@ Cosh, Erf, Exp, + Hypergeometric2F1, Ln, MathFunction, Sin, @@ -112,6 +113,7 @@ cosh, exp, facet_avg, + hyp2f1, ln, sign, sin, @@ -583,6 +585,19 @@ def _(self, o: Expr, fp: Expr) -> Expr: (f,) = o.ufl_operands return fp * (2.0 / sqrt(pi) * exp(-(f**2))) + @process.register(Hypergeometric2F1) + @DAGTraverser.postorder + def _(self, o: Expr, ap: Expr, bp: Expr, cp: Expr, fp: Expr) -> Expr: + """Differentiate a hyp2f1.""" + a, b, c, f = o.ufl_operands + if not all(nup is None or isinstance(nup, Zero) for nup in (ap, bp, cp)): + raise NotImplementedError( + "Differentiation of hypergeometric function w.r.t. (a, b, c) is not supported." + ) + + op = (a * b / c) * hyp2f1(a + 1, b + 1, c + 1, f) + return op * fp + # --- Bessel functions @process.register(BesselJ) diff --git a/ufl/mathfunctions.py b/ufl/mathfunctions.py index 3caba472e..1578c5f40 100644 --- a/ufl/mathfunctions.py +++ b/ufl/mathfunctions.py @@ -388,6 +388,42 @@ def evaluate(self, x, mapping, component, index_values): return math.erf(a) +@ufl_type(is_scalar=True, num_ops=4) +class Hypergeometric2F1(Operator): + """Class for hypergeometric 2F1 functions.""" + + __slots__ = "_name" + + def __init__(self, a, b, c, argument): + """Initialise.""" + if not all(is_true_ufl_scalar(f) for f in (a, b, c)): + raise ValueError("Expecting scalar (a, b, c).") + if not is_true_ufl_scalar(argument): + raise ValueError("Expecting scalar argument.") + + Operator.__init__(self, (a, b, c, argument)) + + self._name = "hyp2f1" + + def evaluate(self, x, mapping, component, index_values): + """Evaluate.""" + try: + import scipy.special # type: ignore + except ImportError: + raise ValueError( + "You must have scipy installed to evaluate hypergeometric functions in python." + ) + func = getattr(scipy.special, self._name) + return func( + *(arg.evaluate(x, mapping, component, index_values) for arg in self.ufl_operands) + ) + + def __str__(self): + """Format as a string.""" + a, b, c, arg = self.ufl_operands + return f"{self._name}({a}, {b}, {c}, {arg})" + + @ufl_type(is_abstract=True, is_scalar=True, num_ops=2) class BesselFunction(Operator): """Base class for all bessel functions.""" diff --git a/ufl/operators.py b/ufl/operators.py index ccbeaac98..a32b2b895 100644 --- a/ufl/operators.py +++ b/ufl/operators.py @@ -51,6 +51,7 @@ Cosh, Erf, Exp, + Hypergeometric2F1, Ln, Sin, Sinh, @@ -700,6 +701,15 @@ def erf(f): return _mathfunction(f, Erf) +def hyp2f1(a, b, c, f): + """The Gaussian hypergeometric function 2F1.""" + a = as_ufl(a) + b = as_ufl(b) + c = as_ufl(c) + f = as_ufl(f) + return Hypergeometric2F1(a, b, c, f) + + def bessel_J(nu, f): """Cylindrical Bessel function of the first kind.""" nu = as_ufl(nu) From d19ea23940194cfaf893976af0ddb4bf14373403 Mon Sep 17 00:00:00 2001 From: Pablo Brubeck Date: Thu, 5 Mar 2026 10:15:58 +0000 Subject: [PATCH 2/2] elliptic intergrals --- ufl/__init__.py | 4 ++++ ufl/mathfunctions.py | 2 +- ufl/operators.py | 21 +++++++++++++++++++++ 3 files changed, 26 insertions(+), 1 deletion(-) diff --git a/ufl/__init__.py b/ufl/__init__.py index c6fc6ec8b..692b1c597 100644 --- a/ufl/__init__.py +++ b/ufl/__init__.py @@ -380,6 +380,8 @@ elem_mult, elem_op, elem_pow, + elliptic_E, + elliptic_K, eq, erf, exp, @@ -596,6 +598,8 @@ "elem_mult", "elem_op", "elem_pow", + "elliptic_E", + "elliptic_K", "energy_norm", "eq", "erf", diff --git a/ufl/mathfunctions.py b/ufl/mathfunctions.py index 1578c5f40..a3edd90dc 100644 --- a/ufl/mathfunctions.py +++ b/ufl/mathfunctions.py @@ -390,7 +390,7 @@ def evaluate(self, x, mapping, component, index_values): @ufl_type(is_scalar=True, num_ops=4) class Hypergeometric2F1(Operator): - """Class for hypergeometric 2F1 functions.""" + """The Gaussian hypergeometric function 2F1.""" __slots__ = "_name" diff --git a/ufl/operators.py b/ufl/operators.py index a32b2b895..28206fb2a 100644 --- a/ufl/operators.py +++ b/ufl/operators.py @@ -16,6 +16,7 @@ import operator import warnings +from math import pi from ufl import sobolevspace from ufl.algebra import Conj, Imag, Real @@ -710,6 +711,26 @@ def hyp2f1(a, b, c, f): return Hypergeometric2F1(a, b, c, f) +def elliptic_K(k): + r"""Complete elliptic integral of the first kind. + + This function is defined as + + .. math:: K(k) = \int_0^{\pi/2} [1 - k^2 sin(t)^2]^{-1/2} dt + """ + return (pi / 2) * hyp2f1(1 / 2, 1 / 2, 1, k**2) + + +def elliptic_E(k): + r"""Complete elliptic integral of the second kind. + + This function is defined as + + .. math:: E(k) = \int_0^{\pi/2} [1 - k^2 sin(t)^2]^{1/2} dt + """ + return (pi / 2) * hyp2f1(1 / 2, -1 / 2, 1, k**2) + + def bessel_J(nu, f): """Cylindrical Bessel function of the first kind.""" nu = as_ufl(nu)