当前位置:编程学习 > Matlab >>

matlab 麻烦帮我在看看这个问题怎么解决吧,谢谢。。。。

我这有个高斯消元法的程序,怎么能把这个程序插入我要运算的程序中并且能运算成功 function x=gaussMethod(A,b) %高斯列主元消去法,要求系数矩阵非奇异的, % n = size(A,1); if abs(det(A))<= 1e-8 error('系数矩阵是奇异的'); return; end % for k=1:n ak = max(abs(A(k:n,k))); index = find(A(:,k)==ak); if length(index) == 0 index = find(A(:,k)==-ak); end %交换列主元 temp = A(index,:); A(index,:) = A(k,:); A(k,:) = temp; temp = b(index);b(index) = b(k); b(k) = temp; %消元过程 for i=k+1:n m=A(i,k)/A(k,k); %消除列元素 A(i,k+1:n)=A(i,k+1:n)-m*A(k,k+1:n); b(i)=b(i)-m*b(k); end end %回代过程 x(n)=b(n)/A(n,n); for k=n-1:-1:1; x(k)=(b(k)-A(k,k+1:n)*x(k+1:n)')/A(k,k); end x=x'; end 下面是我的程序 syms Q1 Q2 Q3 n a P P1 P2 P3 R L1 L2 L3 B k h u pi; F1=(Q1+Q2+Q3)*u/(2*pi*n*k*h)*log(a/((pi*R)))+(Q1+Q2+Q3)*u/(B*k*h)*L1+Q1*u/(B*k*h)*L2-P+P1; F2=-Q1*u/(2*pi*n*k*h)*log(a/((pi*R)))+(Q2+Q3)*u/(B*k*h)*L2+Q2*u/(2*pi*n*k*h)*log(a/((pi*R)))-P1+P2; F3=-Q2*u/(2*pi*n*k*h)*log(a/((pi*R)))+Q3*u/(2*pi*n*k*h)*log(a/((pi*R)))+Q3*u/(B*k*h)*L3-P2+P3; n=16; a=450; P=19.5; P1=7.5;P2=7.5;P3=7.5; R=10; L1=1100; L2=600; L3=600; B=1.12; k=0.5; h=16; u=9;pi=3.1415; F=eval([F1 F2 F3]); F1=F(1);F2=F(2);F3=F(3); for i=1:3 for j=1:3 A(i,j)=eval(['diff(F' num2str(i) ',''Q' num2str(j) ''')']); end end b=[P-P1; P1-P2; P2-P3]; Q=inv(A)*b; >> Q1=Q(1);Q2=Q(2);Q3=Q(3); >> Q1=vpa(Q1),Q2=vpa(Q2),Q3=vpa(Q3) 老师不让求逆运算,唉。。。请高手指教了,谢谢
答案:syms Q1 Q2 Q3 n a P P1 P2 P3 R L1 L2 L3 B k h u pi;
F1=(Q1+Q2+Q3)*u/(2*pi*n*k*h)*log(a/((pi*R)))+(Q1+Q2+Q3)*u/(B*k*h)*L1+Q1*u/(B*k*h)*L2-P+P1;
F2=-Q1*u/(2*pi*n*k*h)*log(a/((pi*R)))+(Q2+Q3)*u/(B*k*h)*L2+Q2*u/(2*pi*n*k*h)*log(a/((pi*R)))-P1+P2;
F3=-Q2*u/(2*pi*n*k*h)*log(a/((pi*R)))+Q3*u/(2*pi*n*k*h)*log(a/((pi*R)))+Q3*u/(B*k*h)*L3-P2+P3;
n=16; a=450; P=19.5; P1=7.5;P2=7.5;P3=7.5; R=10; L1=1100; L2=600; L3=600; B=1.12; k=0.5; h=16; u=9;pi=3.1415;
F=eval([F1 F2 F3]);
F1=F(1);F2=F(2);F3=F(3);
for i=1:3
for j=1:3
A(i,j)=eval(['diff(F' num2str(i) ',''Q' num2str(j) ''')']);
end
end
b=[P-P1; P1-P2;  P2-P3];
A=eval(A)
%%%%%%%%%%%%%%%%%%%%%
n = size(A,1);
if abs(det(A))<= 1e-8
error('系数矩阵是奇异的');
return;
end
for k=1:n
ak = max(abs(A(k:n,k)));
index = find(A(:,k)==ak);
if length(index) == 0
index = find(A(:,k)==-ak);
end
temp = A(index,:);
A(index,:) = A(k,:);
A(k,:) = temp;
temp = b(index);b(index) = b(k); b(k) = temp;
for i=k+1:n
m=A(i,k)/A(k,k);
A(i,k+1:n)=A(i,k+1:n)-m*A(k,k+1:n);
b(i)=b(i)-m*b(k);
end
end
x(n)=b(n)/A(n,n);
for  k=n-1:-1:1;
x(k)=(b(k)-A(k,k+1:n)*x(k+1:n)')/A(k,k);
end
Q=x';
format long e
Q1=Q(1),Q2=Q(2),Q3=Q(3),
若满意记得采纳!^.^

上一个:关于matlab中的一点小问题,麻烦帮忙解决一下哈,谢谢:
下一个:matlab 里面矩阵如何既有字符变量又有数值变量

CopyRight © 2012 站长网 编程知识问答 www.zzzyk.com All Rights Reserved
部份技术文章来自网络,