程序代码：

//在VB中写的，转成c大概是这样
//LU分解，由系数矩阵a得到l和u
//需要先初始化对角阵l和u，n为阶数

...
  n=n-1;
  for(p=0;p<=n;p++)
    { for(m=p;m<=n;m++)
        { sum=0;
      for(k=0;k<=n;k++)
            { if (k!=p) sum=sum+l[p][k]*u[k][m]; }
          u[p][m]=(a[p][m]-sum)/l[p][p];
         }
       for(m=p+1;m<=n;m++)
         { sum=0;
       for(k=0;k<=n;k++)
             { if(k!=p) sum=sum+l[m][k]*u[k][p];}
           l[m][p]=(a[m][p]-sum)/u[p][p];
          }
     }

  sum=0;
  for(k=0;k<=n;k++)
     { if(k!=n) sum=sum+l[n][k]*u[k][n];}
    u[n][n]=(a[n][n]-sum)/l[n][n];

 ...


 //后面解线性方程，很容易自己考虑

程序代码：

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

void main()
  {
    int i,j,n,p,m,k;
    int **a;
    float **l,**u;
    float sum;
    system("cls");
    printf("Please input matrix a's order:");
    scanf("%d",&n);
    printf("Matrix a's order is %d.\n",n);
/*Input n order 2-d array a*/
    printf("Please updata matrix a\n");
    a=(int**)malloc(sizeof(int*)*n);
    for(i=0;i<n;i++)
      a[i]=(int*)malloc(sizeof(int)*n);

    for(i=0;i<=n-1;i++)
    {  for(j=0;j<=n-1;j++)
          {  printf("a[%d][%d]=",i,j);
           scanf(" %d",&a[i][j]);
        }
    }

    printf("Matrix a is:\n");
    for(i=0;i<=n-1;i++)
      { for(j=0;j<=n-1;j++)printf("%d\t",a[i][j]);
        printf("\n");
      }
/*Create n order 2-D array l,u*/
    l=(float**)malloc(sizeof(float*)*n);
    for(i=0;i<n;i++)
      l[i]=(float*)malloc(sizeof(float)*n);

    u=(float**)malloc(sizeof(float*)*n);
    for(i=0;i<n;i++)
      u[i]=(float*)malloc(sizeof(float)*n);
    n--;
/*Initialize l and u */
  for(i=0;i<=n;i++)
    { for(j=0;j<=n;j++)
    { l[i][j]=i<j?0:1;
      u[i][j]=i>j?0:1;
         }
     }
  for(p=0;p<=n;p++)
    { for(m=p;m<=n;m++)
        { sum=0;
      for(k=0;k<=n;k++)
        { if (k!=p) sum+=l[p][k]*u[k][m]; }
          u[p][m]=(a[p][m]-sum)/l[p][p];
         }
       for(m=p+1;m<=n;m++)
         { sum=0;
       for(k=0;k<=n;k++)
         { if(k!=p) sum+=l[m][k]*u[k][p];}
           l[m][p]=(a[m][p]-sum)/u[p][p];
          }
     }

  sum=0;
  for(k=0;k<=n;k++)
     { if(k!=n) sum=sum+l[n][k]*u[k][n];}
    u[n][n]=(a[n][n]-sum)/l[n][n];

  printf("\nMatrix l is:\n");
  for(i=0;i<=n;i++)
    { for (j=0;j<=n;j++)
    { printf("%5.2f\t",l[i][j]); }
      printf("\n");
    }
  printf("\n");
  printf("Matrix u is:\n");
  for(i=0;i<=n;i++)
    { for (j=0;j<=n;j++)
    { printf("%5.2f\t",u[i][j]); }
      printf("\n");
    }
/*Free memory*/
    for(i-0;i<n+1;i++)
    { free(a[i]);free(l[i]);free(u[i]);}
    free(a);free(l);free(u);
   getch();
  }