Path: news.mathworks.com!not-for-mail
From: <HIDDEN>
Newsgroups: comp.soft-sys.matlab
Subject: Re: Integer solutions for a_1*x_1+a_2*x_2+a_3*x_3+ ...+a_4*x_n = k
Date: Tue, 18 Nov 2008 20:43:02 +0000 (UTC)
Organization: Xoran Technologies
Lines: 98
Message-ID: <gfv9cm$55g$1@fred.mathworks.com>
References: <gfun9a$41t$1@fred.mathworks.com>
Reply-To: <HIDDEN>
NNTP-Posting-Host: webapp-03-blr.mathworks.com
Content-Type: text/plain; charset="ISO-8859-1"
Content-Transfer-Encoding: 8bit
X-Trace: fred.mathworks.com 1227040982 5296 172.30.248.38 (18 Nov 2008 20:43:02 GMT)
X-Complaints-To: news@mathworks.com
NNTP-Posting-Date: Tue, 18 Nov 2008 20:43:02 +0000 (UTC)
X-Newsreader: MATLAB Central Newsreader 1440443
Xref: news.mathworks.com comp.soft-sys.matlab:501560


"Oriole " <oriole_ni@hotmail.com> wrote in message <gfun9a$41t$1@fred.mathworks.com>...
> Hello, 
> 
> I am wondering how to write a function for solving
> 
> a_1*x_1+a_2*x_2+a_3*x_3+ ...+a_4*x_n = k
> 
> where
> Inputs: 
> n, k, and vector of non-negaive integer coefficients a=[a_1, a_2, ...a_n].     
> k is some small positive integer, for example, k= 3, or k= 6, n is rather large positve integer, for example,, n = 20, or n= 60 or in some rare cases even n=100.
> 
> Output:
> The unknowns x_1, x_2, ...x_n are non-negative integers.

It wasn't clear to me whether you were looking for one solution or all of them.
If you require only one solution,  my dynamic programming idea is applicable and I've implemented it below. It is quite fast, even for n=100


>> aa=round(rand(1,100)*10+2);

>> tic;  xx=dioph(aa,50); toc  %test for n=100 and k=50
Elapsed time is 0.000800 seconds.


Of course, I don't know if my aa test vector is representative of your application...








###############################

function xx=dioph(aa,kk)
%Finds non-negative integer solution xx to the equation
%
%   aa'*xx=kk
%
% where aa is a strictly positive integer vector and kk is a non-negative
% scalar

if any(aa-round(aa))
    error 'Not integer valued a'
end

if any(aa<=0)
    error 'Should be strictly positive a'
end

nn=length(aa);
states=(-1:kk).';
Nstates=kk+2;
Nactions=kk+1;

Smap=repmat(states,[1,Nactions]);
Dmap=repmat((0:kk),[Nstates,1]);

TerminalCost=states;  TerminalCost(1)=2*kk;

Value=TerminalCost;
StateAction=zeros(Nstates,Nactions);

for ii=1:nn
        
    
    M=Smap-aa(ii)*Dmap;
    M(M<0)=-1;
    M=M+2;
    
 
    [Value,StateAction(:,ii)]=min(Value(M),[],2);
    StateAction(:,ii)=StateAction(:,ii)-1;
    
    
   
    if (Value(end)==0)
        LastStage=ii;
       break; 
    end
    
end

xx=zeros(1,nn);
TheState=kk;

for jj=LastStage:-1:1

 TheAction=StateAction(TheState+2,jj);
 xx(jj)=TheAction;
 TheState=TheState-aa(jj)*TheAction;
 
    
end