计算方法大作业Word格式.docx
《计算方法大作业Word格式.docx》由会员分享，可在线阅读，更多相关《计算方法大作业Word格式.docx（24页珍藏版）》请在冰豆网上搜索。
clearall;
clc;
A=[1111111;
2111111;
3211111;
4321111;
5432111;
6543211;
7654321];
b=[7;
8;
10;
13;
17;
22;
28];
n=length（A）;
IP=eye（n）;
U=zeros（n）;
L=eye（n）;
[m,p]=max（A（:
1））;
C1=b
（1）;
b
（1）=b（p）;
b（p）=C1;
C2（1:
n）=IP（1,:
）;
IP（1,:
）=IP（p,:
IP（p,:
）=C2（1:
n）;
n）=A（1,:
A（1,:
）=A（p,:
A（p,:
U（1,:
）=A（1,:
fork=2:
n
L（k:
n,k-1）=A（k:
n,k-1）/U（k-1,k-1）;
j=1;
fori=k:
S（j）=A（i,k）-L（i,1:
k-1）*U（1:
k-1,k）;
j=j+1;
[m,p]=max（abs（S（:
）））;
C2（1:
n）=A（k,1:
A（k,1:
n）=A（k+p-1,1:
A（k+p-1,1:
7）=C2（1:
n）=IP（k,1:
IP（k,1:
n）=IP（k+p-1,1:
IP（k+p-1,1:
n）=C2（1:
U（k,k:
n）=A（k,k:
n）-L（k,1:
k-1,k:
ifk<
A（k+1:
n,k）=A（k+1:
n,k）-L（k+1:
n,1:
end
disp（'
L='
disp（L）;
U='
disp（U）;
IP='
disp（IP）;
y=L\b;
x=U\y;
y='
disp（y）;
x='
disp（x）;
4、运行结果：
列主元三角分解原理很简单，用MATLAB实现起来却不容易。
对MATLAB还需要深入学习。
实验三SOR迭代法
一、实验内容：
采用SOR方法求解方程
采用0向量为初始向量，迭代以
结束，但迭代次数不操作1000次，松弛因子使用
其中
，给出使用的各松弛因子及对应迭代次数。
二、实验实现：
1、程序代码如下：
A=[10-10;
-110-2;
0-210];
b=[9;
7;
6];
x0=[0;
0;
0];
e=10^（-6）;
n=1000;
x=x0;
j=18;
w=0.1*j;
fork=1:
fori=1:
length（b）
x（i）=w*（b（i）-A（i,:
）*x）/A（i,i）+x（i）;
ifnorm（x-x0）<
e
disp（w）;
disp（x'
disp（k）;
------------'
）
return
x0=x;
2、运行结果如下：
上图显示为松弛因子分别为0.5、1、1.5和1.8时的结果。
实验四多项式插值
下列数据点的插值
x
1
4
9
16
25
36
49
64
y
2
3
5
6
7
8
可以得到平方根函数的近似，在区间[0,64]上作图.（注：
画图时以步长0.1计算出各点插值，再画出641个点的图，若能给出具体解析式，直接画出连续图）
（1）用这九个点作8次多项式插值L8（x）.，作图，从得到的结果看在[0,64]上，哪个插值区间更精确？
（2）使用三次插值（不是三弯矩，每个点采用距离最近的四个结点进行三次插值，如插值点为10，则取1,4,9,16四点做三次插值），作图，从得到的结果看在[0,64]上，哪个插值区间更精确？
（3）比较
（1）
（2）的结果。
x=[01491625364964];
y=[012345678];
x1=0:
0.1:
64;
y1=zeros（641,1）;
y2=zeros（641,1）;
forj=1:
1:
641
t=1;
ifi~=k
t=t*（x1（j）-x（i））/（x（k）-x（i））;
y1（j）=y1（j）+t*y（k）;
x2=0:
41
t=1
t=t*（x2（j）-x（i））/（x（k）-x（i））;
y2（j）=y2（j）+t*y（k）;
forj=42:
91
fork=2:
fori=2:
forj=92:
161
fork=3:
fori=3:
forj=162:
251
fork=4:
fori=4:
forj=252:
361
fork=5:
fori=5:
forj=362:
fork=6:
fori=6:
subplot（2,1,1）
plot（x,y,'
ro'
x1,y1,'
b.'
x1,sqrt（x1）,'
r'
gridon;
axis（[070010]）;
subplot（2,1,2）
x2,y2,'
k.'
3、结果分析
（1）由九个点做八次插值，由上图可知，仅在区间[016]上较精确，其他区间均与真实值有较大差距，区间[2564]上，甚至出现严重波动。
（2）做分段三次插值，发现在各个区间上，都能获得比较精确的插值。
仅在区间[04]上有少许偏差。
（3）两个结果对比说明，由于龙格现象的存在，高次插值虽能在更多点上与被插值函数的函数值相同，但从整体上看不一定能取得很好的逼近效果，甚至出现震荡。
这时选用分段低次插值可以取得较好的效果。
实验五数据拟合
1、浓度变化规律
在化学反应中，为研究某化合物的浓度随时间的变化规律，测得一组数据如下：
t（分）
6.4
8.0
8.4
9.28
9.5
9.7
9.86
9
10
11
12
13
14
15
10.2
10.32
10.42
10.5
10.55
10.58
10.6
计算出y与t的二次、三次多项式拟合函数，并推断20，40分钟的浓度值，
（1）程序代码如下：
t=1:
16;
y=[46.48.08.49.289.59.79.861010.210.3210.4210.510.5510.5810.6];
x=1:
50;
y2=polyfit（t,y,2）;
y3=polyfit（t,y,3）;
y_2=poly2str（y2,'
t'
y_3=poly2str（y3,'
y22=polyval（y2,x）;
y33=polyval（y3,x）;
y2_20=polyval（y2,20）
y2_40=polyval（y2,40）
y3_20=polyval（y3,20）
y3_40=polyval（y3,40）
plot（x,y22,'
x,y33,'
b'
t,y,'
ko'
holdon
plot（20,y2_20,'
b*'
40,y2_40,'
holdon
plot（20,y3_20,'
40,y3_40,'
gtext（'
二次拟合'
三次拟合'
（2）运行结果如下：
二次、三次拟合多项式为：
二次、三次拟合情况下，20分钟和40分钟的浓度值为：
二次、三次拟合的拟合曲线为：
2、经济增长模型（注：
可转变为矛盾方程组求解）
增加生产、发展经济所以靠的主要因素有增加投资、增加劳动力以及技术革新等，在研究国民经济产值与这些因素的数量关系时，由于技术水平不像资金、劳动力那样容易定量化，作为初步的模型，可认为技术水平不变，只讨论产值和资金、劳动力之间的关系。
在科学技术发展不快时，如资本主义经济发展的前期，这种模型是有意义的。
用Q，K，L分别表示产值、资金、劳动力，要寻求的数量关系Q（K,L）。
经过简化假设与分析，在经济学中，推导出一个著名的Cobb-Douglas生产函数
（*）
式中
要经由经济统计α数据确定。
现有美国马萨诸塞州1900-1926年上述三个经济指数的统计数据，见下表，试用数据拟合的方法，求出（*）式中的参数
。
马萨诸塞1900-1926年统计数据
t
Q
K
L
1900
1.05
1.04
1901
1.18
1.06
1.08
1902
1.29
1.16
1903
1.3
1.22
1904
1.27
1.17
1905
1.42
1.37
1906
1.5
1.44
1.39
1907
1.52
1.53
1.47
1908
1.46
1.57
1.31
1909
1.6
2.05
1.43
1910
1.69
2.51
1.58
1911
1.81
2.63
1.59
1912
1.93
2.74
1.66
1913
1.95
2.82
1.68
1914
2.01
3.24
1.65
1915
1.62
1916
2.09
3.61
1.86
1917
1.96
4.1
1918
2.2
4.36
1919
2.12
4.77
1920
2.16
4.75
1.9
1921
2.08
4.54
1922
2.24
1.67
1923
2.56
4.58
1.82
1924
2.34
1925
2.45
1.61
1926
2.58
1.64
（1）程序代码
对（*）式两边取对数并整理得：
令A=[1lnKlnL],
则lnQ=Ax,x=A\lnQ
K=[1.041.061.161.221.271.371.441.531.572.052.512.632.742.823.243.243.614.14.364.774.754.544.544.584.584.584.54]'
;
L=[1.051.081.181.221.171.31.391.471.311.431.581.591.661.681.651.621.861.931.961.951.91.581.671.821.61.611.64]'
Q=[1.051.181.291.31.31.421.51.521.461.61.691.811.931.952.0122.091.962.22.122.162.082.242.562.342.452.58]'
lnK=log（K）;
lnL=log（L）;
lnQ=log（Q）;
I=ones（size（Q））;
A=[IlnKlnL];
x=A\lnQ
a=exp（x
（1））
alpha=x
（2）
beta=x（3）
实验六数值积分
一、实验目的
掌握复化梯形公式和复化Simpson公式，会编写程序，上机进行实例计算。
二、实验内容
计算
（1）编写复化梯形公式和复化Simpson公式通用子程序，，分别采用4,8,16,32,64等分区间计算。
程序代码如下：
formatlong
f=inline（'
log（1+x）./x'
a=10^-8;
b=1;
fori=2:
n=2^i;
h=（b-a）/n;
x=zeros（n-1）;
sum=0;
fori=1:
n-1
x（i）=a+i*h;
sum=sum+f（x（i））;
fa=f（a）;
fb=f（b）;
T=h/2*（fa+2*sum+fb）;
n;
T
运行结果为：
（2）使用Romberg求积公式。
M=1;
h=b-a;
J=0;
err=1;
R=zeros（4,4）;
eps=10^-8;
R（1,1）=h*（f（a）+f（b））/2;
while（err>
eps）
J=J+1;
h=h/2;
x=a+h:
2*h:
b-h;
R（J+1,1）=R（J,1）/2+h*sum（f（x））;
forK=1:
min（3,J）
R（J+1,K+1）=R（J+1,K）+（R（J+1,K）-R（J,K））/（4^K-1）;
if（J>
3）
err=abs（R（J+1,4）-R（J,4））;
S=R（J,4）
实验七数值积分
实验内容
利用改进Euler法，经典四级四阶Rung-Kutta方法求其数值解，：
分别取h=0.2,0.1,0.05时数值解。
（注：
初值问题的精确解
）
1、改进Euler法
x0=0;
y0=1;
xn=2;
4*x./y-x.*y'
H=[0.20.10.05];
forj=1:
h=H（j）;
n=（xn-x0）/h;
x=x0;
y=y0;
y1=y+h*f（x,y）;
y2=y+h*f（x+h,y1）;
y=（y1+y2）/2;
x=x+h;
Y（j,i）=y;
h
Y（j,1:
n）
disp（'
_____________'
2、经典R-K方法
fori=1:
K1=h*f（x,y）;
K2=h*f（x+h/2,y+K1/2）;
K3=h*f（x+h/2,y+K2/2）;
K4=h*f（x+h,y+K3）;
y=y+（K1+2*K2+2*K3+K4）/6;
x=x+h;
升级会员 