MATLAB code for finding certain coefficients

Question

Sehra Sahu 2020-10-6

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/606486-matlab-code-for-finding-certain-coefficients

编辑： Bruno Luong 2020-10-11

采纳的回答： Bruno Luong

I am trying to input

In the summation,

s are distinct odd primes such that

How can it be written in MATLAB? Here

values can be 0 also. All the other symbols involved in the expression are natural numbers and r can be any suitable value. I couldn't get any idea.

9 个评论
显示 7更早的评论隐藏 7更早的评论

Bruno Luong 2020-10-10

One more question, does the order of (l1, l2, ... lr) matter? Meaning do you take all the permutations of {li} in the sum or they are selected said once, in the increasing order?

Sehra Sahu 2020-10-10

Order is not the matter. Just to remove the possible repetitions, we can consider them in increasing or decreasing order.

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Bruno Luong 2020-10-10

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/606486-matlab-code-for-finding-certain-coefficients#answer_510191

编辑：Bruno Luong 2020-10-11

在 MATLAB Online 中打开

Some detail of the formula is not totally clear to me (do we consider li>=0 or li>=1), but it goes like this

(EDIT, assume now the constraints li>=1)

m = 20;
r = 2;
verbose = true;
l = Partition(m, r);
% list all the odd primes
ptab = arrayfun(@(x) oddprime(x), 1:m, 'unif', 0);
% do the sum
s = 0;
for i=1:size(l,1)
    [p, li] = alldistinctp(l(i,:), ptab);
    if ~isempty(p) && size(p,2)==r
         q = li./p;
         if verbose
            T = array2table([li; p; q]);
            T.Properties.RowNames = {'l' 'p' 'l/p'};
            disp(T)
         end
         s = s + 1/prod(factorial(q));
    end
end
fprintf('s(m=%d,r=%d) = %f\n', m, r, s)
function [p, l] = alldistinctp(l, ptab)
l(l==0) = [];
n = size(l,2);
p = cell(1,n);
[p{:}] = ndgrid(ptab{l});
p = cat(n+1,p{:});
p = reshape(p, [], n);
% discard n-tuplets with repeated p
b = any(diff(sort(p,2),1,2)==0,2);
p(b,:) = [];
end
function p = oddprime(x)
p = unique(factor(x));
p(p==2) = [];
end
function v = Partition(n, lgt)
% v = Partition(n)
% INPUT
%   n: non negative integer
%   lgt: optional non negative integer
% OUTPUT:
%   v: (m x lgt) non-negative integer array such as sum(v,2)==n
%       each row of v is descending sorted
%       v contains all possible combinations
%       m = p(n) in case lgt == n, where p is the partition function
%       v is (dictionnary) sorted
% Algorithm:
%    Recursive
% Example:
% >> Partition(5)
% 
% ans =
% 
%      5     0     0     0     0
%      4     1     0     0     0
%      3     2     0     0     0
%      3     1     1     0     0
%      2     2     1     0     0
%      2     1     1     1     0
%      1     1     1     1     1
if nargin < 2
    lgt = n;
end
v = PartR(lgt+1, n, Inf);
end % Partition
%% Recursive engine of integer partition
function v = PartR(n, L1, head)
rcall = isfinite(head);
if rcall
    L1 = L1-head;
end
if n <= 2
    if ~rcall
        v = L1;
    elseif head >= L1
        v = [head, L1];
    else
        v = zeros(0, n, class(L1));
    end
else % recursive call
    j = min(L1,head):-1:0;
    v = arrayfun(@(j) PartR(n-1, L1, j), j(:), 'UniformOutput', false);
    v = cat(1,v{:});
    if rcall
        v = [head+zeros([size(v,1),1], class(head)), v];
    end
end
end % PartR

2 个评论
显示无隐藏无

Sehra Sahu 2020-10-11

You may consider

for time being.

Apparently the code shows some error in line 19, says function definition not permitted.

Bruno Luong 2020-10-11

编辑：Bruno Luong 2020-10-11

在 MATLAB Online 中打开

Easy to fix: If you run on older version MATLAB, save the three functions in separate mfiles from the script.

It works for me it give out this for m=20, r=2

           Var1    Var2
           ____    ____
    l       19      1  
    p       19      1  
    l/p      1      1  
           Var1    Var2
           ____    ____
    l       17      3  
    p       17      3  
    l/p      1      1  
           Var1    Var2
           ____    ____
    l       15      5  
    p        3      5  
    l/p      5      1  
           Var1    Var2
           ____    ____
    l       14      6  
    p        7      3  
    l/p      2      2  
           Var1    Var2
           ____    ____
    l       13      7  
    p       13      7  
    l/p      1      1  
           Var1    Var2
           ____    ____
    l       11      9  
    p       11      3  
    l/p      1      3  
s(m=20,r=2) = 3.425000