/*
 * To change this template, choose Tools | Templates
 * and open the template in the editor.
 */
package apfloatextraclasses;

import org.apfloat.Apfloat;
import org.apfloat.ApfloatMath;

/**
 * Egy függvénygyűjtemény, ami segíti a mátrixokkal való munkát.
 * @author Esta
 */
public class ApfloatMatrixMath {

    /**
     * Készít egy size sorból és oszlopból álló egységmátrixot.
     * @param size a sorok száma.
     * @param precision a pontosság.
     * @return e az egységmátrix.
     */
    public static Apfloat[][] createEMatrix(int size, long precision) {
        //Létrehozzuk a mátrixot.
        Apfloat[][] e = new Apfloat[size][size];
        for (int i = 0; i < size; i++) {
            for (int j = 0; j < size; j++) {
                if (i == j) {
                    //a főátló csupa 1-es
                    e[i][j] = new Apfloat(1.0, precision);
                } else {
                    //a többi koordináta 0
                    e[i][j] = new Apfloat(0.0, precision);
                }
            }
        }
        return e;
    }

    /**
     * Megmondja, hogy a kapott m mátrix determinánsa 0-e
     * @param m a mátrix
     * @return true ha det(m)==0.
     */
    private static boolean determinantIsZero(Apfloat[][] m) {
        Apfloat[] tomb = new Apfloat[m.length];
        for (int i = 0; i < m.length; i++) {
            tomb[i] = m[i][i];
        }
        return (ApfloatMath.product(tomb).compareTo(new Apfloat(0.0)) == 0);
    }

    /**
     * Elvégzi az utolsó simításokat.
     * @param m a felső-háromszög mátrix
     * @param e az egységmátrixból keletkezett mátrix
     */
    private static Apfloat[][] finnishInvert(Apfloat[][] m, Apfloat[][] e) {
        //minden oszlopon végigmegyünk
        for (int i = 0; i < m.length; i++) {
            Apfloat l = m[i][i];
            //A főátlóban álló elemeket 1-re állítjuk.
            for (int k = 0; k < m.length; k++) {
                m[i][k] = m[i][k].divide(l);
                e[i][k] = e[i][k].divide(l);
            }
        }//minden oszlopban
        for (int i = m.length - 1; i > 0; i--) {
            //ha a főátló feletti elemek nem nulla értékűek, kinullázzuk őket
            for (int j = i - 1; j >= 0; j--) {
                if (m[j][i].compareTo(new Apfloat(0.0)) != 0) {
                    Apfloat l = m[j][i].divide(m[i][i]);
                    for (int k = 0; k < m.length; k++) {
                        m[j][k] = m[j][k].subtract(l.multiply(m[i][k]));
                        e[j][k] = e[j][k].subtract(l.multiply(e[i][k]));
                    }
                }
            }
        }
        return e;
    }

    /**
     * Megmondja a kapott kvadratikus mátrixr determinánsát.
     * @param m a mátrix.
     * @return det(m).
     * @throws apfloatextraclasses.ApfloatMatrixException ha a mátrix nem megfelelő.
     */
    public static Apfloat determinant(Apfloat[][] m) throws ApfloatMatrixException {
        if (m == null || m.length != m[0].length) {
            //ha az m nem kvadratikus, vagy üres, akkor nem folytathatjuk.
            throw new ApfloatMatrixException("ERROR_INVERT_NOT_NxN");
        }
        Apfloat[][] e = ApfloatMatrixMath.createEMatrix(m.length, m[0][0].precision());
        try {
            ApfloatMatrixMath.makeTriangular(m, e);
        } catch (ApfloatMatrixException ex) {
            if (ex.getMessage().equalsIgnoreCase("ERROR_INVERT_DETERMINANT_IS_ZERO")) {
                return new Apfloat(0.0, m[0][0].precision());
            } else {
                throw new ApfloatMatrixException(ex.getMessage());
            }
        }
        Apfloat[] tomb = new Apfloat[m.length];
        for (int i = 0; i < m.length; i++) {
            tomb[i] = m[i][i];
        }
        return (ApfloatMath.product(tomb));
    }

    /**
     * Elvégzi a Gauss-eliminációt a kapott kvadratikus mátrixra.
     * @param m a mátrix.
     * @return m^-1 az inverz mátrix
     * @throws apfloatextraclasses.ApfloatMatrixException ha a mátrix nem reguláris.
     */
    public static Apfloat[][] gaussianElimination(Apfloat[][] m) throws ApfloatMatrixException {
        if (m == null || m.length != m[0].length) {
            //ha az m nem kvadratikus, vagy üres, akkor m nem jó
            throw new ApfloatMatrixException("ERROR_DET_NOT_NxN");
        }
        Apfloat[][] e = ApfloatMatrixMath.createEMatrix(m.length, m[0][0].precision());
        try {
            ApfloatMatrixMath.makeTriangular(m, e);
        } catch (ApfloatMatrixException ex) {
            if (ex.getMessage().equalsIgnoreCase("ERROR_INVERT_DETERMINANT_IS_ZERO")) {
                throw new ApfloatMatrixException("ERROR_DET_NOT_REGULAR");
            } else {
                throw new ApfloatMatrixException(ex.getMessage());
            }
        }
        if (ApfloatMatrixMath.determinantIsZero(m)) {
            throw new ApfloatMatrixException("ERROR_DET_NOT_REGULAR");
        }
        e = ApfloatMatrixMath.finnishInvert(m, e);
        return e;
    }

    /**
     * Felső-háromszög mátrixot készít az m mátrixból.
     * @param m egy mátrix.
     * @param e egy egységmátrix.
     * @throws apfloatextraclasses.ApfloatMatrixException ha a mátrix nem megfelelő.
     */
    private static void makeTriangular(Apfloat[][] m, Apfloat[][] e) throws ApfloatMatrixException {
        //végigmegyünk minden oszlopon
        for (int i = 0; i < m.length; i++) {
            //és minden sorban nullázzuk a főátló alatti elemeket
            if (m[i][i].compareTo(new Apfloat(0.0)) == 0) {
                ApfloatMatrixMath.switchRow(m, e, i);
            }
            for (int j = i + 1; j < m.length; j++) {
                //Ha főátló alatti elem nem nulla, akkor nullázzuk.
                if (m[j][i].compareTo(new Apfloat(0.0)) != 0) {
                    Apfloat l = m[j][i].divide(m[i][i]);
                    //minden a sorban álló koordinátából kivonjuk az l-szeresét az i. sorénak
                    for (int k = 0; k < m.length; k++) {
                        m[j][k] = m[j][k].subtract(l.multiply(m[i][k]));
                        e[j][k] = e[j][k].subtract(l.multiply(e[i][k]));
                    }
                }
            }

        }
    }

    /**
     * A mátrixok két sorát felcseréli úgy, hogy az k. oszlopban a főátlóban ne legyen 0.
     * @param m az első mátrix
     * @param e a másik mátrix
     * @param k az oszlop
     * @throws apfloatextraclasses.ApfloatMatrixException ha csupa 0 sort talál.
     */
    private static void switchRow(Apfloat[][] m, Apfloat[][] e, int k) throws ApfloatMatrixException {
        int what = 0;
        for (int i = k; i < m[0].length; i++) {
            if (m[i][k].compareTo(new Apfloat(0.0)) != 0) {
                what = i;
                break;
            }
        }
        if (what == 0) {
            throw new ApfloatMatrixException("ERROR_INVERT_DETERMINANT_IS_ZERO");
        } else {
            for (int i = 0; i < m[0].length; i++) {
                Apfloat temp = m[what][i];
                m[what][i] = m[k][i];
                m[k][i] = temp.multiply(new Apfloat(-1.0,m[0][0].precision()));
                temp = e[what][i];
                e[what][i] = e[k][i];
                e[k][i] = temp.multiply(new Apfloat(-1.0,m[0][0].precision()));
            }
        }
    }

    /**
     * Transzponál egy mátrixot
     * @param m a mátrix
     * @return m^T a transzponált mátrix
     * @throws apfloatextraclasses.ApfloatMatrixException ha a mátrix nem megfelelő.
     */
    public static Apfloat[][] transponate(Apfloat[][] m) throws ApfloatMatrixException {
        //null nem transzponálható
        if (m == null) {
            throw new ApfloatMatrixException("ERROR_TRANSP_NULL");
        }
        //készítünk egy új mátrixot
        Apfloat[][] retValue = new Apfloat[m[0].length][m.length];
        //felcseréljük a sorokat és oszlopokat
        for (int i = 0; i < m.length; i++) {
            for (int j = 0; j < m[0].length; j++) {
                retValue[j][i] = m[i][j];
            }
        }
        return retValue;
    }

    /**
     * Összead két mátrixot.
     * @param a az első mátrix.
     * @param b a második mátrix.
     * @return a+b a két mátrix összege.
     * @throws apfloatextraclasses.ApfloatMatrixException ha a mátrixok mérete nem megfelelő.
     */
    public static Apfloat[][] add(Apfloat[][] a, Apfloat[][] b) throws ApfloatMatrixException {
        if (a == null || b == null || a.length != b.length || a[0].length != b[0].length) {
            throw new ApfloatMatrixException("ERROR_MADD");
        }
        Apfloat[][] c = new Apfloat[a.length][a[0].length];
        for (int i = 0; i < a.length; i++) {
            for (int j = 0; j < a[0].length; j++) {
                c[i][j] = a[i][j].add(b[i][j]);
            }
        }
        return c;
    }

    /**
     * Kivon egymásból két mátrixot.
     * @param a az első mátrix.
     * @param b a második mátrix.
     * @return a-b a két mátrix különbsége.
     * @throws apfloatextraclasses.ApfloatMatrixException ha a mátrixok mérete nem megfelelő.
     */
    public static Apfloat[][] subtract(Apfloat[][] a, Apfloat[][] b) throws ApfloatMatrixException {
        if (a == null || b == null || a.length != b.length || a[0].length != b[0].length) {
            throw new ApfloatMatrixException("ERROR_MSUB");
        }
        Apfloat[][] c = new Apfloat[a.length][a[0].length];
        for (int i = 0; i < a.length; i++) {
            for (int j = 0; j < a[0].length; j++) {
                c[i][j] = a[i][j].subtract(b[i][j]);
            }
        }
        return c;
    }

    /**
     * Összeszoroz két mátrixot.
     * @param a az első mátrix.
     * @param b a második mátrix.
     * @return a*b a két mátrix szorzata.
     * @throws apfloatextraclasses.ApfloatMatrixException ha a mátrixok mérete nem megfelelő.
     */
    public static Apfloat[][] multiply(Apfloat[][] a, Apfloat[][] b) throws ApfloatMatrixException {
        if (a == null || b == null || a[0].length != b.length) {
            throw new ApfloatMatrixException("ERROR_MMUL");
        }
        Apfloat[][] retValue = new Apfloat[a.length][b[0].length];
        for (int i = 0; i < a.length; i++) {
            //végigmegyünk minden soron az A-ban
            for (int j = 0; j < b[0].length; j++) {
                //minden oszlopon a B-ben
                Apfloat temp = new Apfloat(0.0, a[0][0].precision());
                for (int k = 0; k < a[0].length; k++) {
                    //sort oszloppal szorzunk
                    temp = temp.add(a[i][k].multiply(b[k][j]));
                }
                retValue[i][j] = temp;
            }
        }
        return retValue;
    }

    /**
     * Skalárisan összeszoroz két sorvektort.
     * @param a az első vektor.
     * @param b a második vektor.
     * @return a*b a két vektor skaláris szorzata.
     * @throws apfloatextraclasses.ApfloatMatrixException ha a vektorok mérete nem megfelelő.
     */
    public static Apfloat scalarRowMultiply(Apfloat[][] a, Apfloat[][] b) throws ApfloatMatrixException {
        if (a == null || b == null || a.length != b.length || a[0].length != b[0].length || a.length != 1) {
            throw new ApfloatMatrixException("ERROR_SRMUL");
        }
        Apfloat retValue = new Apfloat(0.0, a[0][0].precision());
        for (int i = 0; i < a[0].length; i++) {
            //végigmegyünk minden oszlopon az A-ban és B-ben és elvégezzük a szorzást
            retValue = retValue.add(a[0][i].multiply(b[0][i]));
        }
        return retValue;
    }

    /**
     * Skalárisan összeszoroz két oszlopvektort.
     * @param a az első vektor.
     * @param b a második vektor.
     * @return a*b a két vektor skaláris szorzata.
     * @throws apfloatextraclasses.ApfloatMatrixException ha a vektorok mérete nem megfelelő.
     */
    public static Apfloat scalarColoumnMultiply(Apfloat[][] a, Apfloat[][] b) throws ApfloatMatrixException {
        if (a == null || b == null || a.length != b.length || a[0].length != b[0].length || a[0].length != 1) {
            throw new ApfloatMatrixException("ERROR_SCMUL");
        }
        Apfloat retValue = new Apfloat(0.0, a[0][0].precision());
        for (int i = 0; i < a.length; i++) {
            //végigmegyünk minden soron az A-ban és B-ben és elvégezzük a szorzást
            retValue = retValue.add(a[i][0].multiply(b[i][0]));
        }
        return retValue;
    }
}
