# LCM3GCD2 - Editorial

Author: deep2905
Tester: abhidot
Editorialist: iceknight1093

3034

# PREREQUISITES:

Familiarity with GCD and LCM

# PROBLEM:

There is a hidden array A of length N. Find it using upto 3N queries of the following types:

• Given i and j, you will receive \gcd(A_i, A_j)
• Given i, j, k, you will receive \text{lcm}(A_i, A_j, A_k)

# EXPLANATION:

With 3N operations, we need to use about 3 per index.

Note that constraints have N \geq 4, so we certainly need to be able to solve for N = 4. Letâ€™s attempt to do that first.

### Solving N = 4

There are 10 possible queries we can ask on these four elements â€” 6 pairwise GCDs, and 4 triplet LCMs. Letâ€™s ask them all.

Let g_{ij} = \gcd(A_i, A_j), and L_{ijk} = \text{lcm}(A_i, A_j, A_k).

Recall that for two numbers x and y, we have the relation \gcd(x, y)\cdot \text{lcm}(x, y) = x\cdot y.
However, we now have the lcm of three numbers, so weâ€™d like to generalize this a bit.
With a little understanding of how the above formula comes about, thatâ€™s not very hard.

Details

For three numbers x, y, z, we want some relation between \text{lcm}(x, y, z) and their gcds.
This can be done using the inclusion-exclusion principle.

This is an overestimation, since any factor common to two of them has been multiplied twice.
So, we can divide out \gcd(x, y), \gcd(x, z), and \gcd(y, z).

However, this is now an underestimation: any factor thatâ€™s common to all three has been taken out too many times, so we add it back in, i.e, multiply by \gcd(x, y, z).

This gives us

\text{lcm}(x, y, z) = \frac{x\cdot y\cdot z\cdot \gcd(x, y, z)}{\gcd(x, y)\gcd(x, z)\gcd(y, z)}

Note that this isnâ€™t quite a formal proof: why taking gcds is necessary and sufficient to cover the common factors takes some getting used to.
A more basic proof is to apply this argument to each prime factor of the LCM, then take the product of all the resulting prime powers.

You might notice that the formula we derived is pretty close to what we want.
For example, we have

L_{123} = \frac{A_1\cdot A_2\cdot A_3 \cdot \gcd(A_1, A_2, A_3)}{g_{12} \cdot g_{23} \cdot g_{13}}

We know most of the parts of this equation: everything except A_1, A_2, A_3 in fact (since \gcd(A_1, A_2, A_3) = \gcd(g_{12}, g_{13})).
In particular, this tells us the value of the product A_1\cdot A_2\cdot A_3. Let this value be P_{123}

Similarly, we can obtain P_{124}, P_{134}, P_{234}.

Multiplying all these together, we have

P_{123}\times P_{124}\times P_{134}\times P_{234} = (A_1\cdot A_2\cdot A_3\cdot A_4)^3

So, taking the cube root of their product gives us the value of A_1\cdot A_2\cdot A_3\cdot A_4; from which computing the individual values of A_1, A_2, A_3, A_4 is quite easy, since we also know the products of all their triplets.

This solves N = 4 entirely, so letâ€™s move on to a full solution.

### Solving N \gt 4

Letâ€™s first use the N = 4 solution to find the first 4 elements.

Suppose we want to find A_i, where i \gt 4.
Recall the lcm-gcd relation we had from earlier. In particular, applying it to indices (1, 2, i) gives us

L_{12i} = \frac{A_1\cdot A_2\cdot A_i \cdot \gcd(A_1, A_2, A_i)}{g_{12} \cdot g_{2i} \cdot g_{1i}}

However, we already know A_1 and A_2, so the only new pieces of information we require here are L_{12i}, g_{1i}, and g_{2i}; each of which require one query.
This also allows us to compute \gcd(A_1, A_2, A_i), after which computing A_i is easy since we know every other part of that expression.

This takes 3 queries for each index past the fourth, and 10 queries for the first 4, for a total of 3N - 2 queries.

### Implementation and large numbers

The above solution idea did not go into specifics about the intermediate values, but the observant reader may notice that we deal with some rather large numbers.
In particular, (A_1\cdot A_2\cdot A_3\cdot A_4) can be as large as 10^{72}, which definitely doesnâ€™t fit in the limit of long long.

There are a few ways to get around this.
By far the simplest is to simply use a language with support for large numbers, such as Python or Java.
In C++, the simplest way is probably to use a BigInteger class pulled from somewhere.

If youâ€™re not willing to do either of the above, there is still hope.
First off, notice that while the numbers we deal with can get quite large, they are products of numbers that are \leq 10^6.
This means that all their prime factors are \leq 10^6.

Since multiplication and division can be simulated easily on the prime factorizations of numbers, this gives us one way of dealing with the issue.

Once you know A_1\cdot A_2\cdot A_3, prime factorize it by iterating over all numbers from 1 to 10^6. Note that iterating till \sqrt{A_1\cdot A_2\cdot A_3} as usual might TLE.
Similarly, prime factorize the other three product triples, which allows you to compute the prime factorization of A_1 \cdot A_2\cdot A_3\cdot A_4. After this, each of the four elements can be computed from this prime factorization, and weâ€™re done.

Note that every other value does in fact fit inside 10^{18}, as long as you take care of the order of operations correctly (i.e, divide first before multiplying when possible).

### A small challenge

Itâ€™s possible to solve this problem in C++ without using a BigInteger class and without prime factorizations (which would allow you to solve for, say T \leq 10^4 and sum(N) \leq 10^4).

In particular, only the use of 128-bit integers is necessary, i.e, the __int128_t data type. Do you see a way to do this?

# TIME COMPLEXITY

\mathcal{O}(N) per test case.

# CODE:

Setter's code (C++)
/**
* Author : RDP
* There are no two words in the English language more harmful than "good job".
* 1729 ;)
**/
#include <bits/stdc++.h>

using namespace std;
using ll = long long;

/********** Definations, Macros and Debug Stuff  **********/
void debug_out() { cerr << '\n'; }
string to_string(const string &s) { return s; }
{
cerr << " " << to_string(H);
debug_out(T...);
}

#define endl '\n'
#define debug(...) cerr << "[" << #__VA_ARGS__ << "]: ", debug_out(__VA_ARGS__)
#define GODSPEED                 \
ios::sync_with_stdio(false); \
std::cin.tie(NULL);          \
std::cout.tie(NULL);
#define all(x) (x).begin(), (x).end()
const long double EPS = 5e-8;
#define PI 3.1415926535897932384626433832795
const ll MOD = 1000000007;
/**********************************************************/
/**************** Frequently used functions ***************/
template <typename T>
inline void print_vector(vector<T> &a)
{
for (auto &x : a)
cout << x << ' ';
cout << endl;
}
inline ll binary_pow(ll a, ll b)
{
ll res = 1;
while (b > 0)
{
if (b & 1)
res = res * a;
a = a * a;
b >>= 1;
}
return res;
}
inline ll mod_pow(ll x, ll y, ll m = MOD)
{
ll res = 1;
x = x % m;
if (x == 0)
return 0;
while (y > 0)
{
if (y & 1)
res = (res * x) % m;
y = y >> 1;
x = (x * x) % m;
}
return res;
}
inline ll mod_add(ll a, ll b, ll m = MOD)
{
a = a % m;
b = b % m;
return (((a + b) % m) + m) % m;
}
inline ll mod_mul(ll a, ll b, ll m = MOD)
{
a = a % m;
b = b % m;
return (((a * b) % m) + m) % m;
}
inline ll mod_sub(ll a, ll b, ll m = MOD)
{
a = a % m;
b = b % m;
return (((a - b) % m) + m) % m;
}
inline ll mminvprime(ll a, ll b)
{
return mod_pow(a, b - 2, b);
}
inline ll mod_div(ll a, ll b, ll m = MOD)
{
a = a % m;
b = b % m;
return (mod_mul(a, mminvprime(b, m), m) + m) % m;
}
inline ll ceilf(ll x, ll y)
{
return x % y == 0 ? x / y : x / y + 1;
}
// Use this for randomizing things

set<ll> primes;
vector<bool> is_prime;
void precompute_primes(ll n)
{
is_prime.resize(n + 1, 1);
for (ll p = 2; p * p <= n; p++)
{
if (is_prime[p] == true)
{
for (ll i = p * p; i <= n; i += p)
is_prime[i] = false;
}
}
for (ll p = 2; p <= n; p++)
if (is_prime[p])
primes.insert(p);
}

/**********************************************************/
map<pair<int, int>, ll> gcd_cache;
map<tuple<int, int, int>, ll> lcm_cache;
int numq = 0;

class FracInt
{
public:
map<ll, ll> up, down;
FracInt() { ; }
FracInt(__int128_t x)
{
for (auto p : primes)
{
if (p > x)
break;
while (x % p == 0)
x /= p, up[p]++;
up[1]++;
}
}
void mul(__int128_t x)
{
for (auto p : primes)
{
if (p > x)
break;
while (x % p == 0)
x /= p, up[p]++;
}
}
void div(__int128_t x)
{
for (auto p : primes)
{
if (p > x)
break;
while (x % p == 0)
x /= p, down[p]++;
}
normalize();
}
void mul(FracInt y)
{
for (auto p : y.up)
{
up[p.first] += p.second;
}
for (auto p : y.down)
{
down[p.first] += p.second;
}
normalize();
}
void div(FracInt y)
{
for (auto p : y.up)
{
down[p.first] += p.second;
}
for (auto p : y.down)
{
up[p.first] += p.second;
}
normalize();
}
void normalize()
{
for (auto &p : up)
{
p.second -= down[p.first];
down[p.first] = 0;
}
}
__int128_t get_abs()
{
normalize();
__int128_t val = 1;
for (auto &p : up)
{
for (int i = 0; i < p.second; i++)
val *= p.first;
}
return val;
}
};

FracInt cuberoot(FracInt x)
{
x.normalize();
for (auto &p : x.up)
p.second /= 3;
return x;
}

ll query(int i, int j, int k = -1)
{
i++;
j++;
if (k != -1)
k++;
if (k == -1)
{
if (gcd_cache.count({i, j}))
return gcd_cache[{i, j}];
}
else
{
if (lcm_cache.count({i, j, k}))
return lcm_cache[{i, j, k}];
}
numq++;
if (k == -1)
cout << "1 " << i << " " << j << endl;
else
cout << "2 " << i << " " << j << " " << k << endl;
ll x;
cin >> x;
if (k == -1)
{
gcd_cache[{i, j}] = x;
gcd_cache[{j, i}] = x;
}
else
{
vector<int> tmp = {i, j, k};
sort(all(tmp));
do
{
lcm_cache[{tmp[0], tmp[1], tmp[2]}] = x;
} while (next_permutation(all(tmp)));
}
assert(x != -1);
return x;
}
void terminate(vector<ll> a)
{
cout << "3 ";
print_vector(a);
ll x;
cin >> x;
assert(x == 1);
}
FracInt solve_for_mul(int i, int j, int k)
{
// let a,b,c
__int128_t lcmabc = (query(i, j, k));

__int128_t gcdab = (query(i, j));
__int128_t gcdac = (query(i, k));
__int128_t gcdbc = (query(j, k));
__int128_t gcdabc = __int128_t(gcd(ll(gcdab), ll(gcdbc)));
FracInt abc(1);
abc.mul(lcmabc);
abc.mul(gcdab);
abc.mul(gcdac);
abc.mul(gcdbc);
abc.div(gcdabc);

return abc;
}
vector<ll> solve_for4(int a, int b, int c, int d)
{
// let a,b,c,d be first 4 elements.
FracInt abc = solve_for_mul(a, b, c);
FracInt bcd = solve_for_mul(b, c, d);
FracInt acd = solve_for_mul(a, c, d);
FracInt abd = solve_for_mul(a, b, d);

FracInt abcd3(1);
abcd3.mul(abc);
abcd3.mul(bcd);
abcd3.mul(acd);
abcd3.mul(abd);

vector<ll> ans(4);
FracInt a0 = cuberoot(abcd3);
a0.div(bcd);

ans[0] = ll(a0.get_abs());

FracInt a1 = cuberoot(abcd3);
a1.div(acd);

ans[1] = ll(a1.get_abs());

FracInt a2 = cuberoot(abcd3);
a2.div(abd);

ans[2] = ll(a2.get_abs());

FracInt a3 = cuberoot(abcd3);
a3.div(abc);

ans[3] = ll(a3.get_abs());
return ans;
}
void solve_for_index(vector<ll> &A)
{
// we want c from a,b,c
__int128_t a = A[A.size() - 2];
__int128_t b = A[A.size() - 1];
int i = A.size() - 2, j = A.size() - 1, k = A.size();
__int128_t lcmabc = query(i, j, k);
__int128_t gcdab = query(i, j);
__int128_t gcdac = query(i, k);
__int128_t gcdbc = query(j, k);
__int128_t gcdabc = gcd(ll(gcdab), ll(gcdbc));

ll c = ll((lcmabc * gcdab * gcdac * gcdbc) / (a * b * gcdabc));
A.push_back(c);
return;
}

void test_case()
{
numq = 0;
ll t = 1e6;
gcd_cache.clear();
lcm_cache.clear();
int n;
cin >> n;
auto a = solve_for4(0, 1, 2, 3);
for (int i = 4; i < n; i++)
{
solve_for_index(a);
}
terminate(a);
return;
}

int main()
{
// GODSPEED;
precompute_primes(1e6 + 10);
int t = 1;
cin >> t;
for (int i = 1; i <= t; i++)
{
// cout << "Case #" << i << ": ";
test_case();
}
return 0;
}
/*
Some things to remember when you're stuck:
1. Check for edge cases.
2. Stay Calm.
3. Don't be stupid (search for silly mistakes).
4. Read problem again. Approach solution from different point of view.
5. In case of modulo, check for negative result (add MOD).

Some common C++ pit falls:
1. Don't use inbuilt ceil.
2. Never take inputs as double unless it is necessary.
3. Don't pass INT in accumulate.
*/

Tester's code (C++)


// Problem: Twice GCD Thrice LCM
// Contest: CodeChef - STR84TST
// Memory Limit: 256 MB
// Time Limit: 1500 ms
// Author: abhidot

// #pragma GCC optimize("O3,unroll-loops")
// #pragma GCC target("avx2,bmi,bmi2,lzcnt,popcnt")
#include <bits/stdc++.h>
#include <ext/pb_ds/assoc_container.hpp>
#include <ext/pb_ds/tree_policy.hpp>
#define int long long
#define ll long long
#define IOS std::ios::sync_with_stdio(false); cin.tie(NULL);cout.tie(NULL);
#define pb push_back
#define mod 1000000007
#define mod2 998244353
#define lld long double
#define pii pair<int, int>
#define ff first
#define ss second
#define all(x) (x).begin(), (x).end()
#define uniq(v) (v).erase(unique(all(v)),(v).end())
#define rep(i,x,y) for(int i=x; i<y; i++)
#define fill(a,b) memset(a, b, sizeof(a))
#define vi vector<int>
#define V vector
#define setbits(x) __builtin_popcountll(x)
#define w(x)  int x; cin>>x; while(x--)
using namespace std;
using namespace __gnu_pbds;
template <typename num_t> using ordered_set = tree<num_t, null_type, less<num_t>, rb_tree_tag, tree_order_statistics_node_update>;
const long long N=200005, INF=2000000000000000000, inf = 2e9+5;

int power(int a, int b, int p){
if(a==0)
return 0;
int res=1;
a%=p;
while(b>0)
{
if(b&1)
res=(res*a)%p;
b>>=1;
a=(a*a)%p;
}
return res;
}

void print(bool n){
if(n){
cout<<"YES\n";
}else{
cout<<"NO\n";
}
}

struct custom_hash {
static uint64_t splitmix64(uint64_t x) {
// http://xorshift.di.unimi.it/splitmix64.c
x += 0x9e3779b97f4a7c15;
x = (x ^ (x >> 30)) * 0xbf58476d1ce4e5b9;
x = (x ^ (x >> 27)) * 0x94d049bb133111eb;
return x ^ (x >> 31);
}

size_t operator()(uint64_t x) const {
static const uint64_t FIXED_RANDOM = chrono::steady_clock::now().time_since_epoch().count();
return splitmix64(x + FIXED_RANDOM);
}
};

int query1(int a, int b){
cout<<"1 "<<a<<" "<<b<<endl;
int r; cin>>r;
return r;
}

int query2(int a, int b, int c){
cout<<"2 "<<a<<" "<<b<<" "<<c<<endl;
int r; cin>>r;
return r;
}

map<int, int> factorize(int x){
map<int, int> f;
for(int i=2;i<=x;i++){
while(x%i==0){
f[i]++;
x/=i;
}
}
return f;
}

void merge(auto& a, auto b){
for(auto i: b){
a[i.ff]+=i.ss;
}
}

int32_t main()
{
IOS;
w(T){
int n; cin>>n;
int gab = query1(1, 2);
int gbc = query1(2, 3);
int gcd = query1(3, 4);
int gac = query1(1, 3);
int gbd = query1(2, 4);

int labc = query2(1, 2, 3);
int labd = query2(1, 2, 4);
int lacd = query2(1, 3, 4);
int lbcd = query2(2, 3, 4);

int gabc = __gcd(gab, gbc);
int gabd = __gcd(gab, gbd);
int gacd = __gcd(gac, gcd);
int gbcd = __gcd(gbc, gcd);

int abc = (labc/gabc)*gab*gbc*gac;
int bcd = (lbcd/gbcd)*gbc*gcd*gbd;

map<int, int> p = factorize(abc);
map<int, int> q = factorize(abd);
map<int, int> r = factorize(acd);
map<int, int> s = factorize(bcd);
map<int, int> t = p;
merge(t, q);
merge(t, r);
merge(t, s);
// cout<<abc<<" "<<abd<<" "<<acd<<" "<<bcd<<endl;
int d = 1;
for(auto i: t){
int x = i.ss/3;
x-=p[i.ff];
d*=power(i.ff, x, INF);
}
int c = 1;
for(auto i: t){
int x = i.ss/3;
x-=q[i.ff];
c*=power(i.ff, x, INF);
}
int b = 1;
for(auto i: t){
int x = i.ss/3;
x-=r[i.ff];
b*=power(i.ff, x, INF);
}
int a = 1;
for(auto i: t){
int x = i.ss/3;
x-=s[i.ff];
a*=power(i.ff, x, INF);
}

int v[n+1]={0};
v[1]=a, v[2]=b, v[3]=c, v[4]=d;
for(int i=5;i<=n;i++){
int gax = query1(1, i);
int gbx = query1(2, i);
int labx = query2(1, 2, i);
int gabx = __gcd(gax, gbx);
int abx = (labx/gabx)*gab*gax*gbx;
v[i]=abx/(a*b);
}
cout<<3<<" ";
for(int i=1;i<=n;i++) cout<<v[i]<<" ";
cout<<endl;
int rr; cin>>rr;
}
}

Editorialist's code (Python)
import math

print(1, i, j)
x = int(input())
assert x != -1
return x
print(2, i, j, k)
x = int(input())
assert x != -1
return x

def cbrt(x):
lo, hi = 1, 10**24 + 2
while lo < hi:
mid = (lo + hi)//2
if mid**3 < x: lo = mid+1
else: hi = mid
return lo

for _ in range(int(input())):
n = int(input())
a = [0]*n

# Find a[0], a[1], a[2], a[3]

p123 = l123 * g12 * g23 * g13 // math.gcd(g12, g23)
p124 = l124 * g12 * g24 * g14 // math.gcd(g12, g14)
p134 = l134 * g13 * g34 * g14 // math.gcd(g13, g14)
p234 = l234 * g23 * g34 * g24 // math.gcd(g23, g24)

allprod3 = p123 * p124 * p134 * p234
allprod = cbrt(allprod3)

a[0] = allprod // p234
a[1] = allprod // p134
a[2] = allprod // p124
a[3] = allprod // p123

for i in range(4, n):