Implementando el método de Gauss Seidel en Java

Posted on by Rudeus Greyrat

Gauss-Seidel (también conocido como método de desplazamiento sucesivo) es un método de cálculo matemático utilizado principalmente para encontrar la solución de un Sistema de Álgebra Lineal. ¿Alguna vez has oído hablar del Método Jacobi? Bueno, el Gauss-Seidel no es más que una mejor versión del Jacobi.  

Ejemplo :

Input:
Enter the number of variables in the equation:
2
Enter the augmented matrix:
1 2 3
6 5 4
 
Output: 
6.0 5.0 4.0 
1.0 2.0 3.0 
X0 = {0.6666666666666666 1.1666666666666667 }
X1 = {-0.30555555555555564 1.652777777777778 }
X2 = {-0.7106481481481481 1.855324074074074 }
X3 = {-0.8794367283950617 1.9397183641975309 }
X4 = {-0.9497653034979425 1.9748826517489713 }
X5 = {-0.9790688764574759 1.9895344382287379 }
X6 = {-0.9912786985239483 1.9956393492619742 }
X7 = {-0.9963661243849785 1.9981830621924892 }
X8 = {-0.9984858851604077 1.9992429425802039 }
X9 = {-0.9993691188168363 1.999684559408418 }
X10 = {-0.9997371328403484 1.999868566420174 }
X11 = {-0.9998904720168117 1.9999452360084058 }
X12 = {-0.999954363340338 1.999977181670169 }
X13 = {-0.9999809847251406 1.9999904923625702 }
X14 = {-0.9999920769688085 1.9999960384844042 }
X15 = {-0.9999966987370034 1.9999983493685016 }
X16 = {-0.9999986244737512 1.9999993122368755 }
X17 = {-0.9999994268640631 1.9999997134320315 }
X18 = {-0.9999997611933598 1.9999998805966799 }
X19 = {-0.9999999004972331 1.9999999502486165 }
X20 = {-0.9999999585405137 1.9999999792702567 }
X21 = {-0.999999982725214 1.999999991362607 }
X22 = {-0.9999999928021724 1.9999999964010862 }
X23 = {-0.999999997000905 1.9999999985004524 }
X24 = {-0.999999998750377 1.9999999993751885 }
X25 = {-0.9999999994793237 1.9999999997396618 }
X26 = {-0.9999999997830514 1.9999999998915257 }
X27 = {-0.9999999999096048 1.9999999999548024 }
X28 = {-0.9999999999623352 1.9999999999811675 }
X29 = {-0.9999999999843061 1.999999999992153 }
X30 = {-0.9999999999934606 1.9999999999967302 }
X31 = {-0.9999999999972751 1.9999999999986375 }
X32 = {-0.9999999999988646 1.9999999999994322 }
X33 = {-0.9999999999995268 1.9999999999997633 }
X34 = {-0.9999999999998028 1.9999999999999014 }
X35 = {-0.9999999999999176 1.9999999999999587 }
X36 = {-0.9999999999999656 1.9999999999999827 }
X37 = {-0.9999999999999855 1.9999999999999927 }
X38 = {-0.9999999999999938 1.999999999999997 }
X39 = {-0.9999999999999973 1.9999999999999987 }
X40 = {-0.9999999999999988 1.9999999999999993 }
X41 = {-0.9999999999999993 1.9999999999999996 }

Un algoritmo matemático básico para Gauss Seidel:

Dado un conjunto de n ecuaciones y n incógnitas:

a11x1+a12x2+a13x3+......+a1nxn=c1

a21x1+a22x2+a23x3+......+a2nxn=c2
.                        .
.                        .
.                        .
a1nx1+a2nx2+a3nx3+......+annxn=cn

Si los elementos de la diagonal son distintos de cero, cada ecuación se reescribe para la incógnita correspondiente, es decir, la primera ecuación se reescribe con x 1 en el lado izquierdo, la segunda ecuación se reescribe con x 2 en el lado izquierdo y así de la siguiente manera:

x1=(c1-a12x2-a13x3-....-a1nxn)/a11
x2=(c2-a21x1-a23x3-....-a2nxn)/a22
.
.
.
xn=(cn-an1x1-an2x2-....-an.n-1xn-1)/ann

Estas ecuaciones se pueden reescribir en la forma como:

x1=(c1-nj=1,j!=1∑a1jxj)/a11
x2=(c2-nj=1,j!=2∑a2jxj)/a22
.
.
.
x2=(cn-nj=1,j!=n∑anjxj)/ann

Por lo tanto, para cualquier rowi;

xi=(ci- nj=i,j!=i∑aijxj)/aii , i=1,2,..,n

Nota: Gauss-Seidel es aplicable a definida positiva estrictamente diagonal dominante o simétrica.

En resumen, para un sistema lineal de ecuaciones, digamos: 

Ax=b

Fórmula deducida:

xi(k)=(bi-∑j<i a i j x j(k)-∑j>i a i j x j(k-1))/aii

En términos de arrays, la definición de Gauss-Seidel:

x(k)=(D-L)-1(Ux(k-1)+b)

Aquí, D representa la parte diagonal de la array A, -L parte triangular estrictamente inferior de A, U parte triangular estrictamente superior de A.

Algoritmo de programación:

  1. Tome el número de variables en la ecuación y los valores de la array aumentada 
     (formada agregando las columnas de las dos arrays dadas) como entrada.
  2. Comprobando si la array aumentada dada es diagonalmente dominante.
  3. Intentar transformar la array aumentada dada en diagonalmente dominante. Si se demuestra lo contrario, se devuelve un mensaje apropiado al usuario.
  4. Cálculo e impresión del resultado final así como de las iteraciones.
     

Implementación en Java: 

Java

// Implementing Gauss Seidel Method in Java
import java.io.BufferedReader;
import java.io.IOException;
import java.io.InputStreamReader;
import java.io.PrintWriter;
import java.util.Arrays;
import java.util.StringTokenizer;
 
class GFG {
    // we set a max number of iterations to
    // prevent an infinite loop
    public static final int MAX_ITERATIONS = 100;
    private double[][] M;
    public GFG(double[][] matrix) { M = matrix; }
 
    public void print() // printing
    {
        int n = M.length;
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < n + 1; j++)
                System.out.print(M[i][j] + " ");
            System.out.println();
        }
    }
    // attempting to change a matrix to dominant
    // if proved that it is not
    public boolean transformToDominant(int r, boolean[] V,
                                       int[] R)
    {
        int n = M.length;
        if (r == M.length) {
            double[][] T = new double[n][n + 1];
            for (int i = 0; i < R.length; i++) {
                for (int j = 0; j < n + 1; j++)
                    T[i][j] = M[R[i]][j];
            }
            M = T;
            return true;
        }
        for (int i = 0; i < n; i++) {
            if (V[i])
                continue;
            double sum = 0;
            for (int j = 0; j < n; j++)
                sum += Math.abs(M[i][j]);
            if (2 * Math.abs(M[i][r]) > sum) {
                // diagonally dominant?
                V[i] = true;
                R[r] = i;
                if (transformToDominant(r + 1, V, R))
                    return true;
                V[i] = false;
            }
        }
        return false;
    }
    // method to check whether matrix is
    // diagonally dominant or not
    public boolean makeDominant()
    {
        boolean[] visited = new boolean[M.length];
        int[] rows = new int[M.length];
        Arrays.fill(visited, false);
        return transformToDominant(0, visited, rows);
    }
 
    // method to find the solution of the matrix
    // after all conditions are satisfied
    public void solve()
    {
        int iterations = 0;
        int n = M.length;
        double epsilon = 1e-15;
        double[] X = new double[n]; // Approximations
        double[] P = new double[n]; // Prev
        Arrays.fill(X, 0);
        while (true) {
            for (int i = 0; i < n; i++) {
                double sum = M[i][n]; // b_n
                for (int j = 0; j < n; j++)
                    if (j != i)
                        sum -= M[i][j] * X[j];
                // Update xi to use in the next
                // row calculation
                X[i] = 1 / M[i][i] * sum;
            }
            System.out.print("X" + iterations + " = {");
            for (int i = 0; i < n; i++)
                System.out.print(X[i] + " ");
            System.out.println("}");
            iterations++;
            if (iterations == 1)
                continue;
            boolean stop = true;
            for (int i = 0; i < n && stop; i++)
                if (Math.abs(X[i] - P[i]) > epsilon)
                    stop = false;
            if (stop || iterations == MAX_ITERATIONS)
                break;
            P = (double[])X.clone();
        }
    }
 
    public static void main(String[] args)
        throws IOException
    {
        PrintWriter writer
            = new PrintWriter(System.out, true);
        int n = 2, k = 1;
        double[][] M = new double[n][n + 1];
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < n + 1; j++)
                M[i][j] = k++;
        }
        GFG gausSeidel = new GFG(M);
        if (!gausSeidel.makeDominant()) {
 
            // if it is found that a matrix cannot be
            // changed into diagonally dominant then we
            // return the message to the user
 
            writer.println(
                "The system isn't diagonally dominant: "
                + "The method cannot guarantee convergence.");
        }
        writer.println();
        gausSeidel.print();
        gausSeidel.solve();
    }
}
Producción

The system isn't diagonally dominant: The method cannot guarantee convergence.

1.0 2.0 3.0 
4.0 5.0 6.0 
X0 = {3.0 -1.2000000000000002 }
X1 = {5.4 -3.1200000000000006 }
X2 = {9.240000000000002 -6.192000000000002 }
X3 = {15.384000000000004 -11.107200000000004 }
X4 = {25.21440000000001 -18.97152000000001 }
X5 = {40.94304000000002 -31.554432000000016 }
X6 = {66.10886400000004 -51.68709120000003 }
X7 = {106.37418240000007 -83.89934592000006 }
X8 = {170.79869184000012 -135.4389534720001 }
X9 = {273.8779069440002 -217.90232555520015 }
X10 = {438.8046511104003 -349.84372088832026 }
X11 = {702.6874417766405 -560.9499534213124 }
X12 = {1124.899906842625 -898.7199254740999 }
X13 = {1800.4398509481998 -1439.1518807585599 }
X14 = {2881.3037615171197 -2303.843009213696 }
X15 = {4610.686018427392 -3687.348814741914 }
X16 = {7377.697629483828 -5900.958103587062 }
X17 = {11804.916207174125 -9442.7329657393 }
X18 = {18888.4659314786 -15109.57274518288 }
X19 = {30222.14549036576 -24176.51639229261 }
X20 = {48356.03278458522 -38683.62622766818 }
X21 = {77370.25245533636 -61895.00196426909 }
X22 = {123793.00392853818 -99033.20314283055 }
X23 = {198069.4062856611 -158454.3250285289 }
X24 = {316911.6500570578 -253528.12004564624 }
X25 = {507059.2400912925 -405646.192073034 }
X26 = {811295.384146068 -649035.1073168544 }
X27 = {1298073.2146337088 -1038457.3717069671 }
X28 = {2076917.7434139343 -1661532.9947311475 }
X29 = {3323068.989462295 -2658453.991569836 }
X30 = {5316910.983139672 -4253527.586511738 }
X31 = {8507058.173023475 -6805645.338418781 }
X32 = {1.3611293676837562E7 -1.088903374147005E7 }
X33 = {2.17780704829401E7 -1.742245518635208E7 }
X34 = {3.484491337270416E7 -2.787592949816333E7 }
X35 = {5.575186199632666E7 -4.460148839706133E7 }
X36 = {8.920297979412267E7 -7.136238263529813E7 }
X37 = {1.4272476827059627E8 -1.1417981341647702E8 }
X38 = {2.2835962983295405E8 -1.8268770266636324E8 }
X39 = {3.653754083327265E8 -2.923003254661812E8 }
X40 = {5.846006539323624E8 -4.6768052194588995E8 }
X41 = {9.353610468917799E8 -7.48288836313424E8 }
X42 = {1.496577675626848E9 -1.1972621393014784E9 }
X43 = {2.394524281602957E9 -1.9156194240823655E9 }
X44 = {3.831238851164731E9 -3.064991079731785E9 }
X45 = {6.12998216246357E9 -4.903985728770856E9 }
X46 = {9.807971460541712E9 -7.84637716723337E9 }
X47 = {1.569275433746674E10 -1.2554203468773392E10 }
X48 = {2.5108406940546783E10 -2.0086725551237427E10 }
X49 = {4.017345110547485E10 -3.2138760883179886E10 }
X50 = {6.427752176935977E10 -5.142201741428782E10 }
X51 = {1.0284403483157564E11 -8.227522786406052E10 }
X52 = {1.6455045573112103E11 -1.3164036458369684E11 }
X53 = {2.6328072917039368E11 -2.1062458333511496E11 }
X54 = {4.212491666732299E11 -3.36999333337384E11 }
X55 = {6.73998666677768E11 -5.391989333410144E11 }
X56 = {1.0783978666850288E12 -8.627182933468231E11 }
X57 = {1.7254365866966462E12 -1.3803492693561172E12 }
X58 = {2.7606985387152344E12 -2.208558830970988E12 }
X59 = {4.417117661944976E12 -3.533694129554781E12 }
X60 = {7.067388259112562E12 -5.65391060728885E12 }
X61 = {1.13078212145807E13 -9.04625697166336E12 }
X62 = {1.809251394332972E13 -1.4474011154662576E13 }
X63 = {2.8948022309328152E13 -2.3158417847461324E13 }
X64 = {4.631683569492565E13 -3.705346855593932E13 }
X65 = {7.410693711188164E13 -5.928554968950412E13 }
X66 = {1.1857109937901123E14 -9.48568795032078E13 }
X67 = {1.897137590064186E14 -1.517710072051337E14 }
X68 = {3.035420144102704E14 -2.4283361152821512E14 }
X69 = {4.8566722305643325E14 -3.8853377844514544E14 }
X70 = {7.770675568902939E14 -6.216540455122339E14 }
X71 = {1.2433080910244708E15 -9.946464728195755E14 }
X72 = {1.989292945639154E15 -1.591434356511322E15 }
X73 = {3.182868713022647E15 -2.5462949704181165E15 }
X74 = {5.092589940836236E15 -4.0740719526689875E15 }
X75 = {8.148143905337978E15 -6.518515124270381E15 }
X76 = {1.3037030248540764E16 -1.042962419883261E16 }
X77 = {2.0859248397665224E16 -1.668739871813218E16 }
X78 = {3.3374797436264364E16 -2.6699837949011492E16 }
X79 = {5.3399675898022984E16 -4.2719740718418392E16 }
X80 = {8.5439481436836784E16 -6.8351585149469432E16 }
X81 = {1.36703170298938864E17 -1.09362536239151104E17 }
X82 = {2.18725072478302208E17 -1.74980057982641792E17 }
X83 = {3.4996011596528358E17 -2.7996809277222688E17 }
X84 = {5.5993618554445376E17 -4.4794894843556301E17 }
X85 = {8.9589789687112602E17 -7.1671831749690086E17 }
X86 = {1.43343663499380173E18 -1.14674930799504141E18 }
X87 = {2.29349861599008282E18 -1.8347988927920663E18 }
X88 = {3.6695977855841326E18 -2.9356782284673065E18 }
X89 = {5.871356456934613E18 -4.697085165547691E18 }
X90 = {9.394170331095382E18 -7.5153362648763064E18 }
X91 = {1.5030672529752613E19 -1.2024538023802092E19 }
X92 = {2.4049076047604183E19 -1.9239260838083346E19 }
X93 = {3.847852167616669E19 -3.0782817340933358E19 }
X94 = {6.1565634681866715E19 -4.925250774549338E19 }
X95 = {9.850501549098675E19 -7.880401239278941E19 }
X96 = {1.5760802478557882E20 -1.2608641982846306E20 }
X97 = {2.5217283965692612E20 -2.017382717255409E20 }
X98 = {4.034765434510818E20 -3.227812347608655E20 }
X99 = {6.45562469521731E20 -5.164499756173848E20 }

Complejidad temporal: O(N 2 )

Publicación traducida automáticamente

Artículo escrito por dikshapatro y traducido por Barcelona Geeks. The original can be accessed here. Licence: CCBY-SA

Deja una respuesta

Tu dirección de correo electrónico no será publicada. Los campos obligatorios están marcados con *