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 里面矩阵如何既有字符变量又有数值变量