如何在Matlab中计算双重和，其中第二个和的上限是第一个和的下限？

Question

我在 Matlab 中实现一个调用我在不同 .m 文件中编写的其他函数的函数时遇到很大困难。我被卡住的部分是对另一个函数中输入的不同值求和的部分，其中另一个函数内部也有一个和。第一个和的下限是第二个和的上限。

功能如下：

对于作为向量输入的 n 和作为标量输入的 x，我的 Hh(n,x) 函数可以正常工作。由于 n 的向量化输入，In 函数内部的 Hi 上的和可以通过调用 sum(Hh(0:n,x)) 快速计算出来。

我想对 In 函数做同样的事情，但是因为 In 内部的 n 现在的范围是 0 到 k-1，而外部函数中的 k 范围是 1 到 n，所以我不知道如何计算这个双和，其中内部和以外部和的下限作为上限。我想尽可能有效地评估这个双重和，因为以后我想用这些公式做很多模拟。现在我在 n 次评估函数，将每个值存储在一个向量中，然后求和，这在计算上非常密集......

我的 In 函数的 Matlab 代码是：

function in = In(n,c,alphaa,betaa, delta)
ie = 0:n;
in = -(exp(alphaa*c)/alphaa)...
     .*sum((betaa/alphaa).^(n-ie).*Hh(ie,betaa*c-delta))...
     -(betaa/alphaa).^(n+1)
end

外部函数的Matlab代码，我们暂时称它为函数f：

function f = F(n,a,mu,sigma,eta1,T)
for k = 1: n
    vector(k) =  In(k-1,a-mu*T,-eta1,-1/(sigma*sqrt(T)),-(sigma*eta1*sqrt(T)));
end
f = sum(vector);
end

如何让In里面的n的输入向量化，这样我就不用把所有输入的n值都单独存起来再计算总和，而是直接对输入的向量n计算总和。

感谢任何帮助，因为我现在严重卡住了！提前谢谢！

Answer 1

所以要具体回答你的问题（n 向量化的输入），我会执行以下操作：

% 1.为每个k创建一个从0到k-1的向量ivect，以及与迭代ivect对应的向量kvect：

ivect=[];
kvect=[];
for k=1:n
    ivect = [ivect 0:k-1];
    kvect = [kvect (k-1)*ones(1,k)];
end

% 2. 根据您的迭代 i 和 k 创建两个总和内的函数：

    function in = In2(k,ie,c,alphaa,betaa, delta)
in = -(exp(alphaa*c)/alphaa)*((betaa/alphaa).^(k-ie).*func(ie,betaa*c-delta))-(betaa/alphaa).^(k+1)./(k+1);
end

% 3. 然后求和：

f=sum(In2(kvect,ivect,a-mu*T,-eta1,-1/(sigma*sqrt(T)),-(sigma*eta1*sqrt(T))))

起初我认为这种方法会减少计算总和的计算时间（当然第1点是昂贵的，但可以一次性完成并用于多次模拟而无需重建这些向量），但即使只比较行：

tic
for k = 1: n
    vector(k) =  In(k-1,a-mu*T,-eta1,-1/(sigma*sqrt(T)),-(sigma*eta1*sqrt(T)));
end
f = sum(vector)
toc

与

tic
f=sum(In2(kvect,ivect,a-mu*T,-eta1,-1/(sigma*sqrt(T)),-(sigma*eta1*sqrt(T))))
toc

第一种方法得到1.029899秒，第二种方法得到1.684432秒，n = 5000。但幸运的是，我找到了原因！这是因为这种方法实际上产生了更多的人工操作，如在函数 In2 中，我要求我的代码为每个 k 计算 -(betaa/alphaa).^(k+1)./(k+1) 和i，而该术语仅取决于 k。所以现在，让我们定义函数 In3 和 In4：

function in = In3(k,ie,c,alphaa,betaa, delta)
in = -(exp(alphaa*c)/alphaa)*((betaa/alphaa).^(k-ie).*Hh(ie,betaa*c-delta));    end

function in = In4(k,alphaa,betaa)
    in=-(betaa/alphaa).^(k+1);
end

然后打电话

tic
f=sum(In3(kvect,ivect,a-mu*T,-eta1,-1/(sigma*sqrt(T)),-(sigma*eta1*sqrt(T)))) +sum(In4([1:n],-eta1,-1/(sigma*sqrt(T))))
toc

瞧，对于 N=5000，我在 0.929530 秒内获得了相同的结果。总而言之，使用最后一种方法可以节省一些时间，但公平地说，我会说你有两个缺点：

1. 一开始就 ivect 和 kvect 的时间计算——但我坚持认为这些向量可以一次性构建、存储并用于多次模拟——

2. 内存消耗：ivect 和 kvect 的大小为 n(n+1)/2，对于较大的 n 来说，这很多

希望它对您有所帮助，无论如何感谢您提出的有趣问题，这可能会在我自己的研究中为我服务！