声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 3282|回复: 6

[经典算法] 微分求积法(Differential Quadrature Method)

[复制链接]
发表于 2007-9-9 14:09 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
微分求积法,简称 DQ,是一种求解常微分和偏微分方程(组)的有效方法。将常微分和偏微分方程(组)经过DQ离散成线性或非线性方程组后可通过Newmark 和 Newton-Rapuson 求解,也可以直接使用Runge-Kutta法。但是本人在使用这些方法是遇到了困难,主要的问题是边界条件的处理,之后编程和算法上,因为它的离散过程是简单的!很希望和正在或曾经使用此方法的朋友交流,并虚心请教!
回复
分享到:

使用道具 举报

发表于 2007-10-12 09:39 | 显示全部楼层
参考下面的程序吧

  1. % Determination of the coefficients in DQM

  2. %   19/03/2006 By Dr. Yang J., Modified by Mr. Xiang H.J.

  3. %   Department of BC, City University of Hong Kong, Hong Kong

  4. %   [C,X]=DQC(N,M,KOP)

  5. %   KOP, Keyoption. 1--uniform space pattern; 2,3--cosine space pattern

  6. %   N--Total grid point; M--Order

  7. function [C,X]=DQC(N,M,KOP)

  8. %1. Generation of grid point coordinates in X-axis

  9. for i=1:N

  10.     X1(i)=(i-1)/(N-1);

  11.     X2(i)=0.5*(1-cos(pi*(i-1)/(N-1)));

  12. end

  13. X3(1)=0;X3(2)=0.00000001;%0.00001;

  14. for i=3:N

  15.      X3(i)=0.5*(1-cos(pi*(i-1)/(N-1)));

  16. end

  17. X4(1)=0;X4(2)=0.00001;X4(N-1)=0.99999;X4(N)=1.0;

  18. for i=3:N-2

  19.      X4(i)=0.5*(1-cos(pi*(i-2)/(N-3)));

  20. end

  21. X5(1)=0;X5(2)=0.00001;X5(3)=0.001;

  22. X5(N-2)=0.999;X5(N-1)=0.99999;X5(N)=1.0;

  23. for i=4:N-3

  24.      X5(i)=0.5*(1-cos(pi*(i-3)/(N-5)));

  25. end

  26. if KOP==1

  27.     X=X1;

  28. elseif KOP==2

  29.     X=X2;

  30. elseif KOP==3

  31.     X=X3;

  32. elseif KOP==4

  33.     X=X4;

  34. else

  35.     X=X5;

  36. end

  37. %2.   Producing the weighting coefficients

  38. %2.1 Producing the first-order weighting coefficients

  39. C=zeros(N,N,M);

  40. for i=1:N

  41.     Q(i)=1.0;

  42.     for j=1:N

  43.        if i~=j

  44.           Q(i)=Q(i)*(X(i)-X(j));

  45.        end

  46.     end

  47. end

  48. for i=1:N  

  49.     for j=1:N

  50.         if i~=j

  51.           C(i,j,1)=Q(i)/(X(i)-X(j))/Q(j);

  52.           C(i,i,1)=C(i,i,1)-C(i,j,1);

  53.        end

  54.     end

  55. end

  56. %2.2 Producing the High-order weighting coefficients

  57. if M>1

  58.    for k=2:M

  59.      for i=1:N

  60.       for j=1:N

  61.        if i~=j

  62.           C(i,j,k)=k*(C(i,i,k-1)*C(i,j,1)-C(i,j,k-1)/(X(i)-X(j)));

  63.           C(i,i,k)=C(i,i,k)-C(i,j,k);

  64.        end

  65.       end

  66.      end

  67.    end

  68. end
复制代码


来自:http://hi.baidu.com/hjxiang/blog ... 43902fb93820fb.html

评分

1

查看全部评分

 楼主| 发表于 2007-10-26 16:55 | 显示全部楼层

回复 #2 风花雪月 的帖子

谢谢!这个我已经做过了,而且已经用C重新写了一遍。现在的问题就是求解二阶微分方程了,我要用4阶龙格库塔法!正在思考。其实我倒自己写了一个,但是不是很好用。算起来很慢,不知道为什么!楼主,能不能帮着找一找更好的程序!
发表于 2007-10-26 23:33 | 显示全部楼层

回复 #3 su200703 的帖子

我想问一下,为什么不直接用ode45这些,命令直接求解,这个方法精度更高吗?
发表于 2007-12-6 09:26 | 显示全部楼层
原帖由 su200703 于 2007-10-26 16:55 发表
谢谢!这个我已经做过了,而且已经用C重新写了一遍。现在的问题就是求解二阶微分方程了,我要用4阶龙格库塔法!正在思考。其实我倒自己写了一个,但是不是很好用。算起来很慢,不知道为什么!楼主,能不能帮着找 ...


rk法本来的效率就不高,龙格库塔法的程序网络上有很多,你可以借鉴,本站搜索一下就能找到
发表于 2007-12-6 09:28 | 显示全部楼层
原帖由 无水1324 于 2007-10-26 23:33 发表
我想问一下,为什么不直接用ode45这些,命令直接求解,这个方法精度更高吗?


ode45是变步长的,说不定楼主想用定步长或者用自己学习编程也未可知
发表于 2012-4-15 15:25 | 显示全部楼层
还与人看这个帖子吗?程序看不懂,谁给发个结合实例的?
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-12-23 10:50 , Processed in 0.066933 second(s), 20 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

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