Evaluated Powers of Bases if gcd of 3 bases is 1 without using tons of for loops

Question

I'm trying to make a program that will evaluate the gcd of 3 numbers (the bases), and if the GCD of base1 and base2, base2 and base3, and base3 and base1 are all equal to one, then evaulate the bases to exponents in a range. Basically, what I need to do is figure out if their GCD's are equal to one, then calculate the numbers to powers. Here's what it would look like:

bases = 150
powers = 150
base1 = all numbers 1-bases
base2 = all numbers 1-bases
base3 = all numbers 1-bases
if the GCD of all combinations = 1
    do base1^all numbers 3-powers
    do base2^all numbers 3-powers
    do base2^all numbers 3-powers
    then store all of those in an array

Now, I've tried using the dreaded for loops, but it is very slow and I'm not considering it a solution. It only works quickly if bases and powers are 10 or below. How can I do this without using for loops? Or, if I have to use a for loop, how can I reduce the number that is used? Could I computate what numbers together's 3 GCD combinations equal 1? The for loop that I've tried is below:

for i = 1:numbers
    for j = 1:numbers
        for k = 1:numbers
            if gcd(i,j) == 1 && gcd(i,k) == 1 && gcd(k,j) == 1
                for a = 3:powers
                    for b = 3:powers
                        for c = 3: powers
                            x = i^a
                            y = j^b
                            z = k^c
                        end
                    end
                end
            end
        end
    end
end

Answer 1

I think this problem can be solved in a vecctorized manner as follows

%Generate all possible triplets.
[x1,x2,x3]=ndgrid(1:3,1:3,1:3);
v=[x1(:) x2(:) x3(:)];

Now define an anonymous function.

gcd3_test=@(a,b,c)(gcd(a,b)==1 & gcd(b,c)==1 & gcd(a,c)==1)
gcdTest=gcd3_test(v(:,1),v(:,2),v(:,3));
v=v(gcdTest,:);

Similarly generate the triplets for all your powers.

[x1,x2,x3]=ndgrid(3:10,3:10,3:10);
p=[x1(:) x2(:) x3(:)];

Then I guess you will have to use the for loop (but only one) as:

Important: I assume if you run the for loop for powers 3:150, and you need to store all x,y,z then you will need lots of memory (35034 GB, even with single precision). You are not doing that in your code.

So do not try the following for loop.

%DO NOT RUN
for i=1:size(v,1)
    vReplicated=repmat(v(i,:),size(p,1),1);
    v_RaisedTo_P{i}=vReplicated.^p;
end

Note: Observe that you don't care about the order when conducting your GCD test (i.e. your if condition. So I think you can filter out lot of triplets, but this will affect your power calculation.

allTriplets=sort(allTriplets,2);
allTriplets=unique(allTriplets,'rows'); %27 triplets reduce to 10