MATLAB code for finding certain coefficients

4 views (last 30 days)
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 Comments
Bruno Luong
Bruno Luong on 10 Oct 2020
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
Sehra Sahu on 10 Oct 2020
Order is not the matter. Just to remove the possible repetitions, we can consider them in increasing or decreasing order.

Sign in to comment.

Accepted Answer

Bruno Luong
Bruno Luong on 10 Oct 2020
Edited: Bruno Luong on 11 Oct 2020
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 Comments
Sehra Sahu
Sehra Sahu on 11 Oct 2020
You may consider for time being.
Apparently the code shows some error in line 19, says function definition not permitted.
Bruno Luong
Bruno Luong on 11 Oct 2020
Edited: Bruno Luong on 11 Oct 2020
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

Sign in to comment.

More Answers (0)

Tags

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!