Example:
Data Types:









Output Arguments
collapse all
— Greatest common divisor
real, nonnegative integer values
Greatest common divisor, returned as an array of real nonnegative
integer values. Thus, effectively the Euclidean algorithm, applied to a vector of numbers…
d = min(A);
while
true
r = mod(A,d);
if
~any(r)
break
end
% find the smallest nonzero element of r
r(r == 0) = inf;
d = min(r);
end
When the loop terminates, which it must do as long as all of the elements of A are positive integers, then d is the GCD of the group.
Solution 1:
If loops are not your preference, would you be interested in using recursive functions instead?
iif = @(varargin) varargin{2 * find([varargin{1:2:end}], 1, 'first')}();
gcdrec=@(v,gcdr) iif(length(v)==1,v, ...
v(1)==1,1, ...
length(v)==2,@()gcd(v(1),v(2)), ...
true,@()gcdr([gcd(v(1),v(2)),v(3:end)],gcdr));
mygcd=@(v)gcdrec(v(:)',gcdrec);
A=[0,0,0; 2,4,2;2,0,8];
divisor=mygcd(A);
A=A/divisor;
The function defined by
iif
enables the creation of a recursive function,
gcdrec
, which can be used to compute the greatest common divisor of an array. To accomplish this,
iif
implements an inline conditional construct that tests the first argument to determine if it is
true
. If it is, the function returns the second argument. If not, it progresses to the third argument and repeats the process until it finds an argument that satisfies
true
. It is important to use
@()
to safeguard recursive functions and other quantities used within them to prevent errors.
By implementing
iif
, the recursive function
gcdrec
operates in the following manner:
 It returns the input vector as it is, if it is a scalar.
 In case the initial element of the vector is 1, recovery becomes impossible and the function returns 1 quickly, especially for bigger matrices.

If the length of the input vector is two, the function will calculate the greatest common divisor using
gcd
.  Otherwise, it invokes itself using a truncated vector where the initial two elements are replaced by their highest common factor.
The purpose of
mygcd
is to serve as a userfriendly interface.
The speed should be satisfactory, with the only potential issue being the recursion depth for exceptionally large problems. To compare with @Adriaan’s looping version, I conducted a brief timing evaluation using
A=randi(100,N,N)50
,
N=100
,
N=1000
,
N=5000
,
tic
, and
toc
.

N=100
includes the following:
–  looping 0.008 seconds
 recursive 0.002 seconds

N=1000
includes the following items:
–  looping 0.46 seconds
 recursive 0.04 seconds

N=5000
includes the following items:
–  looping 11.8 seconds
 recursive 0.6 seconds
–
–
–
The update reveals that the recursion limit, which is set at 500 by default, was not exceeded due to the absence of a common divisor in the data. However, doubling a random matrix and applying
N=100
would cause the recursion limit to be reached quickly, making it unsuitable for large matrices. On the other hand, @Adriaan’s solution works well for small matrices.
In an attempt to address the recursion limit issue, I experimented with rewriting the input vector to a smaller size at each recursive step. However, this approach, although effective, resulted in a significant decrease in speed, taking 2 seconds for
N=100
and a whopping 261 seconds for
N=1000
. I am still searching for a solution that strikes a balance between a reasonably large matrix size and an acceptable runtime.
Solution 2:
A = [0,0,0; 2,4,2;2,0,8];
B = 1;
kk = max(abs(A(:))); % start at the end
while B~=0 && kk>=0
tmp = mod(A,kk);
B = sum(tmp(:));
kk = kk  1;
end
kk = kk+1;
To compute the sum of all elements in your matrix after applying the modulus operation,
B
initializes a counter. The
kk
is used to run through integers while
mod(A,kk)
computes the modulus after division for each element in
A
. If all elements in the matrix are divisible by 2,
mod(A,kk)
will return 0 for each element. The resulting modulomatrix is transformed into a single column by
sum(tmp(:))
and summed to obtain a number. If this number is 0, it indicates that all elements in
A
are divisible by
kk
. The loop stops as soon as this happens, and the number in
kk
is the common divisor. Since
kk
is decreased every count, it is one value too low, thus one is added. Even though this might not be the quickest approach, it will suffice for now.
The loop has been modified to run in reverse as you are searching for the largest CD and not the smallest one. In case you have a matrix similar to
[4,8;16,8]
, the loop will halt at
2
instead of
4
. I apologize for the mistake, but this solution is now functional, despite the fact that the other two solutions presented here are much quicker.
One way to perform matrix division is shown below.
divided_matrix = A/kk;
Solution 3:
I concur, loops are not favorable to me either. We should eliminate them.
unqA = unique(abs(A(A~=0))).'; %//'
col_extent = [2:max(unqA)]'; %//'
B = repmat(col_extent,1,numel(unqA));
B(bsxfun(@gt,col_extent,unqA)) = 0;
divisor = find(all(bsxfun(@times,bsxfun(@rem,unqA,B)==0,B),2),1,'first');
if isempty(divisor)
out = A;
else
out = A/divisor;
end
Sample runs
Case #1:
A =
0 0 0
2 4 2
2 0 8
divisor =
2
out =
0 0 0
1 2 1
1 0 4
Case #2:
A =
0 3 0
5 7 6
5 0 21
divisor =
1
out =
0 3 0
5 7 6
5 0 21