🇬🇧 English · 🇷🇺 Русский

formula - Arbitrary-precision formula parser and solver

Development status

Development plan:

• Deploy new version with complex numbers support. (v4.0)
• Deploy prebuilt wheels for Linux, macOS (Intel + Apple Silicon), and Windows via cibuildwheel.
• Add per-function regression tests covering all built-in math functions (see tests/functions/).
• Support sign(), abs() and a Pythonic Number wrapper class.
• Complex numbers tests and examples. (Character i is reserved for this by default.)
• Formatting documentation. Limits and caveats.
• Support rand for random number in [0.0 - 1.0].

This project is built with pybind11 and boost.

Installation

pip install formula

Prebuilt wheels are published to PyPI for Python 3.11 – 3.13 on:

• Linux — manylinux_2_28 x86_64 (glibc ≥ 2.28)
• macOS — x86_64 and arm64 (Apple Silicon)
• Windows — AMD64 and x86

If a wheel is unavailable for your platform, pip will fall back to building from the source distribution; in that case a working C++ toolchain is required. See Development → Building from source.

On Windows

• pip install formula

• if you get an `error: Microsoft Visual C++ 14.0 is required.`

  Install Microsoft Visual C++ Build Tools 14.0 from https://visualstudio.microsoft.com/visual-cpp-build-tools/ and try again. Example of the correct selection to install:

  Windows runtime requirements

  On Windows, the Visual C++ 2015 redistributable packages are a runtime requirement for this project. It can be found here.

  If you use the Anaconda python distribution, you may require the Visual Studio runtime as a platform-dependent runtime requirement for you package:

  requirements:
  build:
  - python
  - setuptools
  - pybind11run:
  - python
  - vs2015_runtime # [win]

Documentation

formula contains case sensitive (by default) string parser. Let's imagine that we have a string expression: "(x^2+y)/sin(a*pi)":

>>>fromformulaimportSolver>>>formula=Solver("(x^2+y)/sin(a*pi)", precision=32)

Then we want to calculate the value of this function in the following point:

>>>point= {"x": "3", "y": "3e-20", "a": "-0.5"}

And it is enough to call the formula object to calculate the value of the expression or the derivative of the expression at this point:

>>>formula(point) # (3^2 + 3e-20)/sin(-pi/2)'-9.00000000000000000003'>>>formula(point, derivative="x") # 2*3/sin(-pi/2)'-6'>>>formula(point, derivative=("y", "a")) # [1/sin(-pi/2),- (3^2 + 3e-20) * cos(-pi/2) / sin(-pi/2)]
['-1', '-1.5633175729821453046351823925394e-47'] # cos(-pi/2) is an epsilon, not exact 0

Simple examples

One plus one =)

>>>fromformulaimportSolver, FmtFlags>>>Solver("1+1", precision=32)()
'2'>>>Solver("1+1", 32)(format_digits=20, format_flags=FmtFlags.showpos)
'+2'>>>Solver("1+1", 32)(format_digits=20, format_flags=FmtFlags.fixed|FmtFlags.showpos)
'+2.00000000000000000000'>>>Solver("1+1", 32)(format_digits=20, format_flags=FmtFlags.scientific|FmtFlags.showpos)
'+2.00000000000000000000e+00'

Find the number of PI using arcsin

Precision = 32

>>>fromformulaimportSolver, FmtFlags>>>Solver("2*asin(x)", precision=32)({"x": "1"})
# just 32 digits:'3.1415926535897932384626433832795'>>>Solver("2*asin(x)", 32)({"x": "1"}, format_digits=32)
# by default format_digits is equal to precision:'3.1415926535897932384626433832795'>>>Solver("2*asin(x)", 32)({"x": "1"}, format_digits=31)
# let's round in accordance with format_digits:'3.14159265358979323846264338328'>>>Solver("2*asin(x)", 32)({"x": "1"}, format_digits=30)
'3.14159265358979323846264338328'>>>Solver("2*asin(x)", 32)({"x": "1"}, format_digits=29)
'3.1415926535897932384626433833'>>>Solver("2*asin(x)", 32)(1, format_digits=28)
'3.141592653589793238462643383'>>>Solver("2*asin(x)", 32)(1, format_digits=2)
'3.1'>>>Solver("2*asin(x)", 32)(1, format_digits=1)
'3'>>>Solver("2*asin(x)", 32)(1, format_digits=0)
# show the entire chunk of memory, including insignificant digits:'3.1415926535897932384626433832795028841971'

Precision = 4096

>>>fromformulaimportSolver, FmtFlags>>>Solver("2*asin(x)", precision=4096)(1) # 4095 digits of pi after the point ;-)
'3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117067982148086513282306647093844609550582231725359408128481117450284102701938521105559644622948954930381964428810975665933446128475648233786783165271201909145648566923460348610454326648213393607260249141273724587006606315588174881520920962829254091715364367892590360011330530548820466521384146951941511609433057270365759591953092186117381932611793105118548074462379962749567351885752724891227938183011949129833673362440656643086021394946395224737190702179860943702770539217176293176752384674818467669405132000568127145263560827785771342757789609173637178721468440901224953430146549585371050792279689258923542019956112129021960864034418159813629774771309960518707211349999998372978049951059731732816096318595024459455346908302642522308253344685035261931188171010003137838752886587533208381420617177669147303598253490428755468731159562863882353787593751957781857780532171226806613001927876611195909216420198938095257201065485863278865936153381827968230301952035301852968995773622599413891249721775283479131515574857242454150695950829533116861727855889075098381754637464939319255060400927701671139009848824012858361603563707660104710181942955596198946767837449448255379774726847104047534646208046684259069491293313677028989152104752162056966024058038150193511253382430035587640247496473263914199272604269922796782354781636009341721641219924586315030286182974555706749838505494588586926995690927210797509302955321165344987202755960236480665499119881834797753566369807426542527862551818417574672890977772793800081647060016145249192173217214772350141441973568548161361157352552133475741849468438523323907394143334547762416862518983569485562099219222184272550254256887671790494601653466804988627232791786085784383827967976681454100953883786360950680064225125205117392984896084128488626945604241965285022210661186306744278622039194945047123713786960956364371917287467764657573962413890865832645995813390478027590099465764078951269468398352595709825822620522489407726719478268482601476990902640136394437455305068203496252451749399651431429809190659250937221696461515709858387410597885959772975498930161753928468138268683868942774155991855925245953959431049972524680845987273644695848653836736222626099124608051243884390451244136549762780797715691435997700129616089441694868555848406353422072225828488648158456028506016842739452267467678895252138522549954666727823986456596116354886230577456498035593634568174324112515076069479451096596094025228879710893145669136867228748940560101503308617928680920874760917824938589009714909675985261365549781893129784821682998948722658804857564014270477555132379641451523746234364542858444795265867821051141354735739523113427166102135969536231442952484937187110145765403590279934403742007310578539062198387447808478489683321445713868751943506430218453191048481005370614680674919278191197939952061419663428754440643745123718192179998391015919561814675142691239748940907186494231961567945208095146550225231603881930142093762137855956638937787083039069792077346722182562599661501421503068038447734549202605414665925201497442850732518666002132434088190710486331734649651453905796268561005508106658796998163574736384052571459102897064140110971206280439039759515677157700420337869936007230558763176359421873125147120532928191826186125867321579198414848829164470609575270695722091756711672291098169091528017350671274858322287183520935396572512108357915136988209144421006751033467110314126711136990865851639831501970165151168517143765761835155650884909989859982387345528331635507647918535893226185489632132933089857064204675259070915481416549859461637180270981994309924488957571282890592323326097299712084433573265489382391193259746366730583604142813883032038249037589852437441702913276561809377344403070746921120191302033038019762110110044929321516084244485963766983895228684783123552658213144957685726243344189303968642624341077322697802807318915441101044682325271620105265227211166039666557309254711055785376346682065310989652691862056476931257058635662018558100729360659876486118'

Other examples

>>>fromformulaimportSolver, FmtFlags>>>Solver("9.99 + 9e-20 + 9e-51", precision=64)(None, None, 50, FmtFlags.scientific)
'9.99000000000000000009000000000000000000000000000001e+00'>>>Solver("9.99 + 9e-20 + 9e-51", 64)(None, None, 51, FmtFlags.scientific)
'9.990000000000000000090000000000000000000000000000009e+00'>>>Solver("0 + 9e-20 + 9e-51", 64)(None, None, 31, FmtFlags.scientific)
'9.0000000000000000000000000000009e-20'

Using Number as an arbitrary-precision Python value

Number is a small Pythonic wrapper around Solver that lets you compose expressions with the standard arithmetic and comparison operators while keeping arbitrary precision throughout:

>>>fromformulaimportNumber>>>a=Number("1.0000000000000000000000000001", precision=64)
>>>b=Number("2", precision=64)
>>>str(a+b)
'3.0000000000000000000000000001'>>>str(a*b)
'2.0000000000000000000000000002'>>>str(abs(Number("-3.14", precision=32)))
'3.14'>>>a>bFalse>>>a==Number("1.0000000000000000000000000001", precision=64)
True

Equality compares the (real, imag) strings from pair — exact complex identities compare equal, but expressions with drift do not:

>>>Number("i*i") ==Number("-1") # exact — no driftTrue>>>Number("i^4") ==Number("1") # drift in complex ^ → FalseFalse

Solver.pair(values, format_digits, format_flags) returns a (real_str, imag_str) tuple formatted consistently — real-only results get (value, "0") so pairs are byte-comparable. Used internally by Number.__eq__ and Number.__hash__.

Supported operators: +, -, *, /, ** (also written as ^ inside expressions), abs(), ==, <, <=, >, >=.

Finding all ray–surface intersections

RaySurface finds every intersection of a ray r(t) = O + t·d with an implicit surface F(x, y, z) = 0. Substituting the ray reduces the surface to a single-variable g(t) = F(O + t·d), so the intersections are exactly the real roots of g on [t_min, t_max]. t is measured in units of |d| (the direction is normalized internally), so for a unit d it is the distance along the ray.

>>>fromformulaimportRaySurface>>>rs=RaySurface("x*x + y*y - 1", precision=24) # unit cylinder>>>ts=rs.intersect((-2, 0, 0), (1, 0, 0), t_max=10) # two crossings: t≈1, t≈3>>>len(ts)
2>>>pts=rs.points((-2, 0, 0), (1, 0, 0), t_max=10) # (x, y, z) on the surface>>>len(pts)
2

Pick the root-finder with method=:

method	Best for	Guarantee
auto	default — picks per surface	Sturm for polynomials, else Chebyshev ∪ subdivision
sturm	algebraic surfaces (quadrics, tori, …)	exact, complete real-root count
chebyshev	smooth, low-oscillation analytic (sin/exp/…)	self-validating fit; captures tangencies
subdivision	general/oscillatory surfaces	derivative-bound exclusion (practically reliable)
sampling	quick, well-separated roots	none — may miss thin/tangent features

auto routes algebraic surfaces to sturm (exact) and other real surfaces to chebyshev reconciled with subdivision as a safety net. Even-multiplicity (tangent) roots that plain sampling steps over are recovered:

>>>rs=RaySurface("x*x + y*y + z*z - 1", precision=24) # sphere>>>len(rs.intersect((-2, 1, 0), (1, 0, 0), t_max=10, method="sturm")) # grazing1

Keyword options: t_min (default 0), plus per-method knobs such as max_degree (sturm), cheb_degree (chebyshev) and m2_samples (subdivision). Complex-valued surfaces are supported by sturm (real intersections are the common roots of Re g and Im g); the other backends are real-only.

See doc/ray-surface-intersections.md for the full guide and limits, and doc/ray-surface-design.md for the design notes. Caveats in brief:

• subdivision is practically reliable rather than formally rigorous — its exclusion test is only as good as the estimated bound on g''.
• chebyshev self-validates the fit, not the root isolation. Its Chebyshev→monomial step is ill-conditioned at high degree, so it can silently miss roots on densely oscillatory surfaces; prefer subdivision/auto there.
• sturm is exact only up to a moderate algebraic degree (default cap 16; raise max_degree with care — equally-spaced interpolation degrades at high degree).

Supported functions

These built-in functions are recognized inside expression strings and covered by the per-function test suite under tests/functions/:

Arity	Functions
1	abs, sign, sqrt, exp, log, sin, cos, tan, asin, acos, atan
2	+ (add), - (sub), * (mul), / (div), ^ / ** (pow)

Constants: pi, e. Imaginary unit: i (configurable via the imaginary_unit argument).

Numerical accuracy and drift

precision=N guarantees N decimal digits in the representation of a value, not that every operation lands on a mathematically exact result. Some expressions accumulate floating-point error at the last digits even though the math identity is exact.

Operations that are exact at the configured precision

• +, -, *, / on real numbers.
• +, -, *, / on complex numbers (including long * chains — i*i*i*i*i*i*i*i evaluates to exactly 1).
• ^ with real base and integer exponent (uses square-and-multiply internally): 2^10, 1.1^7, etc.
• abs(z) when re² + im² is a perfect square (abs(3+4*i) == 5, abs(i) == 1).
• Transcendentals at points with exact algebraic results: sin(0), cos(0), exp(0), log(1), etc.
• Pythagorean identity in real space: sin(pi/4)^2 + cos(pi/4)^2 == 1.

Operations that drift (visible at last digit or beyond)

• ^ with a complex base, even for trivial integer exponents. Complex pow(a, b) is computed as exp(b·log(a)), a transcendental chain — there is no integer-exponent specialization for complex values in boost.multiprecision. Drift grows with the magnitude of the exponent and the result:

  Expression	Mathematical value	Computed at precision=24
  i^4	1	1 + i·1e-24
  i^64	1	1 + i·1.2e-23
  (1+i)^16	256	256 + i·3.78e-21
  (1+i)^64	4294967296	4294967296 + i·2.5e-13 (≈ half of the available digits lost)
  exp(2*i*pi)	1	1 + i·1e-24

  Workaround: if the exponent is a small integer, replace z^n with the explicit product z*z*…*z — multiplication is exact for complex values.

• Transcendental functions (sin, cos, exp, log, asin, acos, atan, tan) when the mathematical result is irrational. log(exp(1)) rounds clean to 1 at the configured precision; the full memory chunk (format_digits=0) exposes the round-trip drift: 1.0000000000000000000000018…. The same expression with +i*0 takes the complex path and drifts the other way: 0.99999999999999999999999959….

Development

Building from source

A working C++ toolchain is required when no wheel matches your platform or when you check the project out for development:

OS	Toolchain
Linux	g++ (or clang++); python3-dev
macOS	Xcode Command Line Tools (xcode-select --install)
Windows	Visual Studio Build Tools (C++ workload)

git clone https://github.com/hozblok/formula.git
cd formula
python -m venv .venv &&source .venv/bin/activate
python -m pip install -e .[test]

The editable install builds the C++ extension and places the resulting _formula.so (or .pyd) directly into the source tree, so pytest from the repo root resolves formula._formula correctly.

Running tests

See doc/tests.md.

License

formula is provided under Apache license that can be found in the LICENSE file. By using, distributing, or contributing to this project, you agree to the terms and conditions of this license.