for文・if文を用​いて、条件を満たすパ​ターンの組み合わせを​出したい

MATLAB初心者で手も足も出ないので、教えてください。

【やりたいこと】

条件を満たすような5つの変数のパターンの組み合わせを出したい。

　任意の値：N＝〇

　変数：Za、Zb、Zc、Zd、Ze　（全て整数で、値は5～150の範囲）

　条件：Zc/Za＞1

　　　　(Zb・Ze)/(Zc・Zd)＞１

　　　　(Zb・Ze)/(Za・Zd)＞１

　　　　X=(Za＋Zc)/N　（整数）

　　　　Y=(Za・Zd＋Zb・Ze)/(N・(ZbとZdの最大公約数))　(整数)

　　　　Z＝(Zc・ZdーZb・Ze)/(N・(ZbとZdの最大公約数))　(整数)

出力イメージとしては、

　（Za、Zb、Zc、Zd、Ze）＝（11、11、34、12、38）

　　　　　　　　　　　　　　（10、32、30、74、72）

　　　　　　　　　　　　　　　　　　　　・

　　　　　　　　　　　　　　　　　　　　・

　　　　　　　　　　　　　　　　　　　　・

よろしくお願いします。

　　　

[Hernia Baby]

Zaなどはそれぞれベクトルですか？また、条件を満たすときにどのような操作をするのでしょうか？

[実香]

ベクトルではなく、数値です。

条件が6つあると思うのですが、

条件をひとつ満たせば、次の条件を満たすかトライ、満たさなければ、Za～Zeの変数を変えてもう一度for文を回していく。といった感じでしょうか。

最終的な出力としては、全ての条件を満たした5つの値のパターンがズラズラと出てきたらいいなと考えています。

本当に間違っていると思うのですが、私はこのように書いてみました↓

clear variables;

N=5; %%ピニオンの歯数

for Za=5:15
    for Zb=5:15
        for Zc=5:40
            for Zd=5:40
                for Ze=5:40
                    
                    U=Zc/Za;
                    V=(Zb*Ze)/(Zc*Zd);
                    W=(Zb*Ze)/(Za*Zd);
                    
                    if U>1
                        continue
                    end
                    
                    if V>1
                        continue
                    end
                    
                    if W>1
                        continue
                    end
                    
                    Q=gcd(Zb,Zd);%ZbとZdの最大公約数
                    
                    X=(Za+Zc)/N;
                    Y=(Za*Zd+Zb*Ze)/(N*Q);
                    Z=(Za*Zd-Zb*Ze)/(N*Q);
                    
                    if  X==floor(X)
                       continue
                    end
                    
                    if  Y==floor(Y)
                        continue
                    end
                    
                    if  Z==floor(Z)
                        continue
                    end
                    
                    A=[Za Zb Zc Zd Ze];          
                 
                    
                end
            end
        end
    end
end

 disp=(num2str(A));

Answer 1

[交感神経優位なあかべぇ]

forとifを使って、素直に実装しようとすると、下記のようになりますかね。ただ、これだといつまで経っても処理が終わらないので、パフォーマンスを改善した実装が必要になると思います。

（ちなみにZ＝(Zc・ZdーZb・Ze)/(N・(ZbとZdの最大公約数))　(整数)の条件のZcはZcでいいですかね？（コメントのコード文ではZcではなくZaだったので、どちらかが書き間違えかなと））

N = 5; % ピニオンの歯数
A = zeros(0, 5);
for Za = 5 : 150
    for Zb = 5 : 150
        for Zc = 5 : 150
            for Zd = 5 : 150
                for Ze = 5 : 150
                    if ~(Zc / Za > 1)
                        continue;
                    end
                    
                    if ~((Zb * Ze) / (Zc * Zd) > 1)
                        continue;
                    end
                    
                    if ~((Zb * Ze) / (Za * Zd) > 1)
                        continue;
                    end
                    X = (Za + Zc) / N;
                    if ~(X == floor(X))
                        continue;
                    end

                    Q = gcd(Zb, Zd);%ZbとZdの最大公約数
                    Y = (Za * Zd + Zb * Ze) / (N * Q);
                    if ~(Y == floor(Y))
                        continue;
                    end
                    
                    Z = (Zc * Zd - Zb * Ze) / (N * Q);                    
                    if ~(Z == floor(Z))
                        continue;
                    end
                    
                    A(end + 1, :) = [Za, Zb, Zc, Zd, Ze];
                end
            end
        end
    end
end
A % 結果の表示

上記コードの途中までの結果の答え

[交感神経優位なあかべぇ]

事前に計算できる個所は計算しておき、for文内での計算はできるだけ最小限にするような形で実装しようとしましたが、答えが100万通りを超えていますので、中断しました。

N = 5; % ピニオンの歯数
A = zeros(1000000, 5);
prodList = (5 : 150)' .* (5 : 150);
ZxList = 5 : 150;
gcdList = zeros(146, 146);
for idx = 1 : 146
    gcdList(idx, :) = gcd(ZxList(idx), 5 : 150);
end
prodQNList = gcdList .* N;

setAIdx = 1;
for ZaIdx = 1 : 146
    for ZcIdx = ZaIdx + 1 : 146 % Zc/Za＞1を満たすため、左式を使用。
        X = (ZxList(ZaIdx) + ZxList(ZcIdx)) / N;
        if ~(X == floor(X))
            continue;
        end

        for ZbIdx = 1 : 146
            for ZdIdx = 1 : 146
                for ZeIdx = 1 : 146
                    if ~(prodList(ZbIdx, ZeIdx) > prodList(ZcIdx, ZdIdx))
                        continue;
                    end
                    
                    if ~(prodList(ZbIdx, ZeIdx) > prodList(ZaIdx * ZdIdx))
                        continue;
                    end

                    prodQN = prodQNList(ZbIdx, ZdIdx); % ZbとZdの最大公約数
                    Y = (prodList(ZaIdx, ZdIdx) + prodList(ZbIdx, ZeIdx)) / prodQN;
                    if ~(Y == floor(Y))
                        continue;
                    end
                    
                    Z = (prodList(ZcIdx, ZdIdx) - prodList(ZbIdx, ZeIdx)) / prodQN;                    
                    if ~(Z == floor(Z))
                        continue;
                    end
                    
                    A(setAIdx, :) = [ZxList([ZaIdx, ZbIdx, ZcIdx, ZdIdx, ZeIdx])];
                    if setAIdx >= 1000000
                        error('変数Aの容量が足りません。');
                    end
                    setAIdx = setAIdx + 1;
                end
            end
        end
    end
end

変数Aの容量が足りません。

A(setAIdx : end, :) = []; % 余ったデータの削除
A % 結果の表示

[実香]

ご回答ありがとうございます！！

＞ちなみにZ＝(Zc・ZdーZb・Ze)/(N・(ZbとZdの最大公約数))　(整数)の条件のZcはZcでいいですかね？

Zcで合ってます！！コードを間違えていました…

for文、if文を用いた場合、どうなるのかとても分かりやすかったです。

ありがとうございます。

値が大きくなりすぎないように数値を見直して、この方法で行ってみたいと思います！

Answer 2

[Hernia Baby]

まずはデータを用意

N = 1;
Z = randi([5 150],5000,5);
Za = Z(:,1);
Zb = Z(:,2);
Zc = Z(:,3);
Zd = Z(:,4);
Ze = Z(:,5);

条件を作る

idx1 = Zc./Za > 1;
idx2 = (Zb.*Ze)./(Zc.*Zd)>1;
idx3 = (Zb.*Ze)./(Za.*Zd)>1;

X = (Za+Zc)./N;
Y = (Za.*Zd+Zb.*Ze)./(N.*gcd(Zb,Zd));
Z = (Zc.*Zd-Zb.*Ze)./(N.*gcd(Zb,Zd));
% 整数か判定
idx4 = X == floor(X);
idx5 = Y == floor(Z);
idx6 = Z == floor(Z);

全てを満たす条件

idx = idx1 & idx2 & idx3 & idx4 & idx5 & idx6;

条件に合うものを抽出

Z(idx,:)

ans =

  0×1 empty double column vector

[実香]

ご回答ありがとうございます。

コードを入力してみましたが、上手く表示されません…

どのようにしたら良いでしょうか？？

また、forとifの使用した場合はどうなりますでしょうか…？

[Hernia Baby]

存在しないのはそのようなZa,Zb,Zc,Zd,Zeがここには存在しないからです。

そして新しいコメントをみたのですが、ここでいう組み合わせとは

a = (1:3)';
b = (1:3)';
C = [a,b]

C = 3×2

1 1 2 2 3 3

でいうと

for ii = 1:height(C)
    a(ii)*b(ii) > 5
end

ans = logical
   0

ans = logical
   0

ans = logical
   1

のように3回テストをまわすのではなく、

for ii = 1:length(a)
    for jj = 1:length(b)
       a(ii)*b(jj) > 5
    end
end

ans = logical
   0

ans = logical
   0

ans = logical
   0

ans = logical
   0

ans = logical
   0

ans = logical
   1

ans = logical
   0

ans = logical
   1

ans = logical
   1

のように3×3=9回分テストを回すということなのでしょうか？

[実香]

コメントありがとうございます！！

後者の認識のつもりでした…

分かりづらくてごめんなさい…

ただ、今回の場合forで回すと計算回数も多いため、Hernia Babyさんが最初に教えてくださった、

ランダムで数値を取り出したのち、行列で計算するのがとてもシンプルだと思いました。

とても勉強になりました！！