Gauss-Jordan matrix inversion: Difference between revisions

{{trans|Nim}}

 

 

F transformToRref(&mat)

L(j) 0 .< mat.len

I augmat[i][j] != Float(i == j)

result[i][j] = augmat[i][mat.len + j]

R result

runTest(m1)

runTest(m2)

 

{{out}}

The COPY file of FORMAT, to format a floating point number, can be found at:

[[360_Assembly_include|Include files 360 Assembly]].

GAUSSJOR CSECT

USING GAUSSJOR,R13 base register

PG DC CL80' ' buffer

REGEQU

{{out}}

<pre>

=={{header|ALGOL 60}}==

{{works with|A60}}

comment Gauss-Jordan matrix inversion - 22/01/2021;

integer n;

show("c",c)

end

{{out}}

<pre>

=={{header|ALGOL 68}}==

{{Trans|ALGOL 60}}

PROC rref = ( REF[,]REAL m )VOID:

BEGIN

c := inverse( b );

show( "c", c )

{{out}}

<pre>

1.000 0.000 -2.000 4.000

4.000 1.000 2.000 2.000

</pre>

 

=={{header|C}}==

{{trans|Fortran}}

gjinv - Invert a matrix, Gauss-Jordan algorithm

A is destroyed. Returns 1 for a singular matrix.

}

return 0;

Test program:

#include <stdio.h>

int main (void)

}

return 0;

Output is the same as the Fortran example.

 

=={{header|C sharp|C#}}==

using System;

 

}

}

using System;

 

}

}

{{out}}<pre>

-0.913043478260869 0.246376811594203 0.0942028985507246 0.413043478260869

 

=={{header|C++}}==

#include <cassert>

#include <iomanip>

print(std::cout, inverse(inv));

return 0;

 

{{out}}

0.0000 -1.0000 2.0000

</pre>

 

=={{header|Factor}}==

Normally you would call <code>recip</code> to calculate the inverse of a matrix, but it uses a different method than Gauss-Jordan, so here's Gauss-Jordan.

{{works with|Factor|0.99 2020-01-23}}

prettyprint sequences ;

 

{ 7 -8 9 1 }

{ 1 -2 1 3 }

{{out}}

<pre>

-13/23 16/69 -2/69 13/23

</pre>

=={{header|Fortran}}==

Note that the Crout algorithm is a more efficient way to invert a matrix.

! gjinv - Invert a matrix, Gauss-Jordan algorithm

! A is destroyed.

 

RETURN

Test program:

{{works with|gfortran|8.3.0}}

PROGRAM TINV

IMPLICIT NONE

END DO

 

{{out}}

<pre>

7.000 -8.000 9.000 1.000

1.000 -2.000 1.000 3.000

</pre>

 

=={{header|Generic}}==

// The Generic Language is a database compiler. This code is compiled into database then executed out of database.

 

 

}

=={{header|Go}}==

{{trans|Kotlin}}

 

import "fmt"

b := matrix{{2, -1, 0}, {-1, 2, -1}, {0, -1, 2}}

b.inverse().print("Inverse of B is:\n")

{{out}}

<pre>

===Alternative===

{{trans|PowerShell}}

import (

b := matrix{{2, -1, 0}, {-1, 2, -1}, {0, -1, 2}}

b.inverse().print("Inverse of B is:\n")

{{out}}

<pre>Same output than above</pre>

 

=={{header|Haskell}}==

 

isSquareMatrix xs = null xs || all ((== (length xs)).length) xs

mapM_ print $ inversion a

putStrLn "\ninversion b ="

{{out}}

<pre>

 

Uses Gauss-Jordan implementation (as described in [[Reduced_row_echelon_form#J]]) to find reduced row echelon form of the matrix after augmenting with an identity matrix.

augmentR_I1=: ,. e.@i.@# NB. augment matrix on the right with its Identity matrix

 

'''Usage:'''

1 2 3

4 1 6

_0.8125 0.125 0.1875

0.125 _0.25 0.125

 

=={{header|Java}}==

 

import java.util.Random;

Matrix.product(m, inv).print();

}

 

 

public class Matrix {

elements[j] = tmp;

}

 

{{out}}

'''Works with gojq, the Go implementation of jq'''

 

# it being understood that:

# matrix(0; _; _) evaluates to []

def matrix(m; n; init):

if m == 0 then []

elif m > 0 then

| [range(0;m) | $row ]

else error("matrix\(m);_;_) invalid")

if $nc <= .lead then ., break $out else . end

| .i = $r

.i += 1

| if $nr == .i

then .i = $r

| .lead += 1

else .

end

| reduce range(0; $nr) as $i ( matrix($nr; $nr; 0);

reduce range($nr; 2 *$nr) as $j (.;

'''The Task'''

[ [ 1, 2, 3],

[ 4, 1, 6],

;

 

{{out}}

<pre>

 

'''Built-in''' LAPACK-based linear solver uses partial-pivoted Gauss elimination):

@show I / A

 

'''Native implementation''':

size(A, 1) == size(A, 2) || throw(ArgumentError("A must be squared"))

n = size(A, 1)

 

A = rand(10, 10)

 

{{out}}

=={{header|Kotlin}}==

This follows the description of Gauss-Jordan elimination in Wikipedia whereby the original square matrix is first augmented to the right by its identity matrix, its reduced row echelon form is then found and finally the identity matrix to the left is removed to leave the inverse of the original square matrix.

 

typealias Matrix = Array<DoubleArray>

)

b.inverse().printf("Inverse of B is :\n")

 

{{out}}

 

=={{header|Lambdatalk}}==

{require lib_matrix}

 

[1.4444444444444444,-0.5555555555555556,0.1111111111111111],

[-0.05555555555555555,0.1111111111111111,-0.05555555555555555]]

 

=={{header|Mathematica}} / {{header|Wolfram Language}}==

elimmat = RowReduce[Join[a, IdentityMatrix[Length[a]], 2]];

{{out}}

<pre>{{0.0463243, -0.0563948, -0.0367573, 0.163646}, {0.0866062, 0.0684794, -0.133938, -0.020141}, {0.0694864, -0.0845921, 0.194864, -0.00453172}, {-0.149043, 0.137966, 0.00956697, -0.0699899}}</pre>

 

=={{header|MATLAB}}==

original = M;

[n,~] = size(M);

0, -1, 2 ];

 

 

{{out}}

=={{header|Nim}}==

We reuse the algorithm of task https://rosettacode.org/wiki/Reduced_row_echelon_form (version using floats).

 

const Eps = 1e-10

runTest(m1)

runTest(m2)

 

{{out}}

-0.69565 0.36232 0.07971 0.19565

-0.56522 0.23188 -0.02899 0.56522</pre>

 

=={{header|Perl}}==

Included code from [[Reduced_row_echelon_form#Perl|Reduced row echelon form]] task.

our @m; local *m = shift;

@m or return;

my @rt = gauss_jordan_invert( @gj );

print "After round-trip:\n" . display(\@rt) . "\n";} . "\n"

{{out}}

<pre>Original Matrix:

{{trans|Kotlin}}

uses ToReducedRowEchelonForm() from [[Reduced_row_echelon_form#Phix]]

<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>

<span style="color: #008080;">function</span> <span style="color: #000000;">ToReducedRowEchelonForm</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">M</span><span style="color: #0000FF;">)</span>

<span style="color: #0000FF;">{</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #0000FF;">-</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">}}</span>

<span style="color: #7060A8;">pp</span><span style="color: #0000FF;">(</span><span style="color: #000000;">inverse</span><span style="color: #0000FF;">(</span><span style="color: #000000;">test</span><span style="color: #0000FF;">),{</span><span style="color: #004600;">pp_Nest</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">})</span>

{{out}}

<pre>

===alternate===

{{trans|PowerShell}}

<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>

<span style="color: #008080;">function</span> <span style="color: #000000;">gauss_jordan_inv</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">a</span><span style="color: #0000FF;">)</span>

<span style="color: #7060A8;">puts</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"inv(a) = "</span><span style="color: #0000FF;">)</span>

<span style="color: #7060A8;">pp</span><span style="color: #0000FF;">(</span><span style="color: #000000;">inva</span><span style="color: #0000FF;">,{</span><span style="color: #004600;">pp_Nest</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #004600;">pp_Indent</span><span style="color: #0000FF;">,</span><span style="color: #000000;">9</span><span style="color: #0000FF;">,</span><span style="color: #004600;">pp_FltFmt</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"%4.2f"</span><span style="color: #0000FF;">})</span>

{{out}}

<pre>

 

=={{header|PowerShell}}==

function gauss-jordan-inv([double[][]]$a) {

$n = $a.count

show $invb

 

<b>Output:</b>

<pre>

Use numpy matrix inversion

 

import numpy as np

from numpy.linalg import inv

print(a)

print(ainv)

 

{{out}}

Using the matrix library that comes with default Racket installation

 

 

(require math/matrix

 

(inverse A)

 

{{out}}

=={{header|Raku}}==

(formerly Perl 6)

Uses bits and pieces from other tasks, [[Reduced_row_echelon_form#Raku|Reduced row echelon form]] primarily.

 

^@m .map: { @m[$_].append: identity(+@m)[$_] };

@m.&rref[*]»[+@m .. *];

sub identity ($n) { [ 1, |(0 xx $n-1) ], *.rotate(-1).Array ... *.tail }

 

sub rref (@m) {

 

for ^$rows -> $r {

}

@m[$i, $r] = @m[$r, $i] if $r != $i;

for ^$rows -> $n {

next if $n == $r;

}

++$lead;

@m

}

sub rat-or-int ($num) {

return $num unless $num ~~ Rat|FatRat;

$num.nude.join: '/';

}

say_it( ' Gauss-Jordan Inverted:', my @invert = gauss-jordan-invert @matrix );

say_it( ' Re-inverted:', gauss-jordan-invert @invert».Array );

 

{{out}}

=={{header|REXX}}==

===version 1===

Parse Arg seed nn

If seed='' Then

Return bs

 

{{out}} Using the defaults for seed and n

<pre>show 1 The given matrix

===version 2===

{{trans|PL/I}}

data= 8 3 7 5 9 12 10 11 6 2 4 13 14 15 16 17 /*the matrix element values. */

call build 4 /*assign data elements to the matrix. */

else _= _ right(format(x.r.c, w, f), w + f + length(.))

end /*c*/; say _

{{out|output|text=&nbsp; when using the internal default input:}}

<pre>

-0.1833 -0.7009 -0.1425 0.6163

-0.1064 -0.0855 0.0883 0.0779

</pre>

 

=={{header|Rust}}==

fn main() {

let mut a: Vec<Vec<f64>> = vec![vec![1.0, 2.0, 3.0],

}

}

{{out}}

<pre>

[0.24999999999999997, 0.49999999999999994, 0.7499999999999999]

</pre>

=={{header|Sidef}}==

Uses the '''rref(M)''' function from [https://rosettacode.org/wiki/Reduced_row_echelon_form#Sidef Reduced row echelon form].

{{trans|Raku}}

 

var I = M.len.of {|i|

say gauss_jordan_invert(A).map {

.map { "%6s" % .as_rat }.join(" ")

{{out}}

<pre>

-13/23 16/69 -2/69 13/23

</pre>

 

=={{header|Wren}}==

{{libheader|Wren-matrix}}

The above module does in fact use the Gauss-Jordan method for matrix inversion.

 

var arrays = [

Fmt.mprint(m.inverse, 9, 6)

System.print()

 

{{out}}

=={{header|zkl}}==

This uses GSL to invert a matrix via LU decomposition, not Gauss-Jordan.

m:=GSL.Matrix(3,3).set(1,2,3, 4,1,6, 7,8,9);

i:=m.invert();

i.format(10,4).println("\n");

{{out}}

<pre>

-0.0000, 0.0000, 1.0000

</pre>

{{out}}

<pre>

0.2500, 0.5000, 0.7500

</pre>