一、模拟大数乘法
输入是用字符串表示的乘数与被乘数,输出也需要转化为字符串。
思路是用数组的每个元素表示每一位的数,乘数与被乘数每一位都相乘,对应值再相加,最后for循环处理进位。
function re = simuMultiply(a,b)
A = a - '0';
B = b - '0';
D = arrayfun(@(k) sum(diag(rot90(A'*B),k)),(1-numel(b)):(numel(a)-1));
for ii = length(D):-1:2
D(ii-1) = D(ii-1) + floor(D(ii)/10);
D(ii) = mod(D(ii),10);
end
re = [int2str(D(1)) cast(D(2:end)+'0','char')];
end
二、卷积算法
卷积的概念这里就不解释了,每一位对应相乘后对应相加的意义恰好与卷积相通。
function re = convMultiply(a,b)
D = conv(a - '0',b-'0');
for ii = length(D):-1:2
D(ii-1) = D(ii-1) + floor(D(ii)/10);
D(ii) = mod(D(ii),10);
end
re = [int2str(D(1)) cast(D(2:end)+'0','char')];
end
三、FFT快速算法
根据频域相乘等于时域卷积的原理,可以用FFT算法实现卷积。事实上,MATLAB中的卷积就是用FFT实现的。
function re = fftMultiply(a,b)
A = a - '0';
B = b - '0';
D = ifft(fft([zeros(size(B)) A]).*fft([zeros(size(A)) B]));
for ii = length(D):-1:2
D(ii-1) = D(ii-1) + floor(D(ii)/10);
D(ii) = mod(D(ii),10);
end
re = [int2str(D(1)) cast(D(2:end-1)+'0','char')];
end
四、速度比较
a = randi([0,9],[1,1e5]);
a = num2str(a,'%d');
b = randi([0,9],[1,1e5]);
b = num2str(b,'%d');
%% 1000位×1000位乘法
tic;simuMultiply(a(1:1000),b(1:1000));toc
% 时间已过 7.239500 秒。
tic;convMulitply(a(1:1000),b(1:1000));toc
% 时间已过 0.002231 秒。
tic;fftMultiply(a(1:1000),b(1:1000));toc
% 时间已过 0.002402 秒。
%% 10000位×10000位乘法
tic;simuMultiply(a,b);toc
% Out of Memory,内存不足,无法计算
tic;convMulitply(a,b);toc
% 时间已过 0.260405 秒。
tic;fftMultiply(a,b);toc
% 时间已过 0.053287 秒。
总结:模拟计算的方法中间结果太多,内存占用大。卷积与FFT的方法占用内存小,速度相近,但总的来说,直接使用FFT更快。
易夕:MATLAB Tricks 专栏目录zhuanlan.zhihu.com