又一些c语言“线性代数”函数(有更新,并有一些疑问向大家请教,永不过期)

linux_niu

用 gdb 跟着跑了一遍, 其实错误很明显的, 是浮点数精度的问题.

当 rdechl_mtr 函数中 j 达到 3 的时候, p 中保存的已经是结果式了, 然后在对 0 进行计数时, linux_niu 兄是用

2. if (p[i][j] == 0) zero++;

复制代码

这样的语句. 如果 p[j] 是整形量, 那么这样完全正确, 但是不幸的是, p[j] 是浮点型, 而且其结果又是经过几轮相减后得到的一个接近 0 的小量, 但是又不是 0, 这就是浮点数误差.

在 j == 3 的时候, p[2][j] 的值是 -5.921e-16, 实际上就是零了, 但是又无法满足 p[2][j] == 0 的判断

这种情况下, 我们可以把此语句改成

2. if (fabs(p[i][j]) < 0.0001) zero++;

复制代码

0.0001 也可以换成一个稍大一点或稍小一点的值, 只要保证 [1 它远小于可能出现在矩阵中的最小的非零绝对值], [2 又远大于可能达到的累积误差], 这两个条件即可.

实际上, 以楼主的目的来看, 用分数来表示一个数可能更合理一些, 不过这要进行一些先期的封装.

P.S. 建议 linux_niu 兄对代码进行以下几点改进:
* 对关键变量的作用进行说明
* 在关键代码处加上注释
* 用空行分隔语义段
* for(r=j,tmp=p[j];r<n;r++)p[r]-=p[k][r]*tmp; 这样的代码写成

2. tmp = p[i][j];
3. for (r = j; r < n; r++)
4.     p[i][r] -= p[k][r] * tmp;

复制代码

的形式更易于阅读
* 统一函数接口风格. 比如将名字都改成 matrix_make, matrix_free, matrix_exch, matrix_rdechl, array_make, 这样一看就知道哪些是 matrix 相关的, 哪些是 array 相关的
* malloc 分配的内存要记得释放, 比如 exch_mtr() 中的 t
* 更根本的, 在一个函数中申请一块内存给外部使用, 比如 rdechl_mtr() 中的 p, 这种方式非常不好, 说到底还是代码结构的问题. 在另一帖中 realtang 等兄弟们已经推荐过一些好的资料了, 楼主不妨参考一下.