高斯消元

输入线性方程组的维数，然后随即生成一定有解的线性方程组的增广矩阵，求出解，然后输出时间和方程组的解，以及和标准解的误差（其实就是和标准解的方差）。

 1 #include <cstdio>
 2 #include <algorithm>
 3 #include <cstdlib>
 4 #include <iostream>
 5 #include <cstring>
 6 #include <cmath>
 7 #include <ctime>
 8 using namespace std;
 9 const double eps=1e-9;
10 const int MAX=20000;
11 double ans[MAX]; 
12 int n;
13 double **inputdata; double *result;
14 /*
15     输出标准解
16 */
17 void printresult() {
18    for(int i=0;i<n;i++)
19        printf("%.4lf ",result[i]);
20    printf("
");
21 }
22 /*
23     打印出生成的的增广矩阵
24 */
25 void printb() {
26      for(int i=0;i<n;i++)
27        printf("%.4lf ",inputdata[i][n]);
28    printf("
"); 
29 }
30 /*
31     随机生成增广矩阵
32 */
33 void pro_inputdata()
34 {
35     int i,j,k;
36     double res=0;
37   printf("please input the number of the element:
");
38     scanf("%d",&n);
39     srand(time(0)); 
40     result=(double*)malloc(n*sizeof(double));
41     inputdata=(double**)malloc(n*sizeof(double));
42     for(i=0;i<n;i++)
43          result[i]=(double)(rand()%10000)/1000.0;
44     for(i=0;i<n;i++) {
45         inputdata[i]=(double*)malloc((n+1)*sizeof(double));
46     res=0;
47       for(j=0;j<n;j++) {
48             k=(rand()%2==0)?1:(-1);
49             inputdata[i][j]=k*(double)(rand()%10000)/1000.0;
50             res=res+inputdata[i][j]*result[j];
51         }
52         inputdata[i][n]=res;
53     }
54 }
55 
56 /*
57     高斯消元函数
58 */
59 inline void solve(double **a, double ans[], int n) {
60     for (int i = 0; i < n; ++i) {
61         for (int j = i; j < n; ++j)
62             if (fabs(a[j][i]) > eps) {
63                 for (int k = i; k <= n; ++k) swap(a[j][k], a[i][k]);
64                 break;
65             }
66         for (int j = i+1; j < n; ++j) {
67             if (fabs(a[j][i]) > eps) {
68                 double tmp = a[j][i]/a[i][i];
69                 for (int k = i; k <= n; ++k) a[j][k] -= tmp * a[i][k];
70             }
71         }
72     }
73     for (int i = n - 1; i >= 0; --i) {
74         double sum = 0.0;
75         for  (int j = n - 1; j > i; --j)
76             sum += ans[j] * a[i][j];
77         ans[i] = (a[i][n] - sum) / a[i][i];
78     }
79     for (int i = 0; i < n; ++i) printf("%.4lf ", ans[i]);
80     printf("
");
81     return;
82 }
83 int main(void) {
84     freopen("in.txt", "r", stdin);
85     clock_t begin, end;
86     begin = clock();
87     pro_inputdata();
88     solve(inputdata, ans, n);
89     end = clock(); 
90     printf("Time is: %lf seconds
",(double)(end-begin)/CLOCKS_PER_SEC);
91     double sum=0;
92     for (int i = 0; i < n; ++i) {
93         sum += (ans[i]-result[i])*(ans[i]-result[i]);
94     }
95     printf("Error is %.4lf
",sqrt(sum));
96 
97     return 0;
98 }

囧，等下在贴代码。＝＝

刚才系统剪切板里面的内容不能粘贴到浏览器里面了，注销了一下就可以了=_=