commit a6c66052d01a7925848f6de02bd5b9e74a772179
parent ab38267ed2c9b057de5be6780c0e9184be28dd0b
Author: Robert Russell <robertrussell.72001@gmail.com>
Date: Sat, 27 Jul 2024 00:41:12 -0700
Organize and simplify geometric algebra mul
Currently there is a bug in mul for GeometricAlgebra.
Also, ToString for GeometricAlgebra orders terms weirdly.
Diffstat:
| M | Main.lean | | | 168 | ++++++++++++++++++++++++++++++------------------------------------------------- |
1 file changed, 63 insertions(+), 105 deletions(-)
diff --git a/Main.lean b/Main.lean
@@ -1,7 +1,7 @@
import Batteries.Data.Rat
--------------------------------------------------------------------------------
--- Util
+-- Fin util
def finSum (i : Fin (m + n)) : Sum (Fin m) (Fin n) :=
-- Commute m and n in i's type to please Fin.subNat.
@@ -33,6 +33,19 @@ def finProd (i : Fin (m * n)) : Prod (Fin m) (Fin n) :=
(q, r)
--------------------------------------------------------------------------------
+-- BitVec util
+
+def bitVecFoldl (f : α -> Bool -> α) (bv : BitVec n) (init : α) : α :=
+ Fin.foldl n (fun acc i => f acc (bv.getLsb i)) init
+
+def bitVecEnuml (n : Nat) (f : α -> BitVec n -> α) (init : α) : α :=
+ let loop acc i := f acc (BitVec.ofFin i)
+ Fin.foldl (2^n) loop init
+
+def bitVecPopCount (bv : BitVec n) : Nat :=
+ bitVecFoldl (fun acc bit => acc + bit.toNat) bv 0
+
+--------------------------------------------------------------------------------
-- Notation
instance : HPow Type Nat Type where
@@ -240,52 +253,9 @@ instance [Field F] : Algebra F F where
-- TODO: More Algebra instances
-
-
-
-def bitVecFoldl (f : α -> Bool -> α) (bv : BitVec n) (init : α) : α :=
- sorry
-
-def bitVecEnuml (n : Nat) (f : α -> BitVec n -> α) (init : α) : α :=
- let loop acc i := f acc (BitVec.ofFin i)
- Fin.foldl (2^n) loop init
-
-
-structure BasisBlade (n : Nat) : Type where
- bits : BitVec n
-
--- TODO: ?
-def BasisBlade.empty : BasisBlade n -> Bool := sorry
-
-def BasisBlade.has (bld : BasisBlade n) (i : Fin n) : Bool :=
- bld.bits.getLsb i
-
--- TODO: delete
--- def BasisBlade.cons (bit : Bool) (bld : BasisBlade i) : BasisBlade (i + 1) :=
--- ⟨BitVec.cons bit bld.bits⟩
--- def BasisBlade.enum (n : Nat) : List (BasisBlade n) :=
--- let enum (i : Fin n) (b : List (BasisBlade i.val)) : List (BasisBlade (i.val + 1)) :=
--- (BasisBlade.cons false <$> b) ++ (BasisBlade.cons true <$> b)
--- Fin.hIterate (List ∘ BasisBlade) [⟨BitVec.nil⟩] enum
-
-def BasisBlade.mulScalar [Field F] [AddGroup V] [VectorSpace F V] {b : FiniteBasis V n}
- (q : OrthoQuadraticForm b) (bld0 bld1 : BasisBlade n) : F :=
- bitVecFoldl x (bld0 &&& bld1) 0
-
-def BasisBlade.foldl (n : Nat) (f : α -> BasisBlade n -> α) (init : α) : α :=
- let loop acc i := f acc ⟨BitVec.ofFin i⟩
- Fin.foldl (2^n) loop init
-
-instance : ToString (BasisBlade n) where
- toString bld :=
- let indexToString i acc := if bld.has i then toString i :: acc else acc
- let indexStrings := Fin.foldr n indexToString []
- let sep := if n > 9 then "," else ""
- "{" ++ String.intercalate sep indexStrings ++ "}"
-
-
-
-
+--------------------------------------------------------------------------------
+-- Geometric algebra generation from vector spaces with a quadratic form and
+-- an orthogonal finite basis
structure GeometricAlgebra [Field F] [AddGroup V] [VectorSpace F V]
{b : FiniteBasis V n} (q : OrthoQuadraticForm b) : Type where
@@ -341,77 +311,65 @@ instance [Field F] [AddGroup V] [VectorSpace F V] {b : FiniteBasis V n}
-- Now, if we want the component of v * w corresponding to blade bld = 010, we
-- consider the following pairs of vbld and wbld:
-- vbld wbld
--- 010 000
--- 011 001
-- 000 010
-- 001 011
--- 110 100
--- 111 101
+-- 010 000
+-- 011 001
-- 100 110
-- 101 111
--- (These are just all 2^3 3-bit numbers, but with one column (vbld in this
+-- 110 100
+-- 111 101
+-- (These are just all 2^3 3-bit numbers, but with one column (wbld in this
-- case) XORed with bld.) Then we multiply the coefficients corresponding to
--- each row and sum the results. We call each row a "term". We build up each
--- row (i.e., (vbld, wbld) pair) incrementally, bit-by-bit, so an
-
--- When vbld and wbld share a bit, it "cancels"
--- and gets replaced by a scalar in accordance with the quadratic form. Also,
-
-
--- 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⟩
+-- each row and sum the results. Additionally, we must calculate a sign from
+-- vbld and wbld (in accordance with the equation ei * ej = - ej * ei for all
+-- basis vectors ei != ej), and we must calculuate a scalar by multiplying the
+-- "squares" of each basis vector shared by vbld and wbld, where "square" means
+-- apply the quadratic form.
instance [Field F] [AddGroup V] [VectorSpace F V] {b : FiniteBasis V n}
{q : OrthoQuadraticForm b} : MulMonoid (GeometricAlgebra q) where
- one := ⟨fun b => if b.bits = 0 then 1 else 0⟩
+ one := ⟨fun bld => if bld = 0 then 1 else 0⟩
mul v w :=
- let coords (bld : BasisBlade n) : F :=
- let addTerm (acc : F) (vbld : BasisBlade n) : F :=
- let wbld : BasisBlade n := ⟨vbld.bits ^^^ bld.bits⟩
- let bldScalar : F :=
- let vCoeff : F := v.coords vbld
- let wCoeff : F := w.coords wbld
- acc + v.coords vbld * w.coords wbld * vbld.mul wbld
- BasisBlade.foldl n addTerm 0
+ let sign (x y : BitVec n) : F :=
+ let addSwaps (acc : Nat) (i : Fin n) : Nat :=
+ acc + bitVecPopCount (x &&& (y >>> (i.val + 1)))
+ let nSwaps := Fin.foldl n addSwaps 0
+ if nSwaps % 2 = 0 then 1 else -1
+
+ let scalar (x y : BitVec n) : F :=
+ let z := x &&& y
+ let mulSquares (acc : F) (i : Fin n) : F :=
+ if z.getLsb i then acc * q.evalBasis i else acc
+ Fin.foldl n mulSquares 1
+
+ let coords (bld : BitVec n) : F :=
+ let addTerm (acc : F) (vbld : BitVec n) : F :=
+ let wbld : BitVec n := vbld ^^^ bld
+ acc + sign vbld wbld * scalar vbld wbld * v.coords vbld * w.coords wbld
+ bitVecEnuml n addTerm 0
⟨coords⟩
- --mul v w :=
- -- 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 coords (b : BasisBlade 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
- -- ⟨coords⟩
---------------------------------------------------------------------------------
+instance [Field F] [AddGroup V] [VectorSpace F V] {b : FiniteBasis V n}
+ {q : OrthoQuadraticForm b} : Algebra F (GeometricAlgebra q) where
--------------------------------------------------------------------------------
+-- Main
-def main : IO Unit :=
+def main : IO Unit := do
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
+ let V : Type := F ^ 3
+ let b : FiniteBasis V 3 := FiniteBasis.pow (FiniteBasis.field F) 3
+ let q : OrthoQuadraticForm b := ⟨
+ fun
+ | 0 => 1
+ | 1 => 1
+ | 2 => 1
+ ⟩
+ let G : Type := GeometricAlgebra q
+
+ let v : G := ⟨fun _ => 1⟩
+ let w : G := v * v
+
+ IO.println v
IO.println w
