❶ 如何用MATLAB編寫的拉格朗日插值演算法的程序、二階龍格-庫塔方法的程序和SOR迭代法的程序
拉格朗日function y=lagrange(x0,y0,x)
n=length(x0);m=length(x);
for i=1:m
z=x(i);
s=0.0;
for k=1:n
p=1.0;
for j=1:n
if j~=k
p=p*(z-x0(j))/(x0(k)-x0(j));
end
end
s=p*y0(k)+s;
end
y(i)=s;
end SOR迭代法的Matlab程序
function [x]=SOR_iterative(A,b)
% 用SOR迭代求解線性方程組,矩陣A是方陣
x0=zeros(1,length(b)); % 賦初值
tol=10^(-2); % 給定誤差界
N=1000; % 給定最大迭代次數
[n,n]=size(A); % 確定矩陣A的階
w=1; % 給定鬆弛因子
k=1;
% 迭代過程
while k<=N
x(1)=(b(1)-A(1,2:n)*x0(2:n)')/A(1,1);
for i=2:n
x(i)=(1-w)*x0(i)+w*(b(i)-A(i,1:i-1)*x(1:i-1)'-A(i,i+1:n)*x0(i+1:n)')/A(i,i);
end
if max(abs(x-x0))<=tol
fid = fopen('SOR_iter_result.txt', 'wt');
fprintf(fid,'\n********用SOR迭代求解線性方程組的輸出結果********\n\n');
fprintf(fid,'迭代次數: %d次\n\n',k);
fprintf(fid,'x的值\n\n');
fprintf(fid, '%12.8f \n', x);
break;
end
k=k+1;
x0=x;
end
if k==N+1
fid = fopen('SOR_iter_result.txt', 'wt');
fprintf(fid,'\n********用SOR迭代求解線性方程組的輸出結果********\n\n');
fprintf(fid,'迭代次數: %d次\n\n',k);
fprintf(fid,'超過最大迭代次數,求解失敗!');
fclose(fid);
end Matlab中龍格-庫塔(Runge-Kutta)方法原理及實現龍格-庫塔(Runge-Kutta)方法是一種在工程上應用廣泛的高精度單步演算法。由於此演算法精度高,採取措施對誤差進行抑制,所以其實現原理也較復雜。該演算法是構建在數學支持的基礎之上的。龍格庫塔方法的理論基礎來源於泰勒公式和使用斜率近似表達微分,它在積分區間多預計算出幾個點的斜率,然後進行加權平均,用做下一點的依據,從而構造出了精度更高的數值積分計算方法。如果預先求兩個點的斜率就是二階龍格庫塔法,如果預先取四個點就是四階龍格庫塔法。一階常微分方程可以寫作:y'=f(x,y),使用差分概念。
(Yn+1-Yn)/h= f(Xn,Yn)推出(近似等於,極限為Yn')
Yn+1=Yn+h*f(Xn,Yn)
另外根據微分中值定理,存在0<t<1,使得
Yn+1=Yn+h*f(Xn+th,Y(Xn+th))
這里K=f(Xn+th,Y(Xn+th))稱為平均斜率,龍格庫塔方法就是求得K的一種演算法。
利用這樣的原理,經過復雜的數學推導(過於繁瑣省略),可以得出截斷誤差為O(h^5)的四階龍格庫塔公式:
K1=f(Xn,Yn);
K2=f(Xn+h/2,Yn+(h/2)*K1);
K3=f(Xn+h/2,Yn+(h/2)*K2);
K4=f(Xn+h,Yn+h*K3);
Yn+1=Yn+h*(K1+2K2+2K3+K4)*(1/6);
所以,為了更好更准確地把握時間關系,應自己在理解龍格庫塔原理的基礎上,編寫定步長的龍格庫塔函數,經過學習其原理,已經完成了一維的龍格庫塔函數。
仔細思考之後,發現其實如果是需要解多個微分方程組,可以想像成多個微分方程並行進行求解,時間,步長都是共同的,首先把預定的初始值給每個微分方程的第一步,然後每走一步,對多個微分方程共同求解。想通之後發現,整個過程其實很直觀,只是不停的逼近計算罷了。編寫的定步長的龍格庫塔計算函數:
function [x,y]=runge_kutta1(ufunc,y0,h,a,b)%參數表順序依次是微分方程組的函數名稱,初始值向量,步長,時間起點,時間終點(參數形式參考了ode45函數)
n=floor((b-a)/h);%求步數
x(1)=a;%時間起點
y(:,1)=y0;%賦初值,可以是向量,但是要注意維數
for ii=1:n
x(ii+1)=x(ii)+h;
k1=ufunc(x(ii),y(:,ii));
k2=ufunc(x(ii)+h/2,y(:,ii)+h*k1/2);
k3=ufunc(x(ii)+h/2,y(:,ii)+h*k2/2);
k4=ufunc(x(ii)+h,y(:,ii)+h*k3);
y(:,ii+1)=y(:,ii)+h*(k1+2*k2+2*k3+k4)/6;
%按照龍格庫塔方法進行數值求解
end
調用的子函數以及其調用語句:
function dy=test_fun(x,y)
dy = zeros(3,1);%初始化列向量
dy(1) = y(2) * y(3);
dy(2) = -y(1) + y(3);
dy(3) = -0.51 * y(1) * y(2);
對該微分方程組用ode45和自編的龍格庫塔函數進行比較,調用如下:
[T,F] = ode45(@test_fun,[0 15],[1 1 3]);
subplot(121)
plot(T,F)%Matlab自帶的ode45函數效果
title('ode45函數效果')
[T1,F1]=runge_kutta1(@test_fun,[1 1 3],0.25,0,15);%測試時改變test_fun的函數維數,別忘記改變初始值的維數
subplot(122)
plot(T1,F1)%自編的龍格庫塔函數效果
title('自編的 龍格庫塔函數')
❷ sor迭代演算法的matlab程序
function [x k]=EqtsSOR(A,b,x0,omiga,eps)
%超鬆弛(SOR,Successive Over-Relaxation)迭代法求解線性方程組Ax=b
%[x k]=EqtsSOR(A,b,x0,eps)
%x:解向量,列向量
%k:迭代次數
%A:系數矩陣
%b:列向量
%x0:迭代初始值,列向量
%omiga:鬆弛因子,可預設,預設值為1,即為GS迭代法
%eps:誤差限,可預設,預設值為0.5e-6
%
%應用舉例:
%A=[4 3 0;3 4 -1;0 -1 4];b=[24;30;-24];x0=[1;1;1];omiga=1.25;
%[x k]=EqtsSOR(A,b,x0,omiga,0.5e-6)
%x=EqtsSOR(A,b,x0)
if nargin==4
eps=0.5e-6;
end
if nargin==3
omiga=1;
eps=0.5e-6;
end
%檢查輸入參數
n=length(b);
if size(A,1) ~= n || n ~= length(x0)
disp('輸入參數有誤!');
x=' ';
k=' ';
return;
end
%迭代求解
k=0;
x=zeros(n,1);
while 1
k=k+1;
for i=1:n
z=0;
for j=1:i-1
z=z+A(i,j)*x(j);
end
for j=i+1:n
z=z+A(i,j)*x0(j);
end
x(i)=(1-omiga)*x0(i)+omiga*(b(i)-z)/A(i,i);
end
if norm(x-x0)<=eps || k==30
break;
end
x0=x;
end
if k==30
disp('迭代次數太多!')
x=' ';
end
return;
❸ 請c/c++高手看一下關於sor程序的問題
錯誤不在sor,在這里:
for(j=0;j<m;j++)
{
for(l=0;l<n;l++)
{
a[j][l]=rand();
b[j][l]=rand();
c[j][l]=rand();
d[j][l]=rand();
e[j][l]=rand();
f[j][l]=rand();
u[j][l]=0.0;
}
}
懷疑你的dmatrix函數寫的有問題,代碼寫得很亂,不知道在干什麼。估計是內存分配有問題。我這里的運行時錯誤是訪問了非法地址0。環境是VS 2005.
我把你的程序稍微修改了一下,可以看得比較清楚錯誤的地方。
#include <stdlib.h>
#include <stdio.h>
//#include "nrutil.h"
#include <time.h>
#include <math.h>
#include <stddef.h>
#define MAXITS 1000
#define NR_END 1
#define FREE_ARG char*
#define EPS 1.0e-5
void nrerror(char error_text[]);
void sor(double **a,double **b,double **c,double **d,double **e,double **f,double **u,int jmax,double rjac)
{
int ipass,j,jsw,l,lsw,n;
double anorm,anormf=0.0,omega=1.0,resid;
printf("enter sor\n");
for(j=2;j<jmax;j++)
for(l=2;l<jmax;l++)
anormf+=fabs(f[j][l]);
for(n=1;n<=MAXITS;n++)
{
anorm=0.0;
jsw=1;
for(ipass=1;ipass<=2;ipass++)
{
lsw=jsw;
for(j=2;j<jmax;j++)
{
for(l=lsw+1;l<jmax;l+=2)
{
resid=a[j][l]*u[j+1][l]+b[j][l]*u[j-1][l]+c[j][l]*u[j][l+1]+d[j][l]*u[j][l-1]+e[j][l]*u[j][l]-f[j][l];
anorm+=fabs(resid);
u[j][l]-=omega*resid/e[j][l];
}
lsw=3-lsw;
}
jsw=3-jsw;
omega=(n==1&&ipass==1?1.0/(1.0-0.5*rjac*rjac):1.0/(1.0-0.25*rjac*rjac*omega));
}
if(anorm<EPS*anormf)
return;
}
nrerror("MAXITS exceeded");
}
void nrerror(char error_text[])
{
fprintf(stderr,"Numerical Recipes run-time error...\n");
fprintf(stderr,"%s\n",error_text);
fprintf(stderr,"...now exiting to system...\n");
exit(1);
}
double **dmatrix(long nrl, long nrh, long ncl, long nch)
{
long i, nrow=nrh-nrl+1,ncol=nch-ncl+1;
double **m;
m=(double **) malloc((size_t)((nrow+NR_END)*sizeof(double*)));
if (!m) nrerror("allocation failure 1 in matrix()");
m += NR_END;
m -= nrl;
m[nrl]=(double *) malloc((size_t)((nrow*ncol+NR_END)*sizeof(double)));
if (!m[nrl]) nrerror("allocation failure 2 in matrix()");
m[nrl] += NR_END;
m[nrl] -= ncl;
for(i=nrl+1;i<=nrh;i++) m[i]=m[i-1]+ncol;
printf("dmatrix succeeded!\n");
return m;
}
main()
{
double **a,**b,**c,**d,**e,**f,**u,rjac;
int j,l,jmax;
double pi=3.1415926;
long m=10,n=10;
a=dmatrix(1,m,1,n);
b=dmatrix(1,m,1,n);
c=dmatrix(1,m,1,n);
d=dmatrix(1,m,1,n);
e=dmatrix(1,m,1,n);
f=dmatrix(1,m,1,n);
u=dmatrix(1,m,1,n);
printf("%x, %x, %x, %x, %x, %x, %x\n", (int)a, (int)b, (int)c, (int)d, (int)e, (int)f, (int)u);
scanf("%d",&jmax);
rjac=cos(pi/jmax);
srand((unsigned)time(NULL));
printf("init begin...\n");
for(j=0;j<m;j++)
{
for(l=0;l<n;l++)
{
a[j][l]=rand();
b[j][l]=rand();
c[j][l]=rand();
d[j][l]=rand();
e[j][l]=rand();
f[j][l]=rand();
u[j][l]=0.0;
}
}
printf("init finished!\n");
sor(a,b,c,d,e,f,u,jmax,rjac);
for(j=0;j<10;j++)
{
for(l=0;l<10;l++)
{
printf("%f\t",u[j][l]);
}
printf("\n");
}
getch();
}
❹ SOR迭代法 求MATLAB程序
迭代法 matlab實現代碼如下
function [x,n] = jacobi(A,b,x0,eps,varargin)
if nargin ==3
eps = 1.0e-6;
M = 200;
elseif nargin<3
disp('輸入參數數目不足3個');
return
elseif nargin ==5
M = varargin{1};
end
D = diag(diag(A)); %%求A的對角矩陣
L = -tril(A,-1); %%求A的下三角矩陣
U = -triu(A,1);%%求A的上三角矩陣
B = D(L+U);
f = D;
x = B*x0+f;
n = 1;%迭代次數
while norm(x-x0)>=eps
x0 = x;
x = B*x0+f
n = n+1;
if(n>=M)
disp('Warning:迭代次數太多,可能不收斂!')
return;
end
end
運行效果如下:
,A對稱正定,且0<ω<2,則解Ax=b的SOR方法收斂。
❺ 求助,MATLAB關於sor超鬆弛演算法
function [x,n]=SOR(A,b,x0,w,eps,M)
if nargin==4
eps= 1.0e-6;
M = 200;
elseif nargin<4
error
return
elseif nargin ==5
M = 200;
end
if(w<=0 || w>=2)
error;
return;
end
D=diag(diag(A)); %求A的對角矩陣
L=-tril(A,-1); %求A的下三角陣
U=-triu(A,1); %求A的上三角陣
B=inv(D-L*w)*((1-w)*D+w*U);
f=w*inv((D-L*w))*b;
x=B*x0+f;
n=1; %迭代次數
while norm(x-x0)>=eps
x0=x;
x =B*x0+f;
n=n+1;
if(n>=M)
disp('Warning: 迭代次數太多,可能不收斂!');
return;
end
end
這個程序無誤。
w=1時,它們倆是等價的。它們得到的結果如果精度要求低可能一樣,但是精度高了後就會不一樣的。
❻ sor是誰提出的
D. M. Young於20世紀70年代提出逐次超鬆弛(Successive Over Relaxation)迭代法,簡稱SOR方法。
它是為了解決大規模系統的線性等式提出來的,在GS法基礎上為提高收斂速度,採用加權平均而得到的新演算法。由於超鬆弛迭代法公式簡單,編製程序容易,很多工程學、計算數學中都會應用超鬆弛迭代方法。使用超鬆弛迭代法的關鍵在於選取合適的鬆弛因子,如果鬆弛因子選取合適,則會大大縮短計算時間。
SOR方法的思路為:如果能夠簡單有效地確定單個樣本加入樣本集後對訓練結果的影響,一方面,出現新的樣本時可以利用原來的一訓練結果而不必重新開始;另一方面,讓訓練樣本逐個進入樣本集可以簡化尋優過程,提高演算法速度。這實際上是將樣本集中的樣本數減少到一個。
對於逐次超鬆弛迭代法,鬆弛因子的選取對演算法的收斂速度有很大影響,通常對於方程組Ax=Y,若A為正定矩陣,則當0<w <2時,逐次超鬆弛迭代恆收斂。
❼ 什麼是SOR
SOR是英文Standard Operation Records的縮寫,它的中文含意是記錄、憑證類文件。
SOR包括了:⑴物料管理記錄;⑵生產操作記錄;⑶質量檢驗與管理記錄;⑷設備運行與管理記錄;⑸衛生操作與管理記錄;⑹銷售記錄;⑺驗證報告與驗證記錄;⑻人員培訓與考核記錄;⑼文件發放。
SOR 基本解釋:
Schele Outlook Report 計劃遠景報告;sorbose 己酮糖;sale or return 銷售或利潤;
SOR 網路解釋:
1. 逐次超鬆弛:選主元素消去法的應用四 矩陣的三角分解五 平方根法及改進的平方根法六 追趕法七 列主元三角色解法 2 線性方程組的迭代解法一 雅可比(Jacobi)迭代法二 高斯-塞德爾(Gauss-seidel)迭代法三 逐次超鬆弛(SOR)迭代法 3 迭代法。
2. 女法師:最佳答案: PET的搭配要看你連的人物而選擇用什麼,其實在乎的主要就是他的光環,在大多數情況下,一般女法師(sor)會選擇噩夢 ACT 2(以後說的都是ACT 2的PET)里的防禦型,他的光環可以凍結敵人,降低他們的移動和攻擊速度。
3.sor:successive over-relaxation; 超鬆弛迭代法。
4.sor:subsystem – out – of – service – request; 子系統中斷服務請求。
5.sor:sales order; 訂單。
6.sor:synchrotron orbital radiation; 同步輻射。
拓展資料:
1.SOR 雙語例句Some of t hese food companies are a t t he forefron t of bes t prac t ice in t he food ins t ry - and t hey can't afford t o have t heir brand t ain t ed by con t amina t ion of anysort - including lubrican t s.
這其中的一些公司都處於食品行業執行最佳規范的前沿-因而,它們不希望它們的品牌因為某種形式的污染事件而受到玷污-包括潤滑油污染。
2. Compared with the result of theSORmethod, the result of this algorithm shows the algorithm can accelerate converg ence.
比較了該演算法和用超鬆弛法計算的結果,表明採用該演算法確實可以加快收斂速度。
3. The numerical results also express that the stability of new algorithm is better than that of theSORmethod of Cottle and Coheen。
數值結果還表明,新演算法比Cottle和Coheen的SOR法穩定。
❽ 急需用SOR法求解方程組的程序要求精度為10的-6次方
SOR迭代法的Matlab程序
function [x]=SOR_iterative(A,b)
% 用SOR迭代求解線性方程組,矩陣A是方陣
x0=zeros(1,length(b)); % 賦初值
tol=10^(-2); % 給定誤差界
N=1000; % 給定最大迭代次數
[n,n]=size(A); % 確定矩陣A的階
w=1; % 給定鬆弛因子
k=1;
% 迭代過程
while k<=N
x(1)=(b(1)-A(1,2:n)*x0(2:n)')/A(1,1);
for i=2:n
x(i)=(1-w)*x0(i)+w*(b(i)-A(i,1:i-1)*x(1:i-1)'-A(i,i+1:n)*x0(i+1:n)')/A(i,i);
end
if max(abs(x-x0))<=tol
fid = fopen('SOR_iter_result.txt', 'wt');
fprintf(fid,'\n********用SOR迭代求解線性方程組的輸出結果********\n\n');
fprintf(fid,'迭代次數: %d次\n\n',k);
fprintf(fid,'x的值\n\n');
fprintf(fid, '%12.8f \n', x);
break;
end
k=k+1;
x0=x;
end
if k==N+1
fid = fopen('SOR_iter_result.txt', 'wt');
fprintf(fid,'\n********用SOR迭代求解線性方程組的輸出結果********\n\n');
fprintf(fid,'迭代次數: %d次\n\n',k);
fprintf(fid,'超過最大迭代次數,求解失敗!');
fclose(fid);
end
❾ SOR迭代 程序
用戶信息
您好:金305
試用期 一級 ( 4 )
• 我的提問 • 個人中心
• 我的回答 • 退出
同類熱點問題
• 過程,模塊,函數,類,他們有什麼不同?
• 誰了解字體數據結構?
• 我用CAD2004作圖,在列印時為什麼有些字體打不出來
• 二叉樹(中序,前序,順序,靜態,傳地址)
• 不做題得高分(進來看看吧)!
更多>>
您想在自己的網站上展示
網路「知道」上的問答嗎?
來獲取免費代碼吧!
--------------------------------------------------------------------------------
如對「知道」有疑問或建議,請點此與我們聯系
訂閱該問題
待解決