Problem193 平方因子を含まない数

問題

正の整数 n が, 任意の素数の2乗によって割り切れないとき,
n を"平方因子を含まない"(squarefree)と呼ぶ.
1, 2, 3, 5, 6, 7, 10, 11 は平方因子を含まないが, 4, 8, 9, 12 は含む.

2の50乗未満で平方因子を含まない数はいくつあるか?

ソース

using System;
using System.Collections.Generic;

class Program
{
    static long[] mHeihouSuuArr;

    //エラトステネスの篩をしてから、素数の平方の配列を作成
    static void DeriveHeihouSuuArr()
    {
        var IsSosuuArr = new System.Collections.BitArray((int)mJyousuu25 + 1);

        for (int I = 2; I <= IsSosuuArr.Count - 1; I++) {
            IsSosuuArr[I] = true;
        }
        for (int I = 2; I * I <= IsSosuuArr.Count - 1; I++) {
            if (I != 2 && I % 2 == 0) continue;

            if (IsSosuuArr[I]) {
                for (int J = I * 2; J <= IsSosuuArr.Count - 1; J += I) {
                    IsSosuuArr[J] = false;
                }
            }
        }
        var HeihouSuuList = new List<long>();
        for (int I = 2; I <= IsSosuuArr.Count - 1; I++)
            if (IsSosuuArr[I]) HeihouSuuList.Add((long)I * I);

        mHeihouSuuArr = HeihouSuuList.ToArray();
    }

    static long mJyousuu25 = 1;
    static long mJyousuu50 = 1;

    static void Main()
    {
        var sw = System.Diagnostics.Stopwatch.StartNew();

        for (int I = 1; I <= 25; I++) mJyousuu25 *= 2;
        for (int I = 1; I <= 50; I++) mJyousuu50 *= 2;
        Console.WriteLine("2の25乗={0}", mJyousuu25);
        Console.WriteLine("2の50乗={0}", mJyousuu50);

        DeriveHeihouSuuArr();

        Console.WriteLine("Answer={0}。経過時間={1}", Solve(mJyousuu50), sw.Elapsed);
    }

    //指定した上限未満の平方因子を含まない数の個数を返す
    static long Solve(long pJyougen)
    {
        //2*2の倍数、または、3*3の倍数、または、5*5の倍数、または、7*7の倍数・・・
        //の集合の要素数を、包除原理で求める
        long WaSyuugouCnt = 0;
        for (int I = 0; I <= mHeihouSuuArr.GetUpperBound(0); I++) {
            int SelectCnt = I + 1;
            bool HasYouso = false;

            foreach (long EachProd in DeriveProdEnum(SelectCnt)) {
                HasYouso = true;
                long CurrCnt = (pJyougen - 1) / EachProd;

                //奇数なら加算し、偶数なら減算
                if (SelectCnt % 2 == 1)
                    WaSyuugouCnt += CurrCnt;
                else WaSyuugouCnt -= CurrCnt;
            }
            if (HasYouso == false) break;
        }
        return pJyougen - 1 - WaSyuugouCnt;
    }

    struct JyoutaiDef
    {
        internal int CurrInd;
        internal int SelectedCnt;
        internal long ProdVal;
    }

    //引数の個数を選んだ素数の平方の積の列挙を返す
    static IEnumerable<long> DeriveProdEnum(int pNeedCnt)
    {
        var Stk = new Stack<JyoutaiDef>();
        JyoutaiDef WillPush;
        WillPush.CurrInd = 0;
        WillPush.SelectedCnt = 0;
        WillPush.ProdVal = 1;
        Stk.Push(WillPush);

        while (Stk.Count > 0) {
            JyoutaiDef Popped = Stk.Pop();

            //クリア判定
            if (Popped.SelectedCnt == pNeedCnt) {
                yield return Popped.ProdVal;
                continue;
            }

            for (int I = Popped.CurrInd; I <= mHeihouSuuArr.GetUpperBound(0); I++) {
                WillPush.CurrInd = I + 1;
                WillPush.SelectedCnt = Popped.SelectedCnt + 1;

                //葉ノードが、2の50乗未満である、必要条件のチェック
                int RestCnt = pNeedCnt - Popped.SelectedCnt;
                long wkSahen = mJyousuu50;
                wkSahen /= Popped.ProdVal;
                for (int J = 1; J <= RestCnt; J++)
                    wkSahen /= mHeihouSuuArr[I];
                if (wkSahen < 1) break;

                WillPush.ProdVal = Popped.ProdVal * mHeihouSuuArr[I];

                //2の50乗未満である、十分条件のチェック
                if (WillPush.ProdVal >= mJyousuu50) break;

                Stk.Push(WillPush);
            }
        }
    }
}

実行結果

2の25乗=33554432
2の50乗=1125899906842624
Answer=684465067343069。経過時間=00:00:03.0597556

解説

   2の50乗 < A * B
⇔ 2の50乗 / A / B < 1
ですので、
2の50乗をAで割って、端数を切り捨て、
さらにBで割って、端数を切り捨てた数が1未満であることは、

2の50乗 < A * B
であるための必要条件であることを使ってます。