MATLAB计算幂与阶乘比值时出现NaN问题的代码修改咨询

代码出现NaN的核心是双精度浮点数溢出：
当循环变量m增大到170以上时，分子(2*d)^m、分母factorial(m)的单独计算值都会超过MATLAB双精度浮点数的最大表示上限（约1.8e308），溢出后两项都会被判定为Inf，而MATLAB中Inf/Inf的运算结果就是NaN。一旦累加过程中出现NaN，最终求和结果就会全程保持NaN。
实际上两个大数的比值(2*d)^m/factorial(m)本身远小于双精度上限，只是单独计算分子、分母时提前溢出才触发异常。

最简便高效的修复方式是改用递推方式计算通项，避免单独计算大数幂、阶乘这类容易溢出的项。观察通项规律：第m项的a1值和第m-1项的比值为(2*d)/m，从m=0的初始值开始递推即可，全程不会出现单独的超大数溢出问题。
修改后的可运行代码：

l=0;
lp=2;
d=200;
a2= ((-sqrt(d))^lp);
a3 =((-sqrt(d))^l);
sum1 = 0;
% m=0时的初始a1值: (2^0*d^0)/factorial(0) = 1
a1 = 1;
sum1 = sum1 + a1*a2*a3;
% 从m=1开始递推计算
for m = 1:500
    a1 = a1 * (2*d)/m;
    mul = a1*a2*a3;
    sum1 = sum1 + mul;
end

补充说明 #

你当前代码实现的求和本质是指数函数的泰勒展开，当循环上限取500时，求和结果已经和解析解几乎完全一致，可以直接用解析解计算获得更高效率，不需要循环：

% 解析解直接计算结果
sum_res = a2*a3*exp(2*d);

如果遇到更复杂的、无法直接递推的通项求和场景，也可以先在对数域计算每一项的对数值，把幂、阶乘运算转换为对数域的加减、数乘运算，算完后再取指数还原，同样可以避免大数溢出问题。

内容的提问来源于stack exchange，提问作者Jasmine