2
\$\begingroup\$

The purpose of this code is to use it to calculate some number theoretic function efficiently for the numbers from 1 to 1000000, but I am not sure if it will be really useful. The main purpose was to program a piece of C++ code because I have almost no experience in this language. So it would be nice if you could give me some feedback not to the algorithm but to its C++ implementation. There are also some things that are unclear to me: should the duplication of the block with the case statement be avoided? Should the usage of the global variables be avoided? Is the usage of an array appropriate or should I use some container datatype?

/*
if a positive integer number n has a prime divisor p>=sqrt(n), then n=p*q and p and q are uniquely defined and (p,q)=1 and q<=sqrt(n).
More general: if n is a positive integer number then there is a prime p such that n=(p^e)*q1*q2, such that (p,q1)=(p,q2)=(q1,q2)=1 and q1<=sqrt(n), q2<=sqrt(n). p,e,q1,q2 are not uniquely define.
this program calculates and stores such factors q1,q2,p^e, for every n between sqrt(limit) and limit and prints some sample output
*/
#include <iostream>
#include <cmath>
#include <cassert> 
const int limit=1000000;
int sqrtlimit;

int primefactor[limit]={0};

int smallfactor1[limit]={0};
int primepowerfactor[limit]={0};
int smallfactor2[limit]={0};


void init_sieve(void){
    /* 
    initializes the array primefactor
    primefactor[n]==n if n is a prime or n==1 or n==0, otherwise
    primefactor[n] is largest prime factor of n that 
    is smaller than sqrtkimit
    */
    sqrtlimit=int(std::pow(limit,.5))+1;
    for (int n=0;n<limit;n++){
        primefactor[n]=n;
    }
    for (int n=2;n<=sqrtlimit;n++){
        if (primefactor[n]==n){
            for (int i=n;i<limit;i+=n){
                primefactor[i]=n;
            }
        }
    }
}

void split_factors(void){
    /* 
    uses the initialized array primefactor
    sets the arrays smallfactor1, smalfactor2, primepowerfactor
    splits a number in three, pairwise relativly prime factors
    where the first and the third are smaller than sqrtlimit and the  second is a prime power
    */
    enum State {
        continue_first_factor,
        in_second_factor,
        continue_third_factor
        };
    for (int n=sqrtlimit; n<limit; n++){
        int current_prod=1;
        int quotient=n;
        int last_prime=0;
        int prime_exponent=0;
        int first_factor=1;
        int second_factor=1;
        int third_factor=1;
        int prime_power=1;
        int p;
        State state=continue_first_factor;
        while(primefactor[quotient]>1){
            p=primefactor[quotient];
            quotient/=p;
            /*
            if last prime power was found, 
            put it to first factor
            and restart calculating the next prime power
            */
            if (p!=last_prime){
                switch (state) {
                    case continue_first_factor:
                        first_factor*=prime_power;
                        break;
                    case in_second_factor:
                        second_factor=prime_power;
                        current_prod=1;
                        state=continue_third_factor;
                        break;
                    case continue_third_factor:
                        third_factor*=prime_power;
                        break;
                    default:
                        assert(0);
                        break;
                }
                last_prime=p;
                prime_power=p;
            }
            else {
                prime_power*=p;
            }
            if (state==continue_first_factor){
                current_prod*=p;
                if (current_prod>sqrtlimit){
                    state=in_second_factor;
                }
            }
            /* process the last prime power */
            if (quotient==1){
                switch (state) {
                    case continue_first_factor:
                        first_factor*=prime_power;
                        break;
                    case in_second_factor:
                        second_factor=prime_power;
                        current_prod=1;
                        state=continue_third_factor;
                        break;
                    case continue_third_factor:
                        third_factor*=prime_power;
                        break;
                    default:
                        assert(0);
                        break;
                }
            }
        }
        smallfactor1[n]=first_factor;
        primepowerfactor[n]=second_factor;
        smallfactor2[n]=third_factor;
    }
}

int main() {
    init_sieve();
    split_factors();  
    // some sample output
    std::cout<<"limit "<<limit<<std::endl;
    std::cout<<"sqrtlimit "<<sqrtlimit<<std::endl;
    std::cout<<std::endl;
    std::cout<<"print some prime factorization"<<std::endl;
    for (int n=limit-20;n<limit; n++){
        int m=n;
        std::cout << m << " = " ;
        bool first=true;
        while (m>1){
            if (first){
                first=false;
            }
            else {
                std::cout << " * ";
            }
            std::cout <<primefactor[m];
            m=m/primefactor[m];
        }
        std::cout <<std::endl;
    }
    std::cout<<std::endl<<"smallfactor * primepower * smallfactor"<<std::endl;
    std::cout<<"\tthe second factor is a prime power"<<std::endl;
    std::cout<<"\t'!' at the end prime power exponent > 1"<<std::endl;

    for (int n=limit-20;n<limit; n++){
        std::cout<<n<<" = "<<smallfactor1[n]<<" * "<<primepowerfactor[n]<<" * "<<smallfactor2[n];
        if (primefactor[primepowerfactor[n]]!=primepowerfactor[n]){
            std::cout<<" !";
        }
        std::cout<<std::endl;
    }
}
\$\endgroup\$

1 Answer 1

3
\$\begingroup\$

Questions

You’ve asked for a review of the C++ specifically, not the algorithm, so that’s what I’ll do. But first, I’ll answer your questions.

Should the duplication of the block with the case statement be avoided?

Yes.

Unless it’s really trivial—we’re talking like two lines, maybe three, max, and even then only if those lines are not complex or “clever”—you should never be shy about breaking repeated chunks of code out into separate functions.

In fact, forget “repeated”; anytime you have a subsection of an algorithm that seems somewhat self-contained, you should consider splitting it out. But that’s obviously especially true for repeated blocks.

There are a lot of good reasons for this:

  • Simplification: It makes any code that repeated block appears in shorter and simpler, and that’s usually a win for understanding and maintenance.
  • Reusability: If you’ve used it twice, there’s a good chance you could use it thrice, and you’d save time in the long run.
  • Modularization: When you break a chunk of code out into a separate function, you can treat function as a separate unit. You can test it separately. You can optimize it separately. You can work with it to improve it—and, by extension, anything that uses it—without worrying about breaking any code that uses it (assuming you’re properly testing it, of course… which you should).

And with any modern compiler, there is likely to be little to no cost to breaking code out into a separate function… especially if you make that function body visible everywhere (which will happen by default if you make it a constexpr function, which you should really try to do).

Should the usage of the global variables be avoided?

Yes.

Now, there are situations where global variables (or something similar, like a function that returns a reference to a static variable) are okay. But those are last resort situations: you should do everything in your power to avoid global variables, and only use them as a last resort, when every other option is worse.

By making your variables global, you have introduced four major problems:

  1. You’ve made your algorithm less efficient. Because of their nature, global variables must be allocated in memory somewhere. (Usually! Technically there is room for compilers to cheat and avoid actually allocating them. In practice, that’s extremely unlikely.) That means that the compiler no longer has the option of eliding them if it figures out a trick where it doesn’t need them. Not only that, they are very likely allocated somewhere in memory that is not close to where you happen to be working, meaning that it’s not going to be hot in cache.
  2. You’ve made your algorithm functionally impossible to use in concurrent code. This is the future, and everything is concurrent these days. You probably can’t get a general purpose consumer CPU with less than 4 cores these days, even on a cheap machine. (Hell, my router has 4 cores!) And that’s not even considering GPGPU stuff, which trades in hundreds or thousands of cores.
  3. You’ve made your algorithm harder to test. Global variables are difficult to mock, and damn near impossible to isolate, making testing a royal pain. You should really take writing testable code seriously; in fact, you should make it a top priority. Any code I see without tests, I consider garbage code, and I won’t allow it anywhere near any serious projects I’m working on.
  4. You’ve made your algorithm buggy. Even if you’re a genius programmer who only ever writes perfect code and never, ever introduces bugs… this code is unreliable. That’s because even if your code is perfect, someone else could fugger with those global variables and completely screw up the results of your algorithm. And you’d probably never know why: this would one of those horror-story bugs you hear about—the kind where the code works perfectly except on Bob’s computer and only then on Tuesdays and only if the Leafs happen to be playing that night.

We don’t say “avoid global variables” because we’re jerks who like to make up rules and force them on programmers so we can feel superior. Global variables are really, really bad. They are inefficient, they are untestable, and they are dangerous. Yes, there sometimes situations where they’re okay (like std::cout; that’s a global variable (well, technically, in a namespace, but close enough)). But you should only use them when you really, really must.

Is the usage of an array appropriate or should I use some container datatype?

This answer is less clear cut than the others, but I’d say that you should use std::vector instead.

First, if you must use arrays, I’d recommend using std::array rather than C arrays. std::array has safety features that make it less dangerous than a C array, as well as some really sweet usability enhancements.

However, the key issue with your usage is that those are some really freaking huge arrays. Take primefactor for example. It has a million elements. And each element is an int, which is probably going to be 4 bytes. So primefactor alone takes up 4 MB (roughly 3.8 MiB). Well, on Windows, last I checked (and, admittedly, I don’t check Windows stuff that frequently), the default stack size, at least for threads, is only 1 MiB. primefactor blows right past that… and you have four arrays that big. There is a decent likelihood that even if this program compiles, it will crash or behave erratically on some platforms, due to stack overflow.

For datasets that big, you should use dynamic allocation. That way, if there just isn’t enough memory, at least you will get predictable behaviour, rather than a crash. (Theoretically. In practice, the platform may over-commit and then start killing processes, and your program may end up crashing anyway. But at least you tried.)

So I’d say no to arrays for these huge datasets, and suggest using std::vector instead. In practice, not a single line of your code will have to change (yanno, other than those 4 array declarations). However… once you start using modern C++ tools like std::array or std::vector, you could make several improvements. I’ll get to that in the code review.

Code review

#include <cassert>
const int limit=1000000;

You should put some space between the headers and the actual code. Whitespace costs nothing, but makes a huge difference for readability.

const int limit=1000000;

const is the wrong word here; you want constexpr. const means “this value will not change”; constexpr means “this value is a constant”. Note that those are not the same thing. const variables can be set any time, even at runtime. The const just means that after it is set, it won’t change. But it’s still a variable, not a constant.

constexpr (when applied to variables) means a value is a true constant.

Also, I’d recommend putting some space around the =. Again, whitespace costs nothing, but makes things so much easier to read.

You might also consider using digit separators to make large numbers more readable:

constexpr auto limit = 1'000'000;

(The auto is optional, but honestly, when I write modern C++, everything is auto these days. Specific types usually don’t matter, interfaces do.)

const int limit=1000000;
int sqrtlimit;

int primefactor[limit]={0};

int smallfactor1[limit]={0};
int primepowerfactor[limit]={0};
int smallfactor2[limit]={0};

Since these are all global variables, they are all zero-initialized. Which means, you technically don’t need to write the ={0}… but it doesn’t hurt.

However, the fact that these are all zero-initialized is not really a good thing… because that’s just wasted cycles, since the first thing you do when the program actually gets going is to really initialize them with their real values.

This wasted initialization won’t be a problem if you don’t use global variables… which is what I’m about to suggest.

void init_sieve(void){

Okay, first, we don’t write void in empty parameter lists. That’s an ancient practice from C. It doesn’t apply to C++, and just clutters things up unnecessarily.

More importantly, though, when you have an “init” function… that’s a code smell. In C++, initialization and cleanup is usually handled automatically… in constructors and destructors. Those make initialization impossible to forget, and they make cleanup happen in every circumstance. And they’re usually no brainers; you don’t even really need to think about them, because they just work automatically.

Here’s what a class for your algorithm might look like, roughly:

class your_algo
{
    enum class factor_split_state
    {
        continue_first_factor,
        in_second_factor,
        continue_third_factor
    };

    static constexpr int _default_limit = 1'000'000;
    // better make sure a million can fit in an int!
    //
    // something like:
    //      static_assert(static_cast<long>(std::numeric_limits<int>::max()) >= 1'000'000L);
    //
    // or you could replace all the "int"s with "value_type",
    // and set "value_type" to "int" or "long", depending
    // on what fits:
    //      using value_type = std::conditional_t<
    //          (static_cast<long>(std::numeric_limits<int>::max()) >= 1'000'000L),
    //          int,
    //          long>;
    //
    //      static constexpr auto _default_limit = value_type(1'000'000L);
    //
    //      value_type limit = _default_limit;
    //      value_type sqrtlimit;
    //      std::vector<value_type> primefactor;
    //      // ... and so on...

public:
    your_algo() :
        your_algo(_default_limit)
    {}

    explicit your_algo(int lim) :
        limit{lim},
        sqrtlimit{int(std::pow(lim, 0.5)) + 1},
        primefactor(lim),
        smallfactor1(lim),
        primepowerfactor(lim),
        smallfactor2(lim)
    {
        // i've just duplicated your algorithm exactly, because you didn't
        // want that critiqued

        init_sieve();
        split_factors();
    }

    // i've just left all the variables public, because i don't understand
    // what you actually *need* to be public
    //
    // this is *not* good practice; you should decide what interface you
    // actually need
    int limit = _default_limit;
    int sqrtlimit;
    std::vector<int> primefactor;
    std::vector<int> smallfactor1;
    std::vector<int> primepowerfactor;
    std::vector<int> smallfactor2;

private:
    constexpr auto init_sieve() -> void
    {
        // for (int n=0;n<limit;n++){
        //     primefactor[n]=n;
        // }
        std::iota(primefactor.begin(), primefactor.end(), 0);

        for (auto n = 2; n <= sqrtlimit; ++n)
        {
            if (primefactor[n] == n)
            {
                for (auto i = n; i < limit; i += n)
                    primefactor[i] = n;
            }
        }
    }

    constexpr auto split_factors() -> void
    {
        // might want to make sure that:
        //      limit == primefactor.size()
        //      limit == smallfactor1.size()
        //      limit == primepowerfactor.size()
        //      limit == smallfactor2.size()

        for (auto n = sqrtlimit; n < limit; ++n)
        {
            auto current_prod = 1;
            auto quotient = n;
            auto last_prime = 0;
            auto prime_exponent = 0;
            auto first_factor = 1;
            auto second_factor = 1;
            auto third_factor = 1;
            auto prime_power = 1;
            auto state = factor_split_state::continue_first_factor;

            while (primefactor[quotient] > 1)
            {
                auto p = primefactor[quotient];
                quotient /= p;

                if (p != last_prime)
                {
                    std::tie(state, first_factor, second_factor, third_factor) =
                        _factor_calculator_thingy(state, prime_power, first_factor, second_factor, third_factor);

                    if (state == factor_split_state::in_second_factor)
                        current_prod = 1;

                    last_prime = p;
                    prime_power = p;
                }
                else {
                    prime_power *= p;
                }

                if (state == factor_split_state::continue_first_factor)
                {
                    current_prod *= p;
                    if (current_prod > sqrtlimit)
                        state = factor_split_state::in_second_factor;
                }

                /* process the last prime power */
                if (quotient == 1)
                {
                    std::tie(state, first_factor, second_factor, third_factor) =
                        _factor_calculator_thingy(state, prime_power, first_factor, second_factor, third_factor);

                    if (state == factor_split_state::in_second_factor)
                        current_prod = 1;
                }
            }

            smallfactor1[n] = first_factor;
            primepowerfactor[n] = second_factor;
            smallfactor2[n] = third_factor;
        }
    }

    // obviously this needs a better name
    static constexpr auto _factor_calculator_thingy(
        factor_split_state state,
        int prime_power,
        int first_factor,
        int second_factor,
        int third_factor
    ) noexcept
        -> std::tuple<factor_split_state, int, int, int>
    {
        switch (state)
        {
        case continue_first_factor:
            return {state, first_factor * prime_power, second_factor, third_factor};
        case in_second_factor:
            return {factor_split_state::continue_third_factor, first_factor, second_factor * prime_power, third_factor};
        case continue_third_factor:
            return {state, first_factor, second_factor, third_factor * prime_power};
        }
    }
};

Now your main() might look like:

auto main() -> int
{
    auto algo = your_algo{};

    std::cout << "limit " << algo.limit << '\n';
    std::cout << "sqrtlimit " << algo.sqrtlimit << '\n';
    std::cout << '\n';

    std::cout << "print some prime factorization\n";
    for (auto n = algo.limit - 20; n < algo.limit; ++n)
    {
        auto m = n;
        std::cout << m << " = " ;

        auto first = true;
        while (m > 1)
        {
            if (not std::exchange(first, false))
                std::cout << " * ";

            std::cout << algo.primefactor[m];
            m = m / algo.primefactor[m];
        }

        std::cout << '\n';
    }

    std::cout << "\nsmallfactor * primepower * smallfactor\n"
                 "\tthe second factor is a prime power\n"
                 "\t'!' at the end prime power exponent > 1\n";

    for (auto n = algo.limit - 20; n < algo.limit; ++n)
    {
        std::cout << n << " = " << algo.smallfactor1[n] << " * "
            << algo.primepowerfactor[n] << " * " << algo.smallfactor2[n];

        if (algo.primefactor[algo.primepowerfactor[n]] != algo.primepowerfactor[n])
            std::cout << " !";

        std::cout << '\n';
    }
}

Okay, back to the review.

sqrtlimit=int(std::pow(limit,.5))+1;

I don’t see the logic of using std::pow(limit, 0.5) rather than std::sqrt(limit). I presume you have a reason?

enum State {
    continue_first_factor,
    in_second_factor,
    continue_third_factor
    };

You should use the modern enum class rather than the old-school enum. It’s a bit more verbose, but much safer.

        int current_prod=1;
        int quotient=n;
        int last_prime=0;
        int prime_exponent=0;
        int first_factor=1;
        int second_factor=1;
        int third_factor=1;
        int prime_power=1;

When I see a long list of variables like this, it really suggests to me that you could use some more types. For example, a type for the factors:

struct factors_t
{
    int first = 1;
    int second = 1;
    int third = 1;
};

Then the list above reduces a bit to:

        auto current_prod = 1;
        auto quotient = n;
        auto last_prime = 0;
        auto prime_exponent = 0;
        auto factors = factors_t{};
        auto prime_power = 1;

That’s an improvement. It also simplifies the helper function:

    static constexpr auto _factor_calculator_thingy(factor_split_state state, int prime_power, factors_t factors) noexcept
        -> std::tuple<factor_split_state, factors_t>
    {
        switch (state)
        {
        case continue_first_factor:
            return {state, {factors.first * prime_power, factors.second, factors.third}};
        case in_second_factor:
            return {factor_split_state::continue_third_factor, {factors.first * prime_power, factors.second, factors.third}};
        case continue_third_factor:
            return {state, {factors.first, factors.second, factors.third * prime_power}};
        }
    }

You can probably bundle some of the other variables into types as well.

        int p;

Now this definitely doesn’t belong.

This would immediately trigger a “nope” from me in a code review, because it is an uninitialized variable. That’s bad.

But “fixing” it by adding an initializer is the wrong move, because the truth is this variable shouldn’t be here at all. It only needs to exist within the loop, so it should be created when it is needed:

        int third_factor=1;
        int prime_power=1;

        State state=continue_first_factor;
        while(primefactor[quotient]>1){
            auto p = primefactor[quotient]; // <-- this is where it needs to be created
            quotient/=p;

All the other variables need to maintain their values across loop iterations, so they have to be instantiated outside of the loop.

switch (state) {
    case continue_first_factor:
        first_factor*=prime_power;
        break;
    case in_second_factor:
        second_factor=prime_power;
        current_prod=1;
        state=continue_third_factor;
        break;
    case continue_third_factor:
        third_factor*=prime_power;
        break;
    default:
        assert(0);
        break;
}

Adding a default case here is perhaps being overly defensive. All modern compilers that I’m familiar with will automatically detect that the switch is exhaustive, because you’ve accounted for every value of the State enum.

So all that’s left is main(). The biggest crime that’s going on in main() is using std::endl. Don’t do that. Yes, I know that tutorials all over the place use it; well, I’m saying they all suck for doing so. Using std::endl is almost always wrong. (Or at the very least, pointlessly wasteful.) In fact, every single usage of std::endl in your code is wrong.

If you want a newline, just use \n.

    bool first=true;
    while (m>1){
        if (first){
            first=false;
        }
        else {
            std::cout << " * ";
        }
        // ... [snip] ...
    }

There’s a little trick for the pattern you’re using above. The key is std::exchange(). std::exchange() takes a variable as its first argument, and a value as its second: std::exchange(variable, value). What it does is:

  1. set the value of variable to value; and
  2. return the old value of variable.

And it does all that maximally efficiently.

So what you do is this:

    auto first = true;
    while (m > 1)
    {
        if (not std::exchange(first, false))
            std::cout << " * ";

        // ... [snip] ...
    }

On the first time through, the exchange() call sets first to false, and returns true (which is the old value of first). Since not true is false, the cout statement never fires.

On the second time through, and every time after, the exchange() sets first to false (pointlessly, because first is already false, but it doesn’t matter on any modern hardware), and returns false (which was the old value of first). not false is true, so the cout statement prints.

This is just one of those patterns you get to know; like for printing a bunch of values in a vector nicely:

auto const vec = std::vector{1, 2, 3, 4};

std::ranges::for_each(vec, [first=true](auto&& v) mutable
    {
        if (not std::exchange(first, false))
            std::cout << ", ";
        std::cout << v;
    });

// prints:
//      1, 2, 3, 4

Now let’s clean up the big printing block:

    std::cout<<std::endl<<"smallfactor * primepower * smallfactor"<<std::endl;
    std::cout<<"\tthe second factor is a prime power"<<std::endl;
    std::cout<<"\t'!' at the end prime power exponent > 1"<<std::endl;

First we remove all the std::endls (and space out the tokens, for readability):

    std::cout << '\n' << "smallfactor * primepower * smallfactor" << '\n';
    std::cout << "\tthe second factor is a prime power" << '\n';
    std::cout << "\t'!' at the end prime power exponent > 1" << '\n';

Next, we merge all the '\n's into the strings:

    std::cout << "\nsmallfactor * primepower * smallfactor\n";
    std::cout << "\tthe second factor is a prime power\n";
    std::cout << "\t'!' at the end prime power exponent > 1\n";

Next, instead of 3 cout statements, we merge them into one:

    std::cout << "\nsmallfactor * primepower * smallfactor\n"
        << "\tthe second factor is a prime power\n"
        << "\t'!' at the end prime power exponent > 1\n";

Next, we take advantage of the fact that C++ will automatically merge adjacent string literals, and turn those 3 operator<<() calls into 1:

    std::cout << "\nsmallfactor * primepower * smallfactor\n"
        "\tthe second factor is a prime power\n"
        "\t'!' at the end prime power exponent > 1\n";

We could leave it at this, which will already be considerably faster than the original code, but there is another option you could consider: raw strings:

    std::cout << R"(
smallfactor * primepower * smallfactor
    the second factor is a prime power
    '!' at the end prime power exponent > 1
)";

That’s about all I can think to recommend without a complete overhaul of the algorithm, and code structure. A starting point for that would be to use a class like I demonstrated above. After that, perhaps the next step would be to look at breaking out parts of the process into independent functions or types—in particular, it probably shouldn’t be necessary to work with 4 arrays. Then, maybe in the long run, look at opportunities for parallelization.

\$\endgroup\$
1
  • \$\begingroup\$ Thanks for the effort you made. I am impressed. Every line of this post is useful to me. I will revise my program and try to reconsider and incorporate the suggestions you have made. \$\endgroup\$
    – miracle173
    Commented Apr 24, 2021 at 9:09

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.