Summation with upper limit itself as a variable

i m trying to optimize the following function

z=-(r)*(L/M)*(1-(((L/M).^x(1))/((gamma(x(1)+1))*symsum(((L/M).^p)/(gamma(p+1)),p,0,x(1)))))+c*x(1);

 with 
syms p
r=2;
L=10;
M=2;
c=1;

 and getting the following error

??? Error using ==> mupadmex
Error in MuPAD command: Illegal range [sum::sum]

can anybody help me how can i set the upper limit as variable in summation??

[Walter Roberson]

You have not shown how you define x ?

[Amir]

it is basically a part of algo where we want to find the best value of function and x(1). range of x(1) is defined within (1,100)

[Amir]

if i give some numeric value instead of x(1), then it runs. according to my observation it takes numeric value in this format. is there anyother way to define summation with upper limit as free variable??

[Walter Roberson]

Please show the statement you used to initialize x and to set the value of x(1)

[Amir]

well! em using the algo available at web

http://www.mathworks.com/matlabcentral/fileexchange/29693-firefly-algorithm/content/fa_ndim.m

and trying to implement it for my objective function.

whole format is as follows;

% ======================================================== % 
% Files of the Matlab programs included in the book:       %
% Xin-She Yang, Nature-Inspired Metaheuristic Algorithms,  %
% Second Edition, Luniver Press, (2010).   www.luniver.com %
% ======================================================== %    

% -------------------------------------------------------- %
% Firefly Algorithm for constrained optimization using     %
% for the design of a spring (benchmark)                   % 
% by Xin-She Yang (Cambridge University) Copyright @2009   %
% -------------------------------------------------------- %

function fa_ndim
% parameters [n N_iteration alpha betamin gamma]
para=[20 50 0.1 0.1 1];

help fa_ndim.m

% Simple bounds/limits for d-dimensional problems
d=1;
Lb=ones(1,d);
Ub=100*ones(1,d);

% Initial random guess
u0=Lb+(Ub-Lb).*rand(1,d);

[u,fval,NumEval]=ffa_mincon(@cost,u0,Lb,Ub,para);

% Display results
bestsolution=u
bestojb=fval
total_number_of_function_evaluations=NumEval

%%%Put your own cost/objective function here --------%%%
%%Cost or Objective function
 function z=cost(x)
% Exact solutions should be (1,1,...,1) 
syms p
r=2;
L=10;
M=2;
c=1;

z=-(r)*(L/M)*(1-(((L/M).^x(1))/((gamma(x(1)+1))*symsum((L/M).^p/(gamma(p+1)),p,0,x(1)))))+c*x(1);

%%%End of the part to be modified -------------------%%%

%%%--------------------------------------------------%%%
%%%Do not modify the following codes unless you want %%%
%%%to improve its performance etc                    %%%
% -------------------------------------------------------
% ===Start of the Firefly Algorithm Implementation ======
%         Lb = lower bounds/limits
%         Ub = upper bounds/limits
%   para == optional (to control the Firefly algorithm)
% Outputs: nbest   = the best solution found so far
%          fbest   = the best objective value
%      NumEval = number of evaluations: n*MaxGeneration
% Optional:
% The alpha can be reduced (as to reduce the randomness)
% ---------------------------------------------------------

% Start FA
function [nbest,fbest,NumEval]...
           =ffa_mincon(fhandle,u0, Lb, Ub, para)
% Check input parameters (otherwise set as default values)
if nargin<5, para=[20 50 0.1 0.1 1]; end
if nargin<4, Ub=[]; end
if nargin<3, Lb=[]; end
if nargin<2,
disp('Usuage: FA_mincon(@cost,u0,Lb,Ub,para)');
end

% n=number of fireflies
% MaxGeneration=number of pseudo time steps
% ------------------------------------------------
% alpha=0.25;      % Randomness 0--1 (highly random)
% betamn=0.20;     % minimum value of beta
% gamma=1;         % Absorption coefficient
% ------------------------------------------------
n=para(1);  MaxGeneration=para(2);
alpha=para(3); betamin=para(4); gamma=para(5);

% Total number of function evaluations
NumEval=n*MaxGeneration;

% Check if the upper bound & lower bound are the same size
if length(Lb) ~=length(Ub),
    disp('Simple bounds/limits are improper!');
    return
end

% Calcualte dimension
d=length(u0);

% Initial values of an array
zn=ones(n,1)*10^100;
% ------------------------------------------------
% generating the initial locations of n fireflies
[ns,Lightn]=init_ffa(n,d,Lb,Ub,u0);

% Iterations or pseudo time marching
for k=1:MaxGeneration,     %%%%%start iterations

% This line of reducing alpha is optional
 alpha=alpha_new(alpha,MaxGeneration);

% Evaluate new solutions (for all n fireflies)
for i=1:n,
   zn(i)=fhandle(ns(i,:));
   Lightn(i)=zn(i);
end

% Ranking fireflies by their light intensity/objectives
[Lightn,Index]=sort(zn);
ns_tmp=ns;
for i=1:n,
 ns(i,:)=ns_tmp(Index(i),:);
end

%%Find the current best
nso=ns; Lighto=Lightn;
nbest=ns(1,:); Lightbest=Lightn(1);

% For output only
fbest=Lightbest;

% Move all fireflies to the better locations
[ns]=ffa_move(n,d,ns,Lightn,nso,Lighto,nbest,...
      Lightbest,alpha,betamin,gamma,Lb,Ub);

end   %%%%%end of iterations

% -------------------------------------------------------
% ----- All the subfunctions are listed here ------------
% The initial locations of n fireflies
function [ns,Lightn]=init_ffa(n,d,Lb,Ub,u0)
  % if there are bounds/limits,
if length(Lb)>0,
   for i=1:n,
   ns(i,:)=Lb+(Ub-Lb).*rand(1,d);
   end
else
   % generate solutions around the random guess
   for i=1:n,
   ns(i,:)=u0+randn(1,d);
   end
end

% initial value before function evaluations
Lightn=ones(n,1)*10^100;

% Move all fireflies toward brighter ones
function [ns]=ffa_move(n,d,ns,Lightn,nso,Lighto,...
             nbest,Lightbest,alpha,betamin,gamma,Lb,Ub)
% Scaling of the system
scale=abs(Ub-Lb);

% Updating fireflies
for i=1:n,
% The attractiveness parameter beta=exp(-gamma*r)
   for j=1:n,
      r=sqrt(sum((ns(i,:)-ns(j,:)).^2));
      % Update moves
if Lightn(i)>Lighto(j), % Brighter and more attractive
   beta0=1; beta=(beta0-betamin)*exp(-gamma*r.^2)+betamin;
   tmpf=alpha.*(rand(1,d)-0.5).*scale;
   ns(i,:)=ns(i,:).*(1-beta)+nso(j,:).*beta+tmpf;
      end
   end % end for j

end % end for i

% Check if the updated solutions/locations are within limits
[ns]=findlimits(n,ns,Lb,Ub);

% This function is optional, as it is not in the original FA
% The idea to reduce randomness is to increase the convergence,
% however, if you reduce randomness too quickly, then premature
% convergence can occur. So use with care.
function alpha=alpha_new(alpha,NGen)
% alpha_n=alpha_0(1-delta)^NGen=10^(-4);
% alpha_0=0.9
delta=1-(10^(-4)/0.9)^(1/NGen);
alpha=(1-delta)*alpha;

% Make sure the fireflies are within the bounds/limits
function [ns]=findlimits(n,ns,Lb,Ub)
for i=1:n,
     % Apply the lower bound
  ns_tmp=ns(i,:);
  I=ns_tmp<Lb;
  ns_tmp(I)=Lb(I);

    % Apply the upper bounds
    J=ns_tmp>Ub;
    ns_tmp(J)=Ub(J);
    % Update this new move
    ns(i,:)=ns_tmp;
  end

%%==== End of Firefly Algorithm implementation ======

[Walter Roberson]

The code you show never has "x" as a symbolic variable, so your upper bound is always numeric. I need to see the version of the code that has the error.

Answer 1

[Roger Stafford]

In the expression

 symsum((L/M)^p/gamma(p+1),p,0,n)

for p an integer ranging from 0 to n it is the sum of the first n+1 terms in the infinite power series expansion for the exponential function exp(L/M). It is extremely doubtful that there exists any known "closed form" for this finite summation in terms of a general n. In integral calculus, there are some integrals for which there are known closed forms (which one can usually find listed in integral tables) but many more for which none are known. In an analogous way there are finite summations for which there are known closed forms such as

 symsum(p^2,p,1,n) = n*(n+1)*(2*n+1)/6

and many more for which none are known. As I say, it is very likely the summation you describe belongs to that latter category. One should therefore not blame Mupad for not finding it.

[Walter Roberson]

Correct, but MuPAD would not return an error about invalid range for that case. MuPAD would create an unevaluted sum() expression (which might get converted to symsum() at the MATLAB level)

Answer 2

[Amir]

i m not getting what is the issue? if i use the same expression without symsum(..), i'm being returned with the optimized values of function 'z' and decision variable 'x', which shows that code is working properly but with the induction of summation with upper bound as decision variable, it generates the error! is there any other way to induce the following summation here???

"∑ p=0 to x {(L/M)^p/p!}"

where xϵ{1,100}

[Roger Stafford]

The issue, Amir, is that (presumably) no-one has yet found out how to express your sum as a compact general formula in terms of an unlimited variable x except by using the type of "symsum" or "sum" notation you already have here, and therefore our computers are also unable to do so. With your limit of 100 it would be possible to do the equivalent of listing all one hundred possibilities, but that is surely not what you are looking for.

[Walter Roberson]

The error you reported earlier has to do with the bounds being invalid, but you have not yet shown us a code example in which the code using a variable as the upper bound.