数值计算——线性方程组的几个例题解决

1、

计算计算向量b的1范数、2范数、无穷范数。

计算A矩阵的1范数、2范数、无穷范数、条件数、F-范数。

clc,clear
 
b = [2
    8
    10];
A = [2 4 -2
    4 9 -3
    -2 -1 7];
 
fprintf('%s\n','向量b的1范数')
fanshub1 = abs(b(1)) + abs(b(2)) + abs(b(3))
fprintf('%s\n','向量b的2范数')
fanshub2 = sqrt((b(1)).^2 + (b(2)).^2 + (b(3)).^2)
fprintf('%s\n','向量b的无穷范数')
fanshub3 = max(abs(b))
 
fprintf('%s\n','矩阵A的1范数')
a = zeros(1,3);
sum = 0;
for i = 1:3
   for j = 1:3
       sum = sum + abs(A(j,i));
   end
   a(i) = sum;
   sum = 0;
end
fanshua1 = max(a)
 
fprintf('%s\n','矩阵A的2范数')
B = A'*A;
[v d] = eig(B);
fanshua2 = sqrt(max(max(d)))
 
fprintf('%s\n','矩阵A的无穷范数')
for i = 1:3
   for j = 1:3
       sum = sum + abs(A(i,j));
   end
   a(i) = sum;
   sum = 0;
end
fanshua3 = max(a)
 
fprintf('%s\n',' 1范数下矩阵A的条件数')
A_ = inv(A);
for i = 1:3
   for j = 1:3
       sum = sum + abs(A_(j,i));
   end
   a(i) = sum;
   sum = 0;
end
conda1 = fanshua1*max(a)
 
fprintf('%s\n',' 2范数下矩阵A的2条件数')
B = A_'*A_;
[v d] = eig(B);
conda2 = sqrt(max(max(d)))*fanshua2
 
fprintf('%s\n','无穷范数下矩阵A的条件数')
for i = 1:3
   for j = 1:3
       sum = sum + abs(A_(i,j));
   end
   a(i) = sum;
   sum = 0;
end
conda3 = max(a)*fanshua3
 
fprintf('%s\n','矩阵A的F范数')
for i = 1:3
   for j = 1:3
       sum = sum + (A(i,j)).^2;
   end
end
fanshuaf = sqrt(sum)
sum = 0;

运行结果如下所示

2、楚列斯基分解

（1）、写出Cholesky分解函数。

用MATLAB代码实现如下所示

Cholesky.m
function [ B ] = Cholesky(A )
%UNTITLED2 此处显示有关此函数的摘要
%   楚列斯基分解
B = A;
n = length(A);
for k = 1:n
   B(k,k) = sqrt(B(k,k));
   for i = k+1:n
      B(i,k) = B(i,k)/B(k,k); 
   end
   for j = k+1:n
      for i = k+1:n
          B(i,j) = B(i,j)-B(i,k)*B(j,k);
      end
   end
end

for i = 1:n-1
    for j = i+1:n
        B(i,j) = 0;
    end
end

end

（2）、

解方程组Hx=b（取n=2，5，10）

用MATLAB编写代码如下

Fiteration.m
function [ x ] = Fiteration( L,b )
%UNTITLED8 此处显示有关此函数的摘要
%   此处显示详细说明
[n1,n]=size(L);
x = zeros(n,1);
if n1==n
    j = 1;
    while j<n+1 && L(j,j)~=0
        x(j,1) = b(j)/L(j,j);
        for i = j+1:n
           b(i) = b(i)-L(i,j)*x(j,1); 
        end
        j=j+1;
    end
end

end

Biteration.m
function [ x ] = Biteration( U,b )
%UNTITLED9 此处显示有关此函数的摘要
%   此处显示详细说明
[n1,n] = size(U);
x = zeros(n,1);
if n1==n
    j = n;
    while j>0 && U(j,j)~=0
        x(j,1) = b(j)/U(j,j);
        for i = 1:j-1
           b(i) = b(i)-U(i,j)*x(j); 
        end
        j = j-1;
    end
    
end

end

test221.m
clc,clear
fprintf('%s\n','n = 2时')
A = [1 1/2
    1/2 1/3];
L = Cholesky(A)
x1 = [1
    1];
b = A*x1
front = Fiteration( L,b );
x = Biteration(L',front)

fprintf('%s\n','n = 5时')
A = zeros(5,5);
for i = 1:5
   for j = 1:5
      A(i,j) = 1/(i+j-1); 
   end
end
L = Cholesky(A)
x1 = ones(5,1);
b = A*x1
front = Fiteration( L,b );
x = Biteration(L',front)

fprintf('%s\n','n = 10时')
A = zeros(10,10);
for i = 1:10
   for j = 1:10
      A(i,j) = 1/(i+j-1); 
   end
end
L = Cholesky(A)
x1 = ones(10,1);
b = A*x1
front = Fiteration( L,b );
x = Biteration(L',front)

运行结果如下所示

（3）、计算H矩阵的1范数、2范数、无穷范数、条件数（取n=2，5，10）

用MATLAB代码实现如下

f222.m
function [  ] = f222(n )
%UNTITLED5 此处显示有关此函数的摘要
%   此处显示详细说明

A = zeros(n,n);
for i = 1:n
   for j = 1:n
      A(i,j) = 1/(i+j-1); 
   end
end

fprintf('%s\n','矩阵的1范数为')
a = zeros(1,n);
sum = 0;
for i = 1:n
   for j = 1:n
       sum = sum + abs(A(j,i));
   end
   a(i) = sum;
   sum = 0;
end
fanshua1 = max(a)

fprintf('%s\n','矩阵A的2范数')
B = A'*A;
[v d] = eig(B);
fanshua2 = sqrt(max(max(d)))

fprintf('%s\n','矩阵A的无穷范数')
for i = 1:n
   for j = 1:n
       sum = sum + abs(A(i,j));
   end
   a(i) = sum;
   sum = 0;
end
fanshua3 = max(a)

fprintf('%s\n','1范数下矩阵A的条件数')
A_ = inv(A);
for i = 1:n
   for j = 1:n
       sum = sum + abs(A_(j,i));
   end
   a(i) = sum;
   sum = 0;
end
conda1 = fanshua1*max(a)

fprintf('%s\n','2范数下矩阵A的2条件数')
B = A_'*A_;
[v d] = eig(B);
conda2 = sqrt(max(max(d)))*fanshua2

fprintf('%s\n','无穷范数下矩阵A的条件数')
for i = 1:n
   for j = 1:n
       sum = sum + abs(A_(i,j));
   end
   a(i) = sum;
   sum = 0;
end
conda3 = max(a)*fanshua3

end
test222.m
clc,clear
fprintf('%s\n','n = 2时')
f222(2)
fprintf('%s\n','n = 5时')
f222(5)
fprintf('%s\n','n = 10时')
f222(10)