output=amicable(N2,varargin)

%Function to compute amicable pairs between two intergers, upper & lower limits.
%
%Amicable numbers are two different numbers so related that the sum of
%the proper divisors of the one is equal to the other, one being considered
%as a proper divisor but not the number itself.
%Example of such a pair is (220, 284); for the proper divisors of
%220 are 1, 2, 4, 5, 10, 11, 20, 22, 44, 55, 110, sum of which is 284;
%and the proper divisors of 284 are 1, 2, 4, 71, 142, sum of which is 220.
%
%Check: http://en.wikipedia.org/wiki/Amicable_number
%
%Put the 'amicable.m' file into your MATLAB path or working directory and
%call the following function syntax
%
% Function syntax:
% output=amicable(N2,N1)
% 
% where, N2 is the upper limit and N1 the lower limit, both positive integers. 
% If only one number is given, it is considered to be the upper limit and 
% the default lower limit is then 220 (as no amicable pair exists below 220).
% If two numbers are given, the upper limit (1st number) must be greater than
% the lower limit (2nd number).
% 
% The 'output' would be a matrix having two columns, for the amicable pairs. 
% The rows show the unique pairs. If no such pairs are found empty matrix is returned.
%
%Written by (C) Abhisek Ukil, [e-mail], 26.09.2008


function output=amicable(N2,varargin)

%Check arguments.
if numel(varargin)==1      
    N1=varargin{1};
else
    N1=220; %default lower limit
end


%error checks
if numel(N1)~=1 || numel(N2)~=1
    error('MATLAB:factor:NonScalarInput','Limits must be scalar.'); 
end

if ischar(N1) || ischar(N2)
    error('MATLAB:factor:InputNotPosInt', 'Limits must be numbers.');
end

if (N1 < 0) || (floor(N1) ~= N1) || (N2 < 0) || (floor(N2) ~= N2)
    error('MATLAB:factor:InputNotPosInt', 'Limits must be positive integers.');
end

if N1>N2,error('MATLAB:factor:Limiterror','Upper limit (1st number) must be greater than lower limit (2nd number).'),end

if N2-N1>1e5,warning('Difference between two limits too high (>1e5). Might be very slow !!!'),end


output=[];

if N2<220, return, end

for i=N1:N2

    n1=i;

    % do not consider repeatations
    if ~isempty(output)
        if ~isempty(intersect(output(:,2),n1,'rows')),continue,end
    end

    % check for amicable pair condition
    f1=allfactor(n1);
    s1=sum(f1);

    n2=s1;
    if n2~=n1
        f2=allfactor(n2);
        s2=sum(f2);

        if s2==n1 & ~isprime(n1) & ~isprime(n2)
            output(end+1,1:2)=[n1 n2];
        end
    end

end

%------------------------------------------------------------

%function to calculate all factors, contribution of 'Jos x@y.z'
function all_factor=allfactor(N)

all_factor = 1:(N/2) ; 
t = N ./ all_factor ; 
all_factor = all_factor(t==fix(t)) ;