Cramer's rule: Difference between revisions

From Rosetta Code
Content added Content deleted
m (promoted to a regular task)
m (→‎{{header|Sidef}}: minor code simplification)
Line 218: Line 218:
var sign = +1
var sign = +1
var pivot = 1
var pivot = 1

a.range.each { |k|
a.range.each { |k|
var r = (k+1 .. a.end)
var r = (k+1 .. a.end)
var previous_pivot = pivot
var previous_pivot = pivot

if ((pivot = a[k][k]) == 0) {
a.swap(r.first_by {|i| a[i][k] != 0 } \\ (return 0), k)
if ((pivot = a[k][k])==0) {
a.swap(r.first_by { |i|
a[i][k] != 0
} \\ (return 0), k)
pivot = a[k][k]
pivot = a[k][k]
sign = -sign
sign = -sign
}
}

r ~X r -> each { |p|
var(i, j) = p...
for i,j in (r ~X r) {
((a[i][j] *= pivot) -= a[i][k]*a[k][j]) /= previous_pivot
((a[i][j] *= pivot)
-= a[i][k]*a[k][j]
) /= previous_pivot
}
}
}
}
Line 254: Line 260:


var free_terms = [-3, -32, -47, 49]
var free_terms = [-3, -32, -47, 49]
var (w, x, y, z) = cramers_rule(matrix, free_terms)...;
var (w, x, y, z) = cramers_rule(matrix, free_terms)...


say "w = #{w}"
say "w = #{w}"

Revision as of 13:57, 5 March 2016

Task
Cramer's rule
You are encouraged to solve this task according to the task description, using any language you may know.

In linear algebra, Cramer's rule is an explicit formula for the solution of a system of linear equations with as many equations as unknowns, valid whenever the system has a unique solution. It expresses the solution in terms of the determinants of the (square) coefficient matrix and of matrices obtained from it by replacing one column by the vector of right hand sides of the equations.


Given

which in matrix format is

Then the values of and can be found as follows:

Task

Given the following system of equations:

solve for , , and , using Cramer's rule.

EchoLisp

<lang scheme> (lib 'matrix) (string-delimiter "") (define (cramer A B (X)) ;; --> vector X (define ∆ (determinant A)) (for/vector [(i (matrix-col-num A))] (set! X (matrix-set-col! (array-copy A) i B)) (// (determinant X) ∆)))

(define (task) (define A (list->array

 	'( 2 -1 5 1 3 2 2 -6 1 3 3 -1 5 -2 -3 3) 4 4))

(define B #(-3 -32 -47 49)) (writeln "Solving A * X = B") (array-print A) (writeln "B = " B) (writeln "X = " (cramer A B))) </lang>

Output:
(task)
Solving A * X = B    
  2   -1   5    1  
  3   2    2    -6 
  1   3    3    -1 
  5   -2   -3   3  
B =      #( -3 -32 -47 49)    
X =      #( 2 -12 -4 1)    

J

Implementation:

<lang J>cramer=:4 :0

 A=. x [ b=. y
 det=. -/ .*
 A %~&det (i.#A) b"_`[`]}&.|:"0 2 A

)</lang>

Task data:

<lang J>A=: _&".;._2]t=: 0 :0

 2 -1  5  1
 3  2  2 -6
 1  3  3 -1
 5 -2 -3  3

)

b=: _3 _32 _47 49</lang>

Task example:

<lang J> A cramer b 2 _12 _4 1</lang>

Perl

<lang perl>use Math::Matrix;

sub cramers_rule {

   my ($A, $terms) = @_;
   my @solutions;
   my $det = $A->determinant;
   foreach my $i (0 .. $#{$A}) {
       my $Ai = $A->clone;
       foreach my $j (0 .. $#{$terms}) {
           $Ai->[$j][$i] = $terms->[$j];
       }
       push @solutions, $Ai->determinant / $det;
   }
   @solutions;

}

my $matrix = Math::Matrix->new(

   [2, -1,  5,  1],
   [3,  2,  2, -6],
   [1,  3,  3, -1],
   [5, -2, -3,  3],

);

my $free_terms = [-3, -32, -47, 49]; my ($w, $x, $y, $z) = cramers_rule($matrix, $free_terms);

print "w = $w\n"; print "x = $x\n"; print "y = $y\n"; print "z = $z\n";</lang>

Output:
w = 2
x = -12
y = -4
z = 1

Perl 6

<lang perl6>sub det(@matrix) {

   my @a = @matrix.map: { [|$_] };
   my $sign = +1;
   my $pivot = 1;
   for ^@a -> $k {
     my @r = ($k+1 .. @a.end);
     my $previous-pivot = $pivot;
     if 0 == ($pivot = @a[$k][$k]) {
       (my $s = @r.first: { @a[$_][$k] != 0 }) // return 0;
       (@a[$s],@a[$k]) = (@a[$k], @a[$s]);
       my $pivot = @a[$k][$k];
       $sign = -$sign;
     }
     for @r X @r -> ($i, $j) {
       ((@a[$i][$j] *= $pivot) -= @a[$i][$k]*@a[$k][$j]) /= $previous-pivot;
     }
   }
   $sign * $pivot

}

sub cramers_rule(@A, @terms) {

   gather for ^@A -> $i {
       my @Ai = @A.map: { [|$_] };
       for ^@terms -> $j {
           @Ai[$j][$i] = @terms[$j];
       }
       take det(@Ai);
   } »/» det(@A);

}

my @matrix = (

   [2, -1,  5,  1],
   [3,  2,  2, -6],
   [1,  3,  3, -1],
   [5, -2, -3,  3],

);

my @free_terms = (-3, -32, -47, 49); my ($w, $x, $y, $z) = |cramers_rule(@matrix, @free_terms);

say "w = $w"; say "x = $x"; say "y = $y"; say "z = $z";</lang>

Output:
w = 2
x = -12
y = -4
z = 1

Ruby

<lang ruby>require 'matrix'

def cramers_rule(a, terms)

   solutions = []
   det = a.determinant
   size = a.row_size;
   size.times do |i|
       cols = []
       (0..i-1).each { |j|
           cols << a.column(j)
       }
       cols << terms
       (i+1...size).each { |j|
           cols << a.column(j)
       }
       ai = Matrix.columns(cols)
       solutions << ai.determinant / det
   end
   solutions

end

matrix = Matrix[

   [2, -1,  5,  1],
   [3,  2,  2, -6],
   [1,  3,  3, -1],
   [5, -2, -3,  3],

]

vector = [-3, -32, -47, 49] puts cramers_rule(matrix, vector)</lang>

Output:
2
-12
-4
1

Sidef

<lang ruby>func det(a) {

   a = a.map{.map{_}}
   var sign = +1
   var pivot = 1
   a.range.each { |k|
     var r = (k+1 .. a.end)
     var previous_pivot = pivot
     if ((pivot = a[k][k])==0) {
       a.swap(r.first_by { |i|
           a[i][k] != 0
       } \\ (return 0), k)
       pivot = a[k][k]
       sign = -sign
     }
     for i,j in (r ~X r) {
       ((a[i][j] *= pivot)
           -= a[i][k]*a[k][j]
       ) /= previous_pivot
     }
   }
   sign * pivot

}

func cramers_rule(A, terms) {

   gather {
       A.each_index { |i|
           var Ai = A.map{.map{_}}
           terms.each_index { |j|
               Ai[j][i] = terms[j]
           }
           take(det(Ai))
       }
   } »/» det(A)

}

var matrix = [

   [2, -1,  5,  1],
   [3,  2,  2, -6],
   [1,  3,  3, -1],
   [5, -2, -3,  3],

]

var free_terms = [-3, -32, -47, 49] var (w, x, y, z) = cramers_rule(matrix, free_terms)...

say "w = #{w}" say "x = #{x}" say "y = #{y}" say "z = #{z}"</lang>

Output:
w = 2
x = -12
y = -4
z = 1