symgeoalg

symbolic geometric algebra calculations
git clone git://git.rr3.xyz/symgeoalg
Log | Files | Refs | README | LICENSE

commit 524507dfa2a332f8f52718c7527f1ed74a15f504
parent a2dd74df3c1f89e01d97074736007a2263186bb8
Author: Robert Russell <robertrussell.72001@gmail.com>
Date:   Wed, 24 Jul 2024 23:06:00 -0700

Add geometric algebras with finite bases

Work in progress, but it seems to mostly work right now.

Diffstat:
MMain.lean | 125++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++-----
1 file changed, 117 insertions(+), 8 deletions(-)

diff --git a/Main.lean b/Main.lean @@ -4,6 +4,18 @@ class SMul (α β : Type) where smul : α -> β -> β infixr:75 " ⋅ " => SMul.smul +def BitVec.foldl (b : BitVec n) (f : α -> Bool -> α) (init : α) := + let rec loop (acc : α) (i : Nat) : α := + if i < n then + let acc' := f acc (BitVec.getLsb b i) + loop acc' (i + 1) + else + acc + loop init 0 + +def BitVec.popCount (b : BitVec n) : Nat := + b.foldl (fun acc bit => acc + Bool.toNat bit) 0 + class AddSemigroup (A : Type) where add : A -> A -> A class AddMonoid (A : Type) extends AddSemigroup A where @@ -12,7 +24,7 @@ class AddGroup (A : Type) extends AddMonoid A where sub : A -> A -> A neg : A -> A := fun x => sub zero x --- Notation for AddXxx +-- Notation for AddXXX instance [AddMonoid A] : OfNat A 0 where ofNat := AddMonoid.zero instance [AddSemigroup A] : Add A where @@ -30,7 +42,7 @@ class MulGroup (M : Type) extends MulMonoid M where div : M -> M -> M inv : M -> M := fun x => div one x --- Notation for MulXxx +-- Notation for MulXXX instance [MulMonoid M] : OfNat M 1 where ofNat := MulMonoid.one instance [MulSemigroup M] : Mul M where @@ -61,7 +73,7 @@ instance [Field F] : VectorSpace F F where smul := MulSemigroup.mul -- Given a type X and a vector space V over F, the type of functions X -> V is a vector space over F. -def FunctionSpace (X V : Type) := X -> V +def FunctionSpace (X V : Type) : Type := X -> V instance [AddGroup V] : AddGroup (FunctionSpace X V) where zero := fun _ => 0 add f g := fun x => f x + g x @@ -69,16 +81,113 @@ instance [AddGroup V] : AddGroup (FunctionSpace X V) where instance [Field F] [AddGroup V] [VectorSpace F V] : VectorSpace F (FunctionSpace X V) where smul a f := fun x => a ⋅ f x -abbrev NSpace (n : Nat) (V : Type) := FunctionSpace (Fin n) V +-- Function spaces with finite domain. +abbrev NSpace (n : Nat) (V : Type) : Type := FunctionSpace (Fin n) V instance [ToString V] : ToString (NSpace n V) where toString f := let coordToString i acc := toString (f i) :: acc s!"({String.intercalate ", " $ Fin.foldr n coordToString []})" - def vec (coords : List V) : NSpace coords.length V := coords.get +def QuadraticForm (F : Type) [Field F] (n : Nat) := Fin n -> F + +structure Blade (n : Nat) : Type where + bits : BitVec n + +def Blade.has (b : Blade n) (i : Fin n) : Bool := b.bits.getLsb i + +def Blade.cons (bit : Bool) (b : Blade i) : Blade (i + 1) := ⟨BitVec.cons bit b.bits⟩ + +def Blade.enum (n : Nat) : List (Blade n) := + let enum (i : Fin n) (b : List (Blade i.val)) : List (Blade (i.val + 1)) := + (Blade.cons false <$> b) ++ (Blade.cons true <$> b) + Fin.hIterate (List ∘ Blade) [⟨BitVec.nil⟩] enum + +instance : ToString (Blade n) where + toString b := + let vectorIndexStrings := Fin.foldl n (fun acc i => if b.has i then toString (n - i) :: acc else acc) [] + if b.bits != 0 then + "e" ++ String.join vectorIndexStrings + else + "" + +/- +def Blade.mulSign (x y : Blade n) : Int := + let rec countInversions (y : BitVec n) (count : Nat) (i : Nat) : Nat := + if i < n then + let y' := y >>> 1 + let count' := count + BitVec.popCount (x.bits &&& y') + countInversions y' count' (i + 1) + else + count + if countInversions y.bits 0 0 % 2 = 0 then 1 else -1 +-/ + +structure GeometricAlgebra (F : Type) [Field F] (n : Nat) (q : QuadraticForm F n) : Type where + fn : FunctionSpace (Blade n) F + +instance [Field F] {q : QuadraticForm F n} : AddGroup (GeometricAlgebra F n q) where + zero := ⟨0⟩ + add u v := ⟨u.fn + v.fn⟩ + sub u v := ⟨u.fn - v.fn⟩ + +instance [Field F] {q : QuadraticForm F n} : VectorSpace F (GeometricAlgebra F n q) where + smul a v := ⟨a ⋅ v.fn⟩ + +-- TODO: Rename variables. Reorder BitVec.foldl args. Separate notion of basis for which we +-- can print coordinates. + +-- Given u,v : GeometricAlgebra F n and b : Blade n, the b-th component of the product u * v +-- is a sum of 2^n terms, where each term is the product of the l-th component of u and the +-- r-th component of v, for all l,r : Blade n that multiply to b (modulo a scalar). The type +-- PartialTerm F i represents part of such a term, specifically the first i bits of l and r +-- and an associated scalar a. We also keep track of the parity of the l bits seen so far, so +-- that we know if we should flip the sign of a when we get the next r bit. This type should +-- be local to the multiplication function, but I don't think structures can be local in Lean +-- (ughh!). +structure PartialTerm (F : Type) [Field F] (i : Nat) where + a : F + lParity : F + l : BitVec i + r : BitVec i + +def PartialTerm.cons [Field F] (a : F) (l r : Bool) (pt : PartialTerm F i) : PartialTerm F (i + 1) := + let a' := a * pt.a * if r then pt.lParity else 1 + let lParity' := pt.lParity * if l then -1 else 1 + ⟨a', lParity', BitVec.cons l pt.l, BitVec.cons r pt.r⟩ + +def PartialTerm.eval [Field F] {q : QuadraticForm F n} (u v : GeometricAlgebra F n q) (pt : PartialTerm F n) : F := + pt.a * u.fn ⟨pt.l⟩ * v.fn ⟨pt.r⟩ + +instance [Field F] {q : QuadraticForm F n} : MulMonoid (GeometricAlgebra F n q) where + one := ⟨fun b => if b.bits = 0 then 1 else 0⟩ + mul u v := + let terms (b : Blade n) (i : Fin n) (acc : List (PartialTerm F i.val)) : List (PartialTerm F (i.val + 1)) := + if b.has i then + (PartialTerm.cons 1 false true <$> acc) ++ (PartialTerm.cons 1 true false <$> acc) + else + (PartialTerm.cons 1 false false <$> acc) ++ (PartialTerm.cons (q i) true true <$> acc) + let fn (b : Blade n) : F := + let ts := Fin.hIterate (List ∘ PartialTerm F) [⟨1, 1, BitVec.nil, BitVec.nil⟩] (terms b) + let x := PartialTerm.eval u v <$> ts + List.foldl AddSemigroup.add 0 x + ⟨fn⟩ + +-- TODO: Generalize this printing functionality to arbitrary FunctionSpaces with +-- enumerable domain? +instance [Field F] [ToString F] {q : QuadraticForm F n} : ToString (GeometricAlgebra F n q) where + toString v := + let termToString (b : Blade n) : String := + toString (v.fn b) ++ if b.bits != 0 then s!" * {toString b}" else "" + String.intercalate " + " (termToString <$> Blade.enum n) + def main : IO Unit := - let u : NSpace 2 Rat := vec [1, 0] - let v : NSpace 2 Rat := vec [0, 1] - let w := u + v + let F : Type := Rat + let n : Nat := 3 + let q : QuadraticForm F 3 := fun + | 0 => 1 + | 1 => 1 + | 2 => 1 + let u : GeometricAlgebra F n q := ⟨fun _ => 1⟩ + let w : GeometricAlgebra F n q := u + u IO.println w