Vector.java
package church.math;

import church.lang.Interval;
import church.lang.operators.Relational;
import church.lang.operators.Relational.$$equal;
import church.lang.operators.Relational.$$neq;
import church.primitives.Integers;

import static church.lang.Array.*;
import static church.lang.ArraySupport.get0;
import static church.lang.Collections.$Size;
import static church.lang.Error.error;
import static church.lang.operators.Arithmetic.*;
import static church.lang.operators.Streams.$$encode;

@SuppressWarnings("unchecked")
public class Vector<T> {
    private static final $$equal<Integer>       $S0 = Integers::$equal;
    private static final $$neq<Integer>         $S1 = Relational.$neq($S0);
    private static final $CommonLength<Integer> $S2 = Vector.commonLength($S1);

    public final T[] elements;

    public Vector(T... elements) {
        this.elements = elements;
    }

    public static <T> Vector<T> vector(T... elements) {
        return new Vector<>(elements);
    }

    public static interface $CreateVector<T> {
        public Vector<T> createVector(int length, java.util.function.IntFunction<T> generator);
    }

    public static <T> $CreateVector<T> createVector($CreateArray<T> $L0) {
        return (length, generator) -> vector($L0.createArray(length, generator));
    }

    public static interface $CreateMatrix<T> {
        public Vector<Vector<T>> createMatrix(int rowCount, int columnCount, church.lang.Functions.Function1<Integer, church.lang.Functions.Function1<Integer, T>> generator);
    }

    public static <T> $CreateMatrix<T> createMatrix($CreateVector<Vector<T>> $L0, $CreateVector<T> $L1) {
        return (rowCount, columnCount, generator) -> $L0.createVector(rowCount, i -> $L1.createVector(columnCount, j -> generator.of(i).of(j)));
    }

    public static <T> T[] elements(Vector<T> $0) {
        return $0.elements;
    }

    public static <T> int size(Vector<T> $0) {
        return $0.elements.length;
    }

    public static <T> T get(Vector<T> $0, int i) {
        return get0($0.elements, i);
    }

    public static interface $CommonLength<T> {
        public T commonLength(T l1, T l2);
    }

    public static <T> $CommonLength<T> commonLength($$neq<T> $L0) {
        return (l1, l2) -> $L0.$neq(l1, l2) ? error("Vector:: incompatible vectors") : l1;
    }

    public static interface $Sum1<T> {
        public T sum1(int n, church.lang.Functions.Function1<Integer, T> f);
    }

    public static <T> $Sum1<T> sum1($Additive_identity<T> $L0, $$sum<T> $L1) {
        return (n, f) -> new Interval(0, n).reduce($L0.additive_identity(), s -> i -> $L1.$sum(s, f.of(i)));
    }

    public static interface $MultiplicativeIdentityMatrix<T> {
        public Vector<Vector<T>> multiplicativeIdentityMatrix(int N);
    }

    public static <T> $MultiplicativeIdentityMatrix<T> multiplicativeIdentityMatrix($CreateMatrix<T> $L0, $Multiplicative_identity<T> $L1, $Additive_identity<T> $L2) {
        return N -> $L0.createMatrix(N, N, i -> j -> i == j ? $L1.multiplicative_identity() : $L2.additive_identity());
    }

    public static interface $Pivot<T> {
        public int pivot(Vector<Vector<T>> A, int row, int col);
    }

    public static <T> $Pivot<T> pivot($Size<Vector<Vector<T>>> $L0, $$neq<T> $L1, $Additive_identity<T> $L2) {
        return new $Pivot<T>() {
            public int pivot(Vector<Vector<T>> A, int row, int col) {
                return row == $L0.size(A) ? error("Division by singular matrix") : $L1.$neq(get(get(A, row), col), $L2.additive_identity()) ? row : pivot(A, (row + 1), col);
            }
        };
    }

    public static interface $GaussJordanStep<T> {
        public church.lang.Functions.Function1<Vector<Vector<T>>, Vector<Vector<T>>> gaussJordanStep(Vector<Vector<T>> A, int j);
    }

    public static <T> $GaussJordanStep<T> gaussJordanStep($Pivot<T> $L0, $CreateVector<Vector<T>> $L1, $Size<Vector<Vector<T>>> $L2, $CreateVector<T> $L3, $Size<Vector<T>> $L4, $$div<T, T> $L5, $$sub<T> $L6, $$prd<T, T> $L7) {
        return (A, j) -> {
            int pivot = $L0.pivot(A, j, j);
            return M -> $L1.createVector($L2.size(A), i -> {
                Vector<T> pivotRow = get(M, pivot);
                if (i == j) {
                    T d = get(get(A, pivot), i);
                    return $L3.createVector($L4.size(pivotRow), k -> $L5.$div(get(pivotRow, k), d));
                }
                int       source    = i == pivot ? j : i;
                Vector<T> sourceRow = get(M, source);
                T         factor    = $L5.$div(get(get(A, source), j), get(get(A, pivot), j));
                return $L3.createVector($L4.size(sourceRow), k -> $L6.$sub(get(sourceRow, k), $L7.$prd(factor, get(pivotRow, k))));
            });
        };
    }

    public static interface $GaussJordanElimination<T> {
        public Vector<Vector<T>> gaussJordanElimination(Vector<Vector<T>> N, Vector<Vector<T>> D, int j);
    }

    public static <T> $GaussJordanElimination<T> gaussJordanElimination($Size<Vector<Vector<T>>> $L0, $GaussJordanStep<T> $L1) {
        return new $GaussJordanElimination<T>() {
            public Vector<Vector<T>> gaussJordanElimination(Vector<Vector<T>> N, Vector<Vector<T>> D, int j) {
                if (j == $L0.size(D)) {
                    return N;
                }
                church.lang.Functions.Function1<Vector<Vector<T>>, Vector<Vector<T>>> op = $L1.gaussJordanStep(D, j);
                return gaussJordanElimination(op.of(N), op.of(D), (j + 1));
            }
        };
    }

    public static interface $BareissStep<T> {
        public church.lang.Functions.Function1<Vector<Vector<T>>, Vector<Vector<T>>> bareissStep(Vector<Vector<T>> A, int j);
    }

    public static <T> $BareissStep<T> bareissStep($Pivot<T> $L0, $Multiplicative_identity<T> $L1, $CreateVector<Vector<T>> $L2, $Size<Vector<Vector<T>>> $L3, $CreateVector<T> $L4, $Size<Vector<T>> $L5, $$dof<T> $L6, $$sub<T> $L7, $$prd<T, T> $L8) {
        return (A, j) -> {
            int pivot = $L0.pivot(A, j, j);
            T   d     = j == 0 ? $L1.multiplicative_identity() : get(get(A, (j - 1)), (j - 1));
            return M -> $L2.createVector($L3.size(A), i -> {
                Vector<T> pivotRow = get(M, pivot);
                if (i == j) {
                    return pivotRow;
                }
                int       source       = i == pivot ? j : i;
                Vector<T> sourceRow    = get(M, source);
                T         sourceFactor = get(get(A, source), j);
                T         pivotFactor  = get(get(A, pivot), j);
                return $L4.createVector($L5.size(sourceRow), k -> $L6.$dof($L7.$sub($L8.$prd(pivotFactor, get(sourceRow, k)), $L8.$prd(sourceFactor, get(pivotRow, k))), d));
            });
        };
    }

    public static interface $BareissElimination<T, U> {
        public U bareissElimination(church.lang.Functions.Function1<Vector<Vector<T>>, church.lang.Functions.Function1<Vector<Vector<T>>, U>> rtn, Vector<Vector<T>> N, Vector<Vector<T>> D, int j);
    }

    public static <T, U> $BareissElimination<T, U> bareissElimination($Size<Vector<Vector<T>>> $L0, $BareissStep<T> $L1) {
        return new $BareissElimination<T, U>() {
            public U bareissElimination(church.lang.Functions.Function1<Vector<Vector<T>>, church.lang.Functions.Function1<Vector<Vector<T>>, U>> rtn, Vector<Vector<T>> N, Vector<Vector<T>> D, int j) {
                if (j == $L0.size(D)) {
                    return rtn.of(N).of(D);
                }
                church.lang.Functions.Function1<Vector<Vector<T>>, Vector<Vector<T>>> op = $L1.bareissStep(D, j);
                return bareissElimination(rtn, op.of(N), op.of(D), (j + 1));
            }
        };
    }

    public static interface $Cofactors<T> {
        public Vector<Vector<T>> cofactors(Vector<Vector<T>> m);
    }

    public static <T> $Cofactors<T> cofactors($BareissElimination<T, Vector<Vector<T>>> $L0, $MultiplicativeIdentityMatrix<T> $L1, $Size<Vector<Vector<T>>> $L2) {
        return m -> $L0.bareissElimination(N -> D -> N, $L1.multiplicativeIdentityMatrix($L2.size(m)), m, 0);
    }

    public static interface $Determinant<T> {
        public T determinant(Vector<Vector<T>> m);
    }

    public static <T> $Determinant<T> determinant($BareissElimination<T, T> $L0, $Size<Vector<Vector<T>>> $L1, $MultiplicativeIdentityMatrix<T> $L2) {
        return m -> $L0.bareissElimination(N -> D -> get(get(D, ($L1.size(m) - 1)), ($L1.size(m) - 1)), $L2.multiplicativeIdentityMatrix($L1.size(m)), m, 0);
    }

    public static <T, U> $$encode<T, Vector<U>> $encode($WriteEmbracedElementsToStream<T, U, String, String> $L0) {
        return (stream, $0) -> $L0.writeEmbracedElementsToStream(stream, $0.elements, "<", ">");
    }

    public static <T> $$equal<Vector<T>> $equal($$equal<T[]> $L0) {
        return ($0, $1) -> $L0.$equal($0.elements, $1.elements);
    }

    public static <T> $$sum<Vector<T>> $sum($CreateVector<T> $L0, $Size<Vector<T>> $L1, $$sum<T> $L2) {
        return (a, b) -> $L0.createVector($S2.commonLength($L1.size(a), $L1.size(b)), i -> $L2.$sum(get(a, i), get(b, i)));
    }

    public static <T> $$neg<Vector<T>> $neg($CreateVector<T> $L0, $Size<Vector<T>> $L1, $$neg<T> $L2) {
        return a -> $L0.createVector($L1.size(a), i -> $L2.$neg(get(a, i)));
    }

    public static interface $DotProduct<T, U> {
        public U dotProduct(Vector<T> a, Vector<U> b);
    }

    public static <T, U> $DotProduct<T, U> dotProduct($Sum1<U> $L0, $Size<Vector<T>> $L1, $Size<Vector<U>> $L2, $$prd<T, U> $L3) {
        return (a, b) -> $L0.sum1($S2.commonLength($L1.size(a), $L2.size(b)), i -> $L3.$prd(get(a, i), get(b, i)));
    }

    public static <T, U> $$prd<Vector<Vector<T>>, Vector<Vector<U>>> $prd($CreateVector<Vector<U>> $L0, $Size<Vector<Vector<T>>> $L1, $CreateVector<U> $L2, $Size<Vector<U>> $L3, $Sum1<U> $L4, $Size<Vector<T>> $L5, $$prd<T, U> $L6) {
        return (m, n) -> $L0.createVector($L1.size(m), i -> $L2.createVector($L3.size(get(n, i)), j -> $L4.sum1($L5.size(get(m, i)), k -> $L6.$prd(get(get(m, i), k), get(get(n, k), j)))));
    }

    public static <T> $$div<Vector<Vector<T>>, Vector<Vector<T>>> $div($GaussJordanElimination<T> $L0) {
        return (n, d) -> $L0.gaussJordanElimination(n, d, 0);
    }

    public static <T> $$rcp<Vector<Vector<T>>> $rcp($GaussJordanElimination<T> $L0, $MultiplicativeIdentityMatrix<T> $L1, $Size<Vector<Vector<T>>> $L2) {
        return d -> $L0.gaussJordanElimination($L1.multiplicativeIdentityMatrix($L2.size(d)), d, 0);
    }

}