标签为'时间换空间'的全部存档

时间换空间/空间换时间？

发表于： 2010年05月24日 分类：数字之美 0条评论 标签: , , , , , .

经常看到网上有人在讨论算法时用时间换空间/空间换时间的说法. 可惜神人们使用的例子都很专业, 再加上一两句晦涩而潇洒的写法, 我等凡人基本上就很难看懂了. 今天在实现数值积分中的龙贝格 (Romberg) 算法时, 算是遇到了一个这方面有点关系的问题, 简单记录一下.

首先上一段800字的原理说明:

……

由于作者心里很清楚读者将会人肉Skip上面一段, 所以这里直接人肉省略其中的797字.

如果执意想看一下, 参考 Wikipedia 好了.

下面是一个很好懂的 MATLAB 实现:

?Download romberg.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63

function [ ] = romberg(aa, bb, kk, iter)
% Romberg Integral Method
% Coded by Nan.Xiao 2010-05-23
% Step 1. Calc TrapeStep
% Step 2. Calc Romberg Method
% Step 3. Extract Diagnal & Output
function [ ] = trapestep(a, b, k)
    format long
%     function y = f(x)
%         y=(2/sqrt(pi))*exp(-x);
%     end
%     Examples Can Also Be:
%     function y = f(x)
%         y=x*sin(x);
%     end
%     Or
    function y = f(x)
        y=x*sqrt(1+x.^2);
    end
k1=2^(k-1);
temp=zeros(k,k1);
mysum=zeros(k,1);
t=zeros(k+1,1);
for i=1:k
    for j=1:(2^(i-1))
    temp(i,j)=f(a+(((2*j-1)*(b-a))/(2^i)));
    end
end
for p=1:k
    mysum(p,1)=sum(temp(p,:));
end
t(1,1)=0.5*(f(a)+f(b));
for q=2:k+1
    t(q,1)=(0.5*t(q-1,1))+(((b-a)*mysum(q-1,1))/(2^(q-1)));
end
end
    format long
trapestep(aa, bb, kk);
diagonal=zeros(kk,1);
tt=zeros(kk,kk);
for v=1:kk
    tt(v,1)=t(v,1);
end
for r=2:kk
for s=2:r
    tt(r,s)=((4^(s-1)*tt(r,s-1))-tt(r-1,s-1))/(4^(s-1)-1);
end
end
save('tt.txt', 'tt', '-ascii', '-double'); % Save the Table
% Extract Table Diagonals
for w=1:kk
    diagonal(w,1)=tt(w,w);
end
disp('Romberg Table Diagonals:')
disp(diagonal);
% Output Final Result
cc=2;
while abs(diagonal(cc,1)-diagonal(cc-1,1))>=iter || abs(diagonal(cc+1,1)-diagonal(cc,1))>=iter
    cc=cc+1;
end
disp('The Final Result is:');
disp(diagonal(cc,1));
end

与网上其他程序不同, 我在这里没有做太多的判断, 将二分的次数k作为一个参数手动设定, 是为了观察计算的具体过程. 同时, 为了方便移植, 仅使用最基本的矩阵操作语句, 没有使用 MATLAB 中的 feval() 函数, 要改变被积函数, 直接修改函数体. 最终目的是, 我要取得计算产生的含有T_m^{(k)}的表格, 以方便完成笔头的作业. 懒到一定程度了是不是?
继续阅读’时间换空间/空间换时间？’