Geometric algebra

Geometric algebra is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

Geometric algebra is an other name for Clifford algebras and it's basically an algebra containing a vector space ${\mathcal {V}}$ and obeying the following axioms:

${\begin{array}{c}(ab)c=a(bc)\\a(b+c)=ab+ac\\(a+b)c=ac+bc\\\forall \mathbf {x} \in {\mathcal {V}},\,\mathbf {x} ^{2}\in \mathbb {R} \end{array}}$

The purpose of this task is to implement such an algebra for a vector space of countable, yet infinite dimension. You'll then pick three random elements of this algebra, along with a random vector, and verify that the axioms are respected.

Perl 6

<lang perl6>module Clifford;

1. Metric signature

our @signature = 1 xx *;

our class MultiVector {...}

sub e(UInt $n?) returns MultiVector is export {

   $n.defined ?? MultiVector.new(:blades(my Real %{UInt} = (1 +< $n) => 1)) !! MultiVector.new

}

my sub grade(UInt $n) is cached { [+] $n.base(2).comb } my sub order(UInt:D $i is copy, UInt:D $j) is cached {

   my $n = 0;
   repeat {

$i +>= 1; $n += [+] ($i +& $j).base(2).comb;

   } until $i == 0;
   return $n +& 1 ?? -1 !! 1;

}

sub metric-product(UInt $i, UInt $j) {

   my $r = order($i, $j);
   my $t = $i +& $j;
   my $k = 0;
   while $t !== 0 {

if $t +& 1 { $r *= @signature[$k]; } $t +>= 1; $k++;

   }
   return $r;

}

class MultiVector {

   has Real %.blades{UInt};
   method clean {

for %!blades {  %!blades{.key} :delete unless .value; }

   }
   multi method gist {

my sub blade-gist($blade) { join( '*', $blade.value, map { "e({$_ - 1})" }, grep +*, ($blade.key.base(2).comb.reverse Z* 1 .. *) ).subst(/<|w>1\*/, ) } if  %!blades == 0 { return '0' } elsif %!blades == 1 { given %!blades.pick { if .key == 0 { return .value.gist; } else { return blade-gist($_); } } } else { return join( ' + ', do for sort *.key, %!blades { .key == 0 ?? .value.gist !! blade-gist($_); } ).subst('+ -','- ', :g); }

   }
   method reals { %!blades.values }
   method max-grade { self.clean(); max map &grade, %!blades.keys }
   method AT-POS($n) {

MultiVector.new: :blades(my Real %{UInt} = %!blades.grep: *.key.&grade == $n)

   }	
   method narrow {

for %!blades { return self if .key > 0 && .value !== 0; } return %!blades{0} // 0;

   }
   method reverse {

[+] do for 0..self.max-grade -> $grade { (-1)**($grade*($grade - 1) div 2) * self[$grade]; }

   }

}

multi infix:<+>(MultiVector $A, MultiVector $B) returns MultiVector is export {

   my Real %blades{UInt} = $A.blades.clone;
   for $B.blades {

%blades{.key} += .value; %blades{.key} :delete unless %blades{.key};

   }
   return MultiVector.new: :%blades;

} multi infix:<+>(Real $s, MultiVector $A) returns MultiVector is export {

   my Real %blades{UInt} = $A.blades.clone;
   %blades{0} += $s;
   %blades{0} :delete unless %blades{0};
   return MultiVector.new: :%blades;

} multi infix:<+>(MultiVector $A, Real $s) returns MultiVector is export { $s + $A } multi infix:<*>(MultiVector $A, MultiVector $B) returns MultiVector is export {

   my Real %blades{UInt};
   for $A.blades -> $a {

for $B.blades -> $b { my $c = $a.key +^ $b.key; %blades{$c} += $a.value * $b.value * metric-product($a.key, $b.key); %blades{$c} :delete unless %blades{$c}; }

   }
   return MultiVector.new: :%blades;

} multi infix:<**>(MultiVector $ , 0) returns MultiVector is export { return MultiVector.new } multi infix:<**>(MultiVector $A, 1) returns MultiVector is export { return $A } multi infix:<**>(MultiVector $A, 2) returns MultiVector is export { return $A * $A } multi infix:<**>(MultiVector $A, UInt $n where $n %% 2) returns MultiVector is export {

   return ($A ** ($n div 2)) ** 2;

} multi infix:<**>(MultiVector $A, UInt $n) returns MultiVector is export {

   return $A * ($A ** ($n div 2)) ** 2;

}

multi infix:<*>(MultiVector $, 0) returns MultiVector is export { MultiVector.new } multi infix:<*>(MultiVector $A, 1) returns MultiVector is export { $A } multi infix:<*>(MultiVector $A, Real $s) returns MultiVector is export {

   return MultiVector.new: :blades(my Real %{UInt} = map { .key => $s * .value }, $A.blades);

} multi infix:<*>(Real $s, MultiVector $A) returns MultiVector is export { $A * $s } multi infix:</>(MultiVector $A, Real $s) returns MultiVector is export { $A * (1/$s) } multi prefix:<->(MultiVector $A) returns MultiVector is export { return -1 * $A } multi infix:<->(MultiVector $A, MultiVector $B) returns MultiVector is export { $A + -$B } multi infix:<->(MultiVector $A, Real $s) returns MultiVector is export { $A + -$s } multi infix:<->(Real $s, MultiVector $A) returns MultiVector is export { $s + -$A }

multi infix:<==>(MultiVector $A, MultiVector $B) returns Bool is export { $A - $B == 0 } multi infix:<==>(Real $x, MultiVector $A) returns Bool is export { $A == $x } multi infix:<==>(MultiVector $A, Real $x) returns Bool is export {

   my $narrowed = $A.narrow;
   $narrowed ~~ Real and $narrowed == $x;

}</lang>

The code required to test the axioms is the following:

<lang perl6>use Clifford; use Test;

plan 4;

sub random {

   [+] map {

Clifford::MultiVector.new: :blades(my Real %{UInt} = $_ => rand.round(.01))

   }, (^32).pick(5);

}

my ($a, $b, $c) = random() xx 3;

ok ($a*$b)*$c == $a*($b*$c), 'associativity'; ok $a*($b + $c) == $a*$b + $a*$c, 'left distributivity'; ok ($a + $b)*$c == $a*$c + $b*$c, 'right distributivity'; my @coeff = (.5 - rand) xx 4; my $v = [+] @coeff Z* map &e, ^4; ok ($v**2).narrow ~~ Real, 'contraction'; </lang>

Output:

1..4
ok 1 - associativity
ok 2 - left distributivity
ok 3 - right distributivity
ok 4 - contraction