## 4.3-1
##
source('utils.R')

if( TRUE ){
    ## Demonstrate Gaussing Elimination with scaled partial pivoting on the matrix from Exercise 4.2-8:
    ##

    W = matrix(data=c(0.21, 0.32, 0.12, 0.31, 0.96, 0.1, 0.15, 0.24, 0.22, 0.71, 0.2, 0.24, 0.46, 0.36, 1.26, 0.61, 0.4, 0.32, 0.2, 1.53), nrow=4, ncol=5, byrow=TRUE)
    n = dim(W)[1]
    A = W[, 1:n]
    b = W[, n+1]
    x_true = solve(A, b)
    results = GE_w_scaled_partial_pivoting(A, b); x_wwx = results$X
    err = sum((x_true - x_wwx)^2)
    print(err)

    ## Generate some random test cases:
    ##
    n = 4
    n = 5
    n = 10
    for(ri in 1:10){
        A = matrix(data=rnorm(n^2), nrow=n, ncol=n, byrow=T)
        b = matrix(data=rnorm(n), nrow=n, ncol=1, byrow=T)
        x_true = solve(A, b)
        results = GE_w_scaled_partial_pivoting(A, b); x_wwx = results$X
        err = sum((x_true - x_wwx)^2)
        print(sprintf('ri= %5d; n= %5d; err= %10.8f', ri, n, err))
    }    
}

if( TRUE ){
    ## Demonstrate Gaussing Elimination with total pivoting on the matrix from Exercise 4.2-8:
    ##

    W = matrix(data=c(0.21, 0.32, 0.12, 0.31, 0.96, 0.1, 0.15, 0.24, 0.22, 0.71, 0.2, 0.24, 0.46, 0.36, 1.26, 0.61, 0.4, 0.32, 0.2, 1.53), nrow=4, ncol=5, byrow=TRUE)
    n = dim(W)[1]
    A = W[, 1:n]
    b = W[, n+1]
    x_true = solve(A, b)
    results = GE_w_total_pivoting(A, b); x_wwx = results$X
    err = sum((x_true - x_wwx)^2)
    print(err)

    ## Generate some random test cases:
    ##
    n = 4
    n = 5
    n = 10
    for(ri in 1:10){
        A = matrix(data=rnorm(n^2), nrow=n, ncol=n, byrow=T)
        b = matrix(data=rnorm(n), nrow=n, ncol=1, byrow=T)
        x_true = solve(A, b)
        results = GE_w_total_pivoting(A, b); x_wwx = results$X
        err = sum((x_true - x_wwx)^2)
        print(sprintf('ri= %5d; n= %5d; err= %10.8f', ri, n, err))
    }
}