lu_decompose.c File Reference

LU decomposition of a square matrix More...

#include <stdio.h>
#include <stdlib.h>
#include <time.h>

Include dependency graph for lu_decompose.c:

Functions
int	lu_decomposition (double **A, double **L, double **U, int mat_size)
	Perform LU decomposition on matrix. More...
void	display (double **A, int N)
	Function to display square matrix.
int	main (int argc, char **argv)
	Main function.

Detailed Description

LU decomposition of a square matrix

Author

Krishna Vedala

Function Documentation

lu_decomposition()

int lu_decomposition	(	double **	A,
		double **	L,
		double **	U,
		int	mat_size
	)

Perform LU decomposition on matrix.

Parameters

[in]	A	matrix to decompose
[out]	L	output L matrix
[out]	U	output U matrix
[in]	mat_size	input square matrix size

{
    int row, col, j;
 
    // regularize each row
    for (row = 0; row < mat_size; row++)
    {
        // Upper triangular matrix
#ifdef _OPENMP
#pragma omp for
#endif
        for (col = row; col < mat_size; col++)
        {
            // Summation of L[i,j] * U[j,k]
            double lu_sum = 0.;
            for (j = 0; j < row; j++) lu_sum += L[row][j] * U[j][col];
 
            // Evaluate U[i,k]
            U[row][col] = A[row][col] - lu_sum;
        }
 
        // Lower triangular matrix
#ifdef _OPENMP
#pragma omp for
#endif
        for (col = row; col < mat_size; col++)
        {
            if (row == col)
            {
                L[row][col] = 1.;
                continue;
            }
 
            // Summation of L[i,j] * U[j,k]
            double lu_sum = 0.;
            for (j = 0; j < row; j++) lu_sum += L[col][j] * U[j][row];
 
            // Evaluate U[i,k]
            L[col][row] = (A[col][row] - lu_sum) / U[row][row];
        }
    }
 
    return 0;
}