数学研发论坛

 找回密码
 欢迎注册
楼主: liangbch

[提问] 征集求解线性方程组的数学软件

[复制链接]
 楼主| 发表于 2010-3-1 18:03:30 | 显示全部楼层
安装了一个matlab 6.5,并尝试写了一个命令行脚本,现分享给大家。

用一个文本编辑器创建一个名为test1.m的文件,保存在C:\MATLAB6p5\work,文件内容如下:
  1. A=[1,1,-5,-1,1;8,5,-1,0,4;3,-7,-5,-1,1;2,1,8,3,-1;0,6,2,7,-2];
  2. b=[9;2;4;-8;15];
  3. xx=A\b;

  4. diary r.txt % 将所有输出保存至r.txt文件
  5. diary on %打开diary开关
  6. prot='A\b is ----------------';

  7. disp('A=');
  8. disp(A);
  9. disp('b=');
  10. disp(b);
  11. disp(prot);
  12. disp(xx);
  13. diary off; %关闭diary开关
  14. exit  % 退出matlab
复制代码
在命令行运行matlab -r test1,则启动matlab,并求解方程组,保存结果到r.txt。

评分

参与人数 1鲜花 +2 收起 理由
wayne + 2 赞一个

查看全部评分

毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2010-3-1 18:50:17 | 显示全部楼层
赞一个先,

不过,我得说我们一般都不是这样用MATLAB的~~
matlab -r test1最终还是会启动MATLAB图形界面的,到不如直接在的GUI界面里面的command windows里面输入三行代码来的爽快:
  1. A=[1,1,-5,-1,1;8,5,-1,0,4;3,-7,-5,-1,1;2,1,8,3,-1;0,6,2,7,-2];
  2. b=[9;2;4;-8;15];
  3. xx=A\b;
复制代码
另外 MATLAB可以和c++进行混合编程的。。。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2010-3-1 19:02:41 | 显示全部楼层
我的解法是消元法,将系数矩阵化为上三角矩阵。
liangbch 发表于 2010-3-1 16:31

一般的计算方法书籍都会讲到线性方程组求解的多种算法,消元法最基础。
我记得还有一个LU分解的方法
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2010-3-1 19:26:20 | 显示全部楼层
我来给一个gsl的版本:
  1. #include <stdio.h>
  2. #include <gsl/gsl_linalg.h>   
  3. int main ()
  4.      {
  5.        double a_data[] = {1,1,-5,-1,1,8,5,-1,0,4,3,-7,-5,-1,1,2,1,8,3,-1,0,6,2,7,-2};   
  6.        double b_data[] = {9,2,4,-8,15};
  7.           
  8.        gsl_matrix_view m = gsl_matrix_view_array (a_data, 5, 5);  
  9.        gsl_vector_view b= gsl_vector_view_array (b_data, 5);   
  10.        gsl_vector *x = gsl_vector_alloc (5);      
  11.        int s;   
  12.        gsl_permutation * p = gsl_permutation_alloc (5);   
  13.        gsl_linalg_LU_decomp (&m.matrix, p, &s);     
  14.        gsl_linalg_LU_solve (&m.matrix, p, &b.vector, x);
  15.        printf ("x = \n");
  16.        gsl_vector_fprintf (stdout, x, "%g");
  17.      
  18.        gsl_permutation_free (p);
  19.        gsl_vector_free (x);
  20.        return 0;
  21.      }
复制代码
编译运行:
H:\tkpng>gcc linear.c -Iinclude -Llib -lgsl -lgslcblas -o gsl

H:\tkpng>gsl.exe
x =
1.5024
1.0006
-2.30719
0.706587
-4.33234
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2010-3-1 22:24:00 | 显示全部楼层
本帖最后由 liangbch 于 2010-3-1 22:31 编辑

我成功地使用VC++6.0跑通了楼上的程序,下面列出我的方法

1)下载WinGsl-Lib-1.4.02.zip from
  http://www6.in.tum.de/~kiss/WinGsl.htm
(2)解压到C:\,WinGsl中的库文件的扩展名是.lib,可用于Visual C++ 6.0 的编译器
(3)设置Visual C++ 6.0编译环境:
  一、将 C:\WinGsl\Bin中的WinGsl.dll和WinGslD.dll复制到C:\Program Files\Microsoft Visual Studio\VC98\Bin下;
  二、打开VC,tools->options->directories,将解压目录下的lib,gsl分别加入到库文件和头文件的搜索路径中,
   即, 将C:\WINGSL加入到Include files 搜索路径,将C:\WINGSL\LIB加入到Library files.

(4) 在新建的工程中将你所使用到的库文件加入到 Link->Object/library modules 中,我在debug模式下添加的是WinGslLib_md.lib,release 模式下添加的是WinGslLib_md.lib,不过好像用错了文件也可以工作。

参考自:DIY部落(http://www.diybl.com/course/3_pr ... 2008326/107387.html)
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2010-3-1 22:25:32 | 显示全部楼层
还有一个问题,14#的数据类型是double型,不知能否使用GMP中的高精度数据类型。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2010-3-1 22:29:03 | 显示全部楼层
12# wayne

11楼的命令行用法 适用求解大量方程组,如果我们有10000个方程组,写一个批处理,采用11楼的的做法很快可以得到答案。如果用GUI则麻烦的多。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2010-3-2 10:15:38 | 显示全部楼层
 在我的应用中,有好多方程组是线性相关的。对于这样的方程组,在消元的过程中,某一个方程的所有系数会同时为0。这时,程序就可以确定方程组无解。
 我特意挑选了一个线性相关的方程组,使用 14# 的代码,用GLS求解,结果发现他依然给出一个解,当然这个解错误的。

 今天凌晨,检查了一下我的程序,在修复了一个bug后发现,整数溢出的概率非常低,再也没有看到整数溢出的报告。所以这个主题贴实际上对我的应用没有多大帮助。不过这个主题贴让我们学到很多,欢迎继续讨论,谢谢所有回复的人。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2010-3-2 10:45:03 | 显示全部楼层
本帖最后由 wayne 于 2010-3-2 10:50 编辑
12# wayne

11楼的命令行用法 适用求解大量方程组,如果我们有10000个方程组,写一个批处理,采用11楼的的做法很快可以得到答案。如果用GUI则麻烦的多。
liangbch 发表于 2010-3-1 22:29


MATLAB是一个很强大的工具,即便是GUI也可以批处理的啊

我自己写了一个批处理的脚本,数据是随机产生的,有1000个方程组;你只需把代码保存到一个test.m文件里,在文件的焦点下按F5运行,或者在command windows里直接输入test也可以运行
  1. %% 产生测试数据:1000个随机的线性方程组 Ax=b, A为5*5方阵,b为5*1列向量
  2. for ii=1:1000
  3. data(ii).A=randi([-2,8],5);
  4. data(ii).b=randi([-2,8],5,1);
  5. end
  6. %% 批量解线性方程组,解统一存在结构表 data的x域里
  7. for ii=1:1000
  8. data(ii).x=data(ii).A\data(ii).b;
  9. end

  10. %% 把数据写入文件leqs.txt
  11. file='leqs.txt'
  12. for ii=1:10
  13. dlmwrite(file, 'A=', 'delimiter', '','-append')
  14. dlmwrite(file, data(ii).A, 'delimiter', '\t','-append')
  15. dlmwrite(file, 'b=', 'delimiter', '','-append')
  16. dlmwrite(file, data(ii).b,'delimiter', '\t', '-append')
  17. dlmwrite(file, 'x=', 'delimiter', '','-append')
  18. dlmwrite(file, data(ii).x, 'delimiter', '\t','-append')
  19. dlmwrite(file, '================================','delimiter', '','-append')
  20. end
复制代码
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2010-3-2 11:40:24 | 显示全部楼层
19# wayne

  你可能没有理解我的意思,我的问题是需要求解的方程组有好多个, 比如说有1000个,我已经用别的程序生成1000个m文件,如果你通过GUI的话,你需要
  1.启动MATLAB
    2.手动输入1000个命令来执行者1000个.m文件,如果我使用批处理方式来处理类似的问题,则可以创建一个名为 all.bat 批处理(内容如下),然后只需要在命令上输入all就可以了。

matlab -r 00001.m
matlab -r 00002.m
matlab -r 00003.m
matlab -r 00004.m
matlab -r 00005.m
matlab -r 00006.m
matlab -r 00006.m
matlab -r 00007.m
matlab -r 00008.m
matlab -r 00009.m
matlab -r 00010.m
...
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
您需要登录后才可以回帖 登录 | 欢迎注册

本版积分规则

小黑屋|手机版|数学研发网 ( 苏ICP备07505100号 )

GMT+8, 2022-5-20 10:02 , Processed in 0.073520 second(s), 15 queries .

Powered by Discuz! X3.4

© 2001-2017 Comsenz Inc.

快速回复 返回顶部 返回列表