How to find a sequence like this?-CodePudding

The problem is as follows.

You are given a natural number n (2 <= n <= 10⁹) and you have to build a sequence with its divisors (d₁, d₂, d₃ ...) so that for every pair of subsequent divisors one of these two rules is complied with. If d_i < d_{i 1} then d_{i 1}%d_i must be 0 and d_{i 1}/d_i must be a prime number and the other way around If d_i > d_{i 1} then d_i%d_{i 1} must be 0 and d_i/d_{i 1} must be a prime number.

Let me give you an example. n = 20 so {1, 2, 4, 20, 10, 5} can be an answer because

1 < 2 && 2%1 == 0 && 2/1 = 2 which is a prime number

2 < 4 && 4%2 == 0 && 4/2 = 2 which is a prime number

4 < 20 && 20%4 == 0 && 20/4 = 5 which is a prime number

20 > 10 && 20 == 0 && 20/10 = 2 which is a prime number

10 > 5 && 10%5 == 0 && 10/5 = 2 which is a prime number

#include <vector>
#include <bitset>
#include <iostream>
#include <algorithm>
#include <cassert>
#include <iterator>

std::size_t constexpr sz_max{ 1500 };
std::vector<int> div_list, ans;
std::vector<std::vector<std::size_t>> edge_list;
std::bitset<sz_max> is_used;

bool is_prime(int src)noexcept {
    if (src == 2 || src == 3) {
        return true;
    }
    if (src <= 1 || (src & 1) == 0 || src % 3 == 0) {
        return false;
    }
    for (int div{ 5 }; div * div <= src; div  = 6) {
        if (src % div == 0 || src % (div   2) == 0) {
            return false;
        }
    }
    return true;
}

bool is_valid(int lhs, int rhs)noexcept {
    return (lhs < rhs&& rhs% lhs == 0 && is_prime(rhs / lhs)) ||
        (lhs > rhs && lhs % rhs == 0 && is_prime(lhs / rhs));
}

bool is_sorted(std::vector<int> const& src)noexcept {
    for (std::size_t i{ 0 }; i   1 != src.size();   i) {
        if (!is_valid(src[i], src[i   1])) {
            return false;
        }
    }
    return true;
}

void DFS(std::size_t i)noexcept {
    ans.emplace_back(div_list[i]);
    is_used.set(i, true);
    for (std::size_t curr : edge_list[i]) {
        if (!is_used.test(curr)) {
            DFS(curr);
        }
    }
    if (ans.size() != div_list.size()) {
        ans.pop_back();
        is_used.set(i, false);
    }
}

int main() {
    // 1
    int n;
    std::cin >> n;
    for (int div{ 1 }; div * div <= n;   div) {
        if (n % div == 0) {
            div_list.emplace_back(div);
            if (div * div != n) {
                div_list.emplace_back(n / div);
            }
        }
    }
    // 2
    std::sort(div_list.begin(), div_list.end());
    // 3
    edge_list.resize(div_list.size());
    for (std::size_t i{ 0 }; i != div_list.size();   i) {
        for (std::size_t j{ 0 }; j != div_list.size();   j) {
            if (i != j && is_valid(div_list[i], div_list[j])) {
                edge_list[i].emplace_back(j);
            }
        }
    }
    // 4
    DFS(0);
    // 5
    assert(is_sorted(ans));
    // 6
    std::copy(ans.cbegin(), std::prev(ans.cend()),
        std::ostream_iterator<int>{ std::cout, " " });
    std::cout << ans.back() << '\n';
    return 0;
}

In //1 I find out the divisors of n and save them in div_list.

In //2 I sort the divisors in ascending order.

In //3 I build the edge list of a graph consisting of the divisors of n as vertices. edge_list[i] is a vector which keeps the indices of the viable neighbours of the divisor div_list[i].

In //4 I build the sequence ans using DFS algorithm which uses the bitset is_used to mark the visited vertices. In this block

    if (ans.size() != div_list.size()) {
        ans.pop_back();
        is_used.set(i, false);
    }

I'm implementing backtracking in case of running into dead end. The constant sz_max is high enough I think to cover the cases where n has too many divisors.

In //5 I double check the sequence ans.

In //6 I output ans.

The problem is that DFS tries every possible path until it finds the correct one. For some numbers the algorithm runs very fast (2162160, 1000000000) but for others it runs for ages (not because It is wrong but because of the recursive calls It has to make). For example 223092870, 555710400 and 707691600.

Question 1: Is is somehow possible to improve the algorithm ?

Question 2: If not and I have to select different strategy which one It has to be ?

CodePudding user response：

This is indeed a search of Hamiltonian path in the lattice of divisors of n. We can prove that such a Hamiltonian path exists, because the lattice of divisors is a product of line graphs for each divisor (see, for example, this answer)

Therefore, you can use this proof to find a path as follows.

Find all prime divisors d_i of n, and their powers p_i.
Start with 1.
If you have constructed a path P using only prime divisors d_j, j<i, you can construct a path additionally using prime divisor d_i: P, d_i * P', d_i^2 * P, d_i^3 * P', ..., d_i^p_i * P. Here I denote by P' a reversed path P, and by k * P a path where each element is multiplied by k.
Repeat step 3 for all prime divisors.

To clarify meaning of step 3. Suppose that you have a path P = 1,2,4, and your next d_i = 5, p_i = 4. Then you get a path (1,2,4), (20,10,5), (25,50,100), (500,250,125). I enclosed each part of the new path in brackets, first part is simply P, second - 5 * P', third - 25 * P, fourth - 125 * P'.

This algorithm is linear in number of divisors.

Here is the main part of solution (I haven't written C in ages, excuse me for its style):

ans.emplace_back(1);
for (int i=0; i<prime_divs.size(); i  ) {
    int d = prime_divs[i];
    int n = ans.size();
    int pos = n-1;
    for (int j=0; j<prime_pows[i]; j  ) {
        for (int k=pos; k>=pos-n 1; k--) {
            ans.emplace_back(ans[k] * d);
        }
        pos  = n;
    }
}

I tested it on your examples, and it outputs a list with all divisors in the required order.

CodePudding user response：

As @stef correctly points out in a comment, the exercise is basically asking you to find a Hamiltonian path through the divisibility lattice of a given number n:

A Hamiltonian path is a path through a graph which visits every vertex exactly once. (If there is also an edge from the last vertex in the path to the first vertex in the path, then it's a Hamiltonian cycle.)
The divisibility lattice is a partial ordering using divisibility, so that a precedes b if a is a divisor of b. But I don't think the lattice actually helps here.

We conceptually construct the graph of valid successors between divisors, which is the same graph explored by the backtracking solution in the OP, using the prime factorization of n. That's a list of primes p₁, p₂, p₃, … p_k and corresponding exponents a₁, a₂, a₃, … a_k such that n = p₁^a₁ × p₂^a₂ × p₃^a₃ × … × p_k^a_k. (To minimize k, we eliminate primes which are not factors of n, so all the exponents are at least 1.)

Every divisor of n can be represented as a product using the same set of primes with corresponding exponents in the inclusive range [0, a_i], and every such list of exponents represents a divisor. Since the list of primes is the same for every vertex in the graph, we just label each vertex with its list of exponents. Then we draw an edge between any two vertices whose exponent lists differ by exactly 1 in exactly one position, which corresponds to the requirement that the ratio between successive divisors in the sequence is a prime.

The resulting graph is, in effect, a hyper-rectangular grid of dimension k without holes, so the number of vertices is the product of the size of each dimension, which is one more than the maximum exponent of that dimension (because the range of exponents also includes 0).

Every such grid has a Hamiltonian path. That's easy to prove, and the proof provides a good hint as to how to construct the sequence. The proof is by induction over k:

If k = 0 (in other words, n = 1), then the only possible path is [1], which is trivially valid.
For the induction step, we assume that there is a valid path for any graph with dimension i < k. We then construct a valid path for n by partitioning its prime factorization into two disjoint non-empty sets of primes, representing a pair of relatively prime integers q and r whose product is n. Evidently, the divisor graphs for q and r both have fewer than k dimensions, so by the induction hypothesis we can create valid sequences Q and R. Also note that every divisor of n is the product of some divisor of q and some divisor of r. We can now produce a valid sequence for n by producing the Cartesian product of Q and R by iterating over R and alternating between Q and it's reverse Q^R. The resulting sequence is:

Q₁ × R₁, Q₂ × R₁, Q₃ × R₁, … Q_|Q| × R₁, Q_|Q| × R₂, Q_|Q|−1 × R₂, Q_|Q|−2 × R₂, … Q₁ × R₂, …

In practice, we don't actually need to construct the graph; the valid sequence described above can be generated directly. How we choose to split the prime factorization of n is irrelevant to the above proof, and the end result will have the same size. While it's possible that a divide-and-conquer bisection will be more efficient for large k, I think the additional complexity just complicates code understanding, so in the sample code below, I just split off R as the last prime, p_k and Q as all the rest (p₁, p₂, … p_k−1.) Furthermore, I forgo recursion, although it might have been more elegant, to construct the sequence left to right, augmenting it as I discover prime factors of n. (The code takes advantage of the fact that reversing a sequence twice returns to the original sequence and that multiplying by a base by a list of powers of a prime can be done by repeatedly multiplying the previous sequence by the prime, so the powers are implicit. I hope that's clear in the code.)

#include <vector>

class DivisorSequence {
    public:
    using value_type = unsigned long;
    
    DivisorSequence() {}
    
    std::vector<value_type>
    divisor_sequence(value_type n) {
        // For simplicity, initialise to the trivial sequence [1]
        std::vector<value_type> divisors{1};
        
        if (n % 2 == 0) n = handle_prime(divisors, n / 2, 2);
        value_type p = 3;
        while (p * p <= n) {
            if (n % p == 0)
                n = handle_prime(divisors, n / p, p);
            p  = 2;
        }
        if (n != 1) handle_prime(divisors, 1, n);
        return divisors;
    }
    
    private:
    // Appends the alternately reversed sequence in divisors with
    // each successive power of p in n. The parameter n is implicit;
    // On entry, n must have been divisible by p, and next_n must be n / p.
    value_type
    handle_prime(std::vector<value_type>& divisors,
                 value_type next_n,
                 value_type p) {
        size_t prev_end = divisors.size();
        size_t prev_len = prev_end;
        // We know that there is at least one divisor.
        while(1) {
            for (size_t i = 1; i <= prev_len;   i)
                divisors.push_back(p * divisors[prev_end - i]);
            if (next_n % p != 0) break;
            next_n /= p;
            prev_end  = prev_len;
        }
        return next_n;
    }
};

There are lots of potential optimisations which could be added to that code, starting with special-casing small primes (particularly 2) and using a more efficient factorization algorithm, perhaps based on a sieve. But it's certainly fast enough for the problematic cases you list.