算法题整理——数论

发表于 | 分类于 |

素数筛

Eratosthenes筛

找出小于n的所有素数，方法是这样的：

首先，从2开始遍历每个小于n的数i，对于i的倍数，我们知道肯定不是素数，所以标记它为合数。

然后优化，假设j<i，那么i的j倍一定在遍历j的倍数的时候就已经标记过了，那么我们在遍历i的倍数的时候，直接从i的i倍开始。

一直遍历到i的j倍不再小于等于n结束。

时间复杂度为O(nloglogn)，结论怎么来的笔者不知道，猜测和自然数中素数的密度有关。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
#include <cstdio>
#include <cstring>
using namespace std;
const int N = 1e7 + 5;
bool vis[N];
int main(){
	memset(vis, 0, sizeof(vis));
	int n;
	scanf("%d", &n);
	for(int i = 2; i <= n; ++i){
		if(vis[i])continue;
		printf("%d is a prime.\n", i);
		for(int j = i; j <= n / i; ++j) vis[i * j] = 1;
	}
	return 0;
}

线性筛

Eratosthenes筛仍然会标记重复的合数，如标记12时，会在扫描2，扫描3的时候各标记一遍。

为此，再次优化算法，得到了时间复杂度为O(N)的线性筛。

1、以此考虑2~N的每个数i

2、若vi=0，说明i是质数，将其保留下来

3、扫描不大于vi的每个质数p，令vi×p=p

这样，每个合数只会被其最小质因子筛一次。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
#include <cstdio>

using namespace std;

int v[10000];
int p[10000];
int m;

int main(){
	int n;
	scanf("%d", &n);
	m = 0;
	for(int i = 2; i <= n; ++i){
		if(!v[i]){
			v[i] = i;
			p[m++] = i;
		}
		for(int j = 0; j < m; ++j){
			if(v[i] < p[j]) break;
			if(i * p[j] > n) break;
			v[i * p[j]] = p[j];
		}
	}
	
	for(int i = 0; i < m; ++i){
		printf("%d\n", p[i]);
	}
	return 0;
}

CH3101_阶乘分解

描述

给定整数 N(1≤N≤10^6)，试把阶乘N!分解质因数，按照算术基本定理的形式输出分解结果中的p_i和c_i即可。

思路

思路1（自己想的）：计算出1~1000的所有质数，然后设1~1000000这么大的数组，下标和值相等，用1~100000挨个除这些质数并统计个数。然后检查得到的商，这时商如果不是1，必定是大于1000的质数，所以再排序，统计1000以上的质数个数。复杂度O(N^1.5)

思路2（蓝书）：因为N！中含有某个质因子的数目就是1~N这些数中含有这个质因子的数目之和。因此，N/p下取整就是所有至少含一个p的整数数目，N/p/p再下取整就是至少含两个p的整数数目。要计算所有p的数目，应计算N/p+N/p/p+……。因为1~N连乘必然含有小于等于N的所有素数，所以要求出1~N的所有素数才行。复杂度O(NlogN)

因为题目数据较弱，我的思路也通过了。但由于此题对于阶乘分解的方法具有普遍意义，我按蓝书上的思路又重写了代码。

解决

思路1

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
#include <cstdio>
#include <cstring>
#include <algorithm>

using namespace std;

const int N = 1e3 + 1;

int num[N * N];

int p[N];

int n;

int con;

void init(){
	con = 0;
	int v[N];
	memset(v, 0, sizeof(v));
	for(int i = 2; i < N; ++i){
		if(v[i] == 0){
			v[i] = i;
			p[con++] = i;
		}
		for(int j = 0; j < con; ++j){
			if(i * p[j] >= N || p[j] > v[i]) break;
			v[i * p[j]] = v[i];
		}
	}
}

int main(){
	init();
	scanf("%d", &n);
	for(int i = 1; i <= n; ++i){
		num[i] = i;
	}
	for(int i = 0; i < con; ++i){
		int c = 0;
		for(int j = 1; j <= n; ++j){
			while(num[j] % p[i] == 0){
				num[j] /= p[i];
				++c;
			}
		}
		if(c){
			printf("%d %d\n", p[i], c);
		}
	}
	sort(num + 1, num + n + 1);
	num[n + 1] = 1000009;
	int c = 1;
	for(int i = 2; i <= n + 1; ++i){
		if(num[i - 1] == 1){
			c = 1;
		}else{
			if(num[i - 1] == num[i]) ++c;
			else{
				printf("%d %d\n", num[i - 1], c);
				c = 1;
			}
		}
	}
	return 0;
}

思路2

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
#include <cstdio>
#include <cstring>

using namespace std;

const int N = 1e6 + 5;

int n;

int con = 0;

int p[N];
int v[N];

void init(){	
	memset(v, 0, sizeof(v));
	for(int i = 2; i <= n; ++i){
		if(!v[i]){
			v[i] = i;
			p[con++] = i;
		}
		for(int j = 0; j < con; ++j){
			if(p[j] * i > n || p[j] > v[i]) break;
			v[p[j] * i] = v[i];
		}
	}
}

int main(){
	scanf("%d", &n);
	init();
	for(int i = 0; i < con; ++i){
		int c = 0;
		int t = n;
		int b = p[i];
		while(t >= b){
			c += t / b;
			t /= b;
		}
		printf("%d %d\n", b, c);
	}
	return 0;
}

POJ2689_Prime Distance

描述

给定两个数l和u，要求找出l到u之间距离最小和最大的相邻素数。

范围

l和u都是int型正整数，u-l小于等于1e6

思路

由于u-l比较小，所以可以筛出根号u以内的所有素数，用它们挨个检验u到l之间的数。因为任何合数都会至少有一个小于等于它算术平方根的因子。所以这样检测不到的一定是素数。

解决

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
#include <cstdio>
#include <cstring>
#include <cmath>

using namespace std;

const long long N = 1e6 + 5;
const long long M = 50000;

long long vis[N];
long long v[M];
long long p[M];

long long l, u, m = 0;

void init(){
	memset(v, 0, sizeof(v));
	for(long long i = 2; i < M; ++i){
		if(!v[i]){
			v[i] = i;
			p[m++] = i;
		}
		for(int j = 0; j < m; ++j){
			if(i * p[j] >= M || v[i] < p[j]) break;
			v[i * p[j]] = p[j];
		}
	}
}

int main(){
	//freopen("input.txt", "r", stdin);
	//freopen("myout.txt", "w", stdout);
	init();
	while(scanf("%lld %lld", &l, &u) != EOF){
		if(l == 1) ++l;
		memset(vis, 0, sizeof(vis));
		for(long long i = 0; i < m; ++i){
			long long time = (l - 1) / p[i] + 1;
			if(time == 1) ++time;
			//printf("%lld %lld\n", p[i], time);
			//if(time < p[i]) break;
			for(; time * p[i] <= u; ++time){
				vis[time * p[i] - l] = 1;
			}
		}
		/*
		for(int i = l; i <= u; ++i){
			printf("%lld %lld\n", vis[i - l], i);
		}
		*/
		long long maxlast, maxcur, maxdlt;
		long long minlast, mincur, mindlt;
		long long flag = 0;
		long long lastp = -1, curp = -1, dlt;
		for(long long i = l; i <= u; ++i){
			if(vis[i - l]) continue;
			else{
				if(lastp == -1){
					lastp = i;
				}else if(curp == -1){
					curp = i;
					dlt = curp - lastp;
					maxlast = minlast = lastp;
					maxcur = mincur = curp;
					maxdlt = mindlt = dlt;
					flag = 1;
				}else{
					lastp = curp;
					curp = i;
					dlt = curp - lastp;
					if(maxdlt < dlt){
						maxlast = lastp;
						maxcur = curp;
						maxdlt = dlt;
					}
					if(mindlt > dlt){
						minlast = lastp;
						mincur = curp;
						mindlt = dlt;
					}
				}
			}
		}
		if(flag){
			printf("%lld,%lld are closest, %lld,%lld are most distant.\n", minlast, mincur, maxlast, maxcur);
		}else{
			printf("There are no adjacent primes.\n");
		}
	}
	return 0;
}

反思和坑

计算的中间结果会超过int的最大范围，最好的办法是全部使用long long。

注意l、u的取值范围包括1，而1是不算素数的。

最后发现一直错误的原因是自作聪明的剪枝：

检查l到u之间的数有没有质因子时，没有把2~50000中所有的素数都检查一遍。

约数

试除法

求某个数的全部约数，只需要扫描1~N^(1/2)中所有数，尝试该数能否整除N即可。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
#include <cstdio>

using namespace std;

int rets[1000];
int con = 0;
int n;

int main(){
	int i;
	scanf("%d", &n);
	for(i = 1; i * i < n; ++i){
		if(n % i == 0){
			rets[con++] = i;
			rets[con++] = n / i;	
		}
	}
	if(n % i == 0){
		rets[con++] = i;
	}
	for(i = 0; i < con; ++i){
		printf("%d\n", rets[i]);
	}
	return 0;
}

试除法的推论

任何正整数N的因数个数不会超过2N^(1/2)

倍数法

求1~N的全部正整数各自的全部因数。可以考虑对1~N的所有正整数遍历其倍数，复杂度为NlogN

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
#include <cstdio>
#include <vector> 
using namespace std;

vector<int> rets[1000];
int con = 0;
int n;

int main(){
	int i;
	scanf("%d", &n);
	for(i = 1; i <= n; ++i){
		for(int j = 1; i * j <= n; ++j){
			rets[i * j].push_back(i);
		}
	}
	for(i = 0; i <= n; ++i){
		for(vector<int>::iterator j = rets[i].begin(); j != rets[i].end(); ++j){
			printf("%d\n", *j);
		}
		printf("\n");
	}
	return 0;
}

倍数法的推论

1~N的全部正整数的因子个数和大约为NlogN

BZOJ1053_反素数

描述

设一个数x的约数个数为h(x)，那么如果任取y < x，有h(y) < h(x)，称x为反素数。

1，2，4，6都是反素数

现给出N，求N以内（包括N）的最大反素数

范围

N <= 2e9

思路

这个题应该是我做不出来的，因为涉及到了一些推导

首先，如果两个数因数个数相同，那么反素数只可能是小的那个

其次，比2e9小的数质因子不会超过9个，这是因为前10个质数的乘积为6,469,693,230，已经大于2e9了

最后，N以内的反素数必定表示为(2^p1)×(3^p2)×(5^p3)×(7^p4)×(11^p5)×(13^p6)×(17^p7)×(19^p8)×(23^p9)×(29^p10)，这其中必定有p1到p10递减

经由以上推导，我们就可以使用一个简单的dfs来确定p1到p10的具体数目

这里还用到了一个东西，如果x的算术基本定理分解式为：

x = (p1^a1)×(p2^a2)×……×(pn^an)

那么x的因数个数是(a1+1)×(a2+1)×……×(an+1)，这结论稍微一动脑子就能想出来

这可以用在dfs里

解决

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
#include <cstdio>

using namespace std;

long long n;
long long p[] = {2, 3, 5, 7, 11, 13, 17, 19, 23};

long long max_num = 0;
long long max_ret = 0;

void dfs(int layer, int max_time, long long num, long long ret){
	//printf("%d %d %lld %lld\n", layer, max_time, num, ret);
	
	if(layer == 9){
		if(max_ret < ret){
			max_ret = ret;
			max_num = num;
		}else if(max_ret == ret && max_num > num){
			max_num = num;
		}
		return;
	}else{
		dfs(layer + 1, max_time, num, ret);
	}
	
	long long new_num = num;
	long long new_ret;
	
	for(int i = 1; i <= max_time; ++i){
		new_num *= p[layer];
		if(new_num > n) break;
		new_ret = ret * (long long)(i + 1);
		dfs(layer + 1, i, new_num, new_ret);
	}
}

int main(){
	scanf("%lld", &n);
	dfs(0, 0x7fffffff, 1, 1);
	printf("%lld\n", max_num);
	return 0;
}

BZOJ1257_余数之和

描述

给出正整数n和k，计算j(n, k)=k mod 1 + k mod 2 + k mod 3 + … + k mod n的值，其中k mod i表示k除以i的余数。
例如j(5, 3)=3 mod 1 + 3 mod 2 + 3 mod 3 + 3 mod 4 + 3 mod 5=0+1+0+3+3=7

范围

1 <= n, k <= 1e9

思路

肯定不能直接算，但怎么简化呢？

可以考虑n > k时的情况，此时余数必为k，总和为k × (n - k)，接下来集中解决n <= k时的情况

令k除以n，考虑：

1
2
3
n = k		余数是0， 商是1
n = k - 1	余数是1， 商是1
n = k - 2	余数是2， 商是1

余数逐个增大1，这种情况在什么时候停止呢？在n小于k的一半的时候，再之后，余数又会逐个增大2，以此类推

1
2
n = k / 2	  余数是0，商是2
n = k / 2 - 1 余数是2，商是2

可见，余数构成等差数列，直到商重置时为止

这样，我们可以计算商为1时的等差数列余数和，然后商为2，一直到商为k

但这样又导致复杂度变成了O(k)，还是太大

我们观察商比较大时的情况，举个例子，假设k = 36，商为9，此时余数可能有0一种，然后商为10，没有可能，商为11，又没有可能

可见，商比较大时，可能有余数的区间更为稀疏

因此，我们只计算商为1到1000的区间，然后剩下的直接用取余式求出k除以1到k除以(k/1000)的余数即可

解决

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
#include <cstdio>

using namespace std;

int main(){
	long long n, k;
	scanf("%lld %lld", &n, &k);
	long long ret = 0;
	if(n > k) ret = k * (n - k), n = k;
	if(k > 1000){
		for(long long div = k / n; div <= 1000; ++div){ 
			long long end = k / div;
			long long start = k / (div + 1) + 1;
			
			if(k - div * end >= end || start > end) continue;
			
			//printf("%lld %lld %lld\n", div, end, start);
			end = k - div * end;
			start = k - div * start;
			ret += (start + end) * ((start - end) / div + 1) / 2;
			//printf("%lld %lld %lld\n", ret, end, start);
		}
		for(long long i = 2; i <= k / 1001; ++i){
			ret += (k % i);
		}
	}else{
		for(long long i = 1; i <= n; ++i){
			ret += (k % i);
		}
	}
	printf("%lld\n", ret);
	return 0;
}

欧拉函数

POJ3090_Visible Lattice Points

描述

从平面直角坐标系上，如果从某个点A出发到另一个点B连成的直线中间没有第三个点，那么可以说在A点能看到B点

现在，给出x和y坐标为在[0，N]之间整数的点，求从原点看去，有多少个可见的点

数据范围

1 < N < 1001

思路

首先，我们知道有三条特殊的线，x轴，y轴和y=x，这三条线上所有点都会被最前面的挡住

然后，如何判断某个点是否被前面的点挡住了呢？很简单，如果这个点的x和y坐标互质，就一定不会有点挡在这个点的前面

问题转换成了求多少个x和y坐标互质的点

然后，我们可以固定x坐标为p，那么显然这个数字就是p的欧拉函数phi(p)

因此答案是∑phi(n) + 3

解决

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
#include <cstdio>
#include <cstring>

using namespace std;

const int maxN = 1005;

int N, n;
int p[maxN];
int v[maxN];
int ret[maxN];
int ans[maxN];
int cnt;

void prime(){
	memset(v, -1, sizeof(v));
	for(int i = 2; i < maxN; i++){
		if(v[i] == -1){
			v[i] = i;
			p[cnt++] = i;
			
		}
		for(int j = 0; p[j] <= v[i] && j < cnt && i * p[j] < maxN; ++j){
			v[i * p[j]] = p[j];
		}
	}
	ret[1] = 1;
	for(int i = 2; i < maxN; ++i){
		ret[i] = ((i / v[i]) % v[i] ? v[i] - 1 : v[i]) * ret[i / v[i]];
	}
	ret[1] = 3;
	for(int i = 2; i < maxN; ++i){
		ret[i] <<= 1;
		ret[i] += ret[i - 1];
	}
}


int main(){
	prime();
	scanf("%d", &N);
	for(int con = 0; con < N; ++con){
		scanf("%d", &n);
		printf("%d %d %d\n", con + 1, n, ret[n]);
	}
	return 0;
}

POJ3696_The Luckiest number

描述

给出一个数L，寻找最短的由8组成的数字，使该数字是L的倍数

输出这个数字有几个8

数据范围

L < 2e9

思路

可以把由8构成的数看作8(10x-1)/9。问题转换成了求最小的x值

形如10x-1的数显然是不含有因子2和5的

所以8(10x-1)若能被L整除，L只能含有不超过3个因子2，并且不含有因子5

我们去掉L中的因子2后，问题转为10x-1能被9L整除

此时L与10互质，根据欧拉定理，我们直到x是phi(9L)的因数

那么，我们只需要列举因数，利用快速幂找出最小的符合条件的数即可

解决

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
#include <cstdio>
#include <cstring>


using namespace std;

const int maxN = 2e6 + 5;

long long p[maxN];
int v[maxN];
int cnt;

void prime(){
	cnt = 0;
	memset(v, -1, sizeof(v));
	for(int i = 2; i < maxN; ++i){
		if(v[i] == -1) v[i] = i, p[cnt++] = i;
		for(int j = 0; j < cnt && p[j] <= v[i] && i * p[j] < maxN; ++j){
			v[i * p[j]] = p[j];
		}
	}
}

long long f(long long n){
	long long ret = 1;
	for(int i = 0; i < cnt && p[i] <= n; ++i){
		if(n % p[i]) continue;
		ret *= (p[i] - 1);
		n /= p[i];
		while(n % p[i] == 0) n /= p[i], ret *= p[i];
		if(n == 1) break;
	}
	if(n != 1) ret *= (n - 1);
	return ret;
}

bool judje(long long i, long long n){
	long long x = 1;
	long long base = 10;
	while(i){
		if(i & 1) x *= base, x %= n;
		base *= base, base %= n;
		i >>= 1;
	}
	return x % n == 1;
}

int main(){
	long long n;
	prime();
	int id = 0;
	while(true){
		scanf("%lld", &n);
		if(!n) break;
		++id;
		printf("Case %d: ", id);
		if(n % 16 == 0 || n % 5 == 0) printf("0\n");
		else{
			while(n % 2 == 0) n /= 2;
			long long phi = f(n * 9);
			long long ret = phi;
			long long con = 1;
			while(n % 3 == 0) n /= 3, con *= 3;
			for(long long i = 1; i * i <= phi; ++i){
				if(ret < i) break;
				if(phi % i) continue;
				if(i % con) continue;
				if(judje(i, n)) ret = i;
				if(ret < phi / i || phi / i % con) continue;
				if(judje(phi / i, n)) ret = phi / i;
			}
			printf("%lld\n", ret);
		}
	}
	return 0;
}

欧几里得算法和扩展

CH3301_同余方程

描述

求关于 x的同余方程 ax ≡ 1(mod b) 的最小正整数解。

数据范围

对于 40% 的数据 2 ≤b≤1,000
对于 60% 的数据 2 ≤b≤50,000,000
对于 100%的数据 2 ≤a, b≤2,000,000,000

思路

先用拓展欧几里得算法求出通解

然后对b取模

解决

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
#include <cstdio>

using namespace std;

void exgcd(long long a, long long b, long long c, long long &x, long long &y){
	if(a == 1){
		x = c;
		y = 0;
	}else{
		long long x1, y1;
		exgcd(b, a % b, c, x1, y1);
		x = y1;
		y = x1 - y1 * (a / b);
	}
}

int main(){
	long long a, b, x, y, f = 0;
	scanf("%lld %lld", &a, &b);
	if(a < b){
		exgcd(b, a, 1, x, y);
		f = 1;
	}else exgcd(a, b, 1, x, y);
	if(f){
		printf("%lld\n", (y % b + b) % b);
	}
	else printf("%lld\n", (x % b + b) % b);
	return 0;
}

CH3201_Hankson的趣味题

描述

已知正整数a0,a1,b0,b1，设某未知正整数x满足：

1
2
1、  x和a0的最大公约数是a1；
2、  x和b0的最小公倍数是b1。

Hankson的“逆问题”就是求出满足条件的正整数x。但稍加思索之后，他发现这样的x并不唯一，甚至可能不存在。因此他转而开始考虑如何求解满足条件的x的个数。请你帮助他编程求解这个问题。

保证a0是a1的倍数，b1是b0的倍数

输入格式和范围

先输入一个整数n表示问题的个数

然后有n行，每行四个数a0，a1，b0，b1

n不大于2000，每次的四个整数不大于2e9

思路

把a0，a1，b0，b1都分解质因数

从题目条件给的倍数关系，可以知道，凡是a1有的质因数，a0一定也有，并且数目只多不少

举例说明：假设a1分解出了两个2，那么a0作为a1的倍数，也至少能分解出两个2

此外，作为a1和x的最大公因数，a0含有的质因子个数为a1和x含有质因子个数的最小值

举例说明：假设a1含有两个2，x含有一个2，那么a0作为最大公因数，一定含有一个2

我们遍历a0的质因数p，假设a0含有m0个p，a1含有m1个p，那么a1中的含p情况有两种：

1
2
3
1.m1 = m0，说明a0是含p较少的那个，所以x至少也要含有m0个p

2.m1 < m0，说明x是含p较少的那个，所以x必须含有m1个p

这样，我们可以把x的含各种质因数的情况列出来

然后，对于b0和b1，有类似的结论

为了加速，我们把2e5以内的质数先筛出来

解决

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
#include <cstdio>
#include <cstring>
#include <map>

using namespace std;

const int maxn = 5e4;
int v[maxn];
int p[maxn];
int cnt;

void initation(){
	cnt = 0;
	memset(v, 0, sizeof(v));
	for(int i = 2; i < maxn; ++i){
		if(!v[i]) p[cnt++] = v[i] = i;
		
		for(int j = 0; j < cnt; ++j){
			if(p[j] > v[i] || i * p[j] >= maxn) break;
			
			v[i * p[j]] = p[j];
		}
	}
}

void resolve(int x, map<int, int> &mpx){
	for(int i = 0; i < cnt; ++i){
		if(x == 1) break;
		else if(x % p[i]) continue;
			
		while(x % p[i] == 0){
			++mpx[p[i]];
			x /= p[i];
		}
	}
	
	if(x != 1) mpx[x] = 1;
}

int main(){
	int n, a0, a1, b0, b1;
	scanf("%d", &n);
	initation();
	while(n--){
		scanf("%d %d %d %d", &a0, &a1, &b0, &b1);
		
		map<int, int> mpa0;
		map<int, int> mpa1;
		map<int, int> mpb0;
		map<int, int> mpb1;
		
		resolve(a0, mpa0);
		resolve(a1, mpa1);
		resolve(b0, mpb0);
		resolve(b1, mpb1);
		/*
		for(map<int, int>::iterator i = mpa0.begin(); i != mpa0.end(); ++i){
			printf("%d %d\n", i->first, i->second);
		}
		*/
		map<int, pair<int, int> > mpx;
		
		for(map<int, int>::iterator i = mpa0.begin(); i != mpa0.end(); ++i){
			int key = i->first;
			int value = i->second;
			if(mpa1[key] == value){
				mpx[key].first = -1;
				mpx[key].second = value;
			}else{
				mpx[key].second = mpx[key].first = mpa1[key];
			}
		}
		/*
		for(map<int, pair<int, int> >::iterator i = mpx.begin(); i != mpx.end(); ++i){
			printf("%d %d %d\n", i->first, i->second.first, i->second.second);
		}
		*/
		int flag = 1;
		for(map<int, int>::iterator i = mpb1.begin(); i != mpb1.end(); ++i){
			int key = i->first;
			int value = i->second;
			if(mpa0[key] == 0){
				if(mpb0[key] == value){
					mpx[key].first = value;
					mpx[key].second = 0;
				}else{
					mpx[key].first = mpx[key].second = value;
				}
			}else{
				if(mpb0[key] == value){
					if(mpx[key].second > value){
						flag = 0;
						break;
					}else if(mpx[key].first == -1){
						mpx[key].first = value;
					}
				}else{
					if(mpx[key].second > value){
						flag = 0;
						break;
					}else if(mpx[key].first == -1){
						mpx[key].first = mpx[key].second = value;
					}else if(mpx[key].first != value){
						flag = 0;
						break;
					}
				}
			}
		}
		/*
		for(map<int, pair<int, int> >::iterator i = mpx.begin(); i != mpx.end(); ++i){
			printf("%d %d %d\n", i->first, i->second.first, i->second.second);
		}
		*/
		if(flag){
			long long ret = 1;
			for(map<int, pair<int, int> >::iterator i = mpx.begin(); i != mpx.end(); ++i){
				ret *= (long long)(i->second.first - i->second.second + 1);
			}
			printf("%lld\n", ret);
		}else{
			printf("0\n");
		}
	}
	return 0;
}