CCPC-Wannafly Winter Camp Day3 (Div2, onsite)
B. 集合
Description:
最近放假了,可怜计划和她的好朋友 Sylvia 去公园玩。
我们可以把地图抽象成一个二维平面,可怜最开始的位置是点 x ,Sylvia 最开始的位置是点 y,公园是一个以点 o 为圆心,半径为 r 的圆 O 。
因为进入公园要买票,因此她们约好了先在不进入公园的情况下,在公园的边界上集合,然后在一起买票进公园。她们想要选择一个最优的集合地点使得双方要走的距离长度和最小。
换句话说,令 d(s,t)d(s,t) 为在平面上,不经过圆 O 的内部,从 s 走到 t 的最短距离,可怜希望在圆 O 的边界上找到点 m,使得 d(x,m)+d(y,m) 尽可能地小。
Input:
输入第一行一个整数 t(1≤t≤106) 表示数据组数。
对于每组数据,输入一行 7 个整数,分别是 x 的横纵坐标, y 的横纵坐标, o 的横纵坐标和半径 r(1≤r≤103) ,其中所有坐标的绝对值不超过 103。
数据保证 x 和 y 都不在圆 O 的内部,即 x,y 和 o 的欧几里得距离都大于等于 r 。
Output:
对于每组数据,输出一行一个整数,表示最小的总距离,保留三位小数。
输入保证每组数据答案的第四位小数都不是 4 或者 5 。
Sample Input:
3
0 1 2 0 1 0 1
0 0 2 0 3 0 1
0 0 2 0 1 3 1
Sample Output:
2.571
2.000
4.472
题目链接
求圆外两点到圆上一点的最短路。
有两种情况:
- 当圆外两点不能以直线到达圆上同一点时两点做圆的切线,此时两点分别到切点加上切点中间弧线的距离是最短距离
- 当圆外两点能以直线到达圆上同一点时在圆上找到一点使得此点与圆外两点连线的线段互相为关于圆的反射线,此时两段线段的距离即为最短距离
圆上反射点的查找可以在圆外两点与圆心连线两线段与圆的两交点之间进行二分查找
AC代码:
#include <bits/stdc++.h>
using namespace std;
typedef long double db;
const db eps = 1e-8;
int Sgn(db Key) {
if (fabs(Key) < eps) return 0;
return Key < 0 ? -1 : 1;
}
struct Point {db X, Y;};
typedef Point Vector;
bool operator == (Point Key1, Point Key2) {return Sgn(Key1.X - Key2.X) == 0 && Sgn(Key1.Y - Key2.Y) == 0;}
Vector operator + (Vector Key1, Vector Key2) {return (Vector){Key1.X + Key2.X, Key1.Y + Key2.Y};}
Vector operator - (Vector Key1, Vector Key2) {return (Vector){Key1.X - Key2.X, Key1.Y - Key2.Y};}
db operator * (Vector Key1, Vector Key2) {return Key1.X * Key2.X + Key1.Y * Key2.Y;}
db operator ^ (Vector Key1, Vector Key2) {return Key1.X * Key2.Y - Key1.Y * Key2.X;}
Vector operator * (Vector Key1, db Key2) {return (Vector){Key1.X * Key2, Key1.Y * Key2};}
Vector operator / (Vector Key1, db Key2) {return (Vector){Key1.X / Key2, Key1.Y / Key2};}
Point Rotate(Point Key, db Angle) {
return (Point){Key.X * cos(Angle) - Key.Y * sin(Angle),
Key.X * sin(Angle) + Key.Y * cos(Angle)};
}
db DisPointToPoint(Point Key1, Point Key2) {return sqrt((Key1 - Key2) * (Key1 - Key2));}
db GetAngle(Vector Key1, Vector Key2) {return fabs(atan2(Key1 ^ Key2, Key1 * Key2));}
struct PolarPoint {db Angle, Length;};
Point PolarPointToPoint(PolarPoint Key) {
return (Point){Key.Length * cos(Key.Angle),
Key.Length * sin(Key.Angle)};
}
struct Line {Point S, T;};
typedef Line Segment;
db Length(Segment Key) {return DisPointToPoint(Key.S, Key.T);}
db DisPointToLine(Point Key1, Line Key2) {return fabs((Key1 - Key2.S) ^ (Key2.T - Key2.S)) / Length(Key2);}
db DisPointToSegment(Point Key1, Segment Key2) {
if (Sgn((Key1 - Key2.S) * (Key2.T - Key2.S)) < 0 || Sgn((Key1 - Key2.T) * (Key2.S - Key2.T)) < 0) {
return min(DisPointToPoint(Key1, Key2.S), DisPointToPoint(Key1, Key2.T));
}
return DisPointToLine(Key1, Key2);
}
struct Circle {Point Center;db Radius;};
int T;
Point Start, End;
Circle Park;
db DisOS, DisOE;
db Angle;
db Ans;
int main(int argc, char *argv[]) {
scanf("%d",&T);
for (int Case = 1; Case <= T; ++Case) {
scanf("%Lf%Lf", &Start.X, &Start.Y);
scanf("%Lf%Lf", &End.X, &End.Y);
scanf("%Lf%Lf%Lf", &Park.Center.X, &Park.Center.Y, &Park.Radius);
Ans = 0.0;
if (Sgn(DisPointToSegment(Park.Center, (Segment){Start, End}) - Park.Radius) >= 0) {
if (Sgn(DisPointToPoint(Park.Center, Start) - Park.Radius) == 0 ||
Sgn(DisPointToPoint(Park.Center, End) - Park.Radius) == 0) {
printf("%.3Lf\n", DisPointToPoint(Start, End));
continue;
}
Angle = GetAngle((Vector){Start - Park.Center}, (Vector){End - Park.Center});
if (Sgn(Angle) == 0) {
printf("%.3Lf\n", DisPointToPoint(Start, Park.Center) + DisPointToPoint(End, Park.Center) - 2 * Park.Radius);
continue;
}
Start = Start - Park.Center; End = End - Park.Center;
Park.Center = Park.Center - Park.Center;
if (Sgn(Start ^ End) < 0) swap(Start, End);
db Left = 0.0, Right = Angle; int Cnt = 0;
while (true) {
db MidAngle = (Left + Right) / 2.0;
Point MidPoint = Rotate(Start, MidAngle) * (Park.Radius / DisPointToPoint(Start, Park.Center));
if (Sgn(GetAngle(MidPoint - Park.Center, Start - MidPoint) - GetAngle(MidPoint - Park.Center, End - MidPoint)) == 0) {
printf("%.3Lf\n", DisPointToPoint(Start, MidPoint) + DisPointToPoint(End, MidPoint));
break;
}
else if (Sgn(GetAngle(MidPoint - Park.Center, Start - MidPoint) - GetAngle(MidPoint - Park.Center, End - MidPoint)) > 0) Right = MidAngle;
else Left = MidAngle;
}
}
else {
DisOS = DisPointToPoint(Park.Center, Start);
DisOE = DisPointToPoint(Park.Center, End);
db Temp1 = Park.Radius / DisOS, Temp2 = Park.Radius / DisOE;
db Angle = GetAngle(Start - Park.Center, End - Park.Center);
if (Sgn(Temp1 - 1.0) < 0) Angle -= acos(Temp1);
if (Sgn(Temp2 - 1.0) < 0) Angle -= acos(Temp2);
if (Sgn(DisOS - Park.Radius) > 0) Ans += sqrt(DisOS * DisOS - Park.Radius * Park.Radius);
if (Sgn(DisOE - Park.Radius) > 0) Ans += sqrt(DisOE * DisOE - Park.Radius * Park.Radius);
Ans += Angle * Park.Radius;
printf("%.3Lf\n", Ans);
}
}
return 0;
}
F. 小清新数论*
Description:
这是一道比较基础的数论题。
给出一个整数 n ,计算 ∑i=1n∑j=1nμ(gcd(i,j)) 。
其中 gcd(i,j) 表示 i,j 的最大公约数, μ(i) 表示莫比乌斯函数,它的一个等价定义是 μ(1)=1 , μ(n)=−∑d<n,d∣nμ(d) 当 n>1 时。
Input:
输入一行包含一个整数 n(1≤n≤107) 。
Output:
输出一行一个整数,表示答案。答案可能很大,请对 998244353 取模后输出。
Sample Input:
5
Sample Output:
14
Sample Input:
100
Sample Output:
3631
题目链接
题意为对于 n 求出 i=1∑nj=1∑nμ(gcd(i,j))
令 gcd(i,j)=k ,考虑枚举 k
k=1∑nμ(k)i=1∑nj=1∑n[gcd(i,j)=k]
化简该式
k=1∑nμ(k)i=1∑⌊kn⌋j=1∑⌊kn⌋[gcd(i,j)=1]
利用莫比乌斯函数的性质 d∣n∑μ(d)=[n=1] 可得
k=1∑nμ(k)i=1∑⌊kn⌋j=1∑⌊kn⌋d∣gcd(i,j)∑μ(d)
考虑枚举 d
k=1∑nμ(k)d=1∑⌊kn⌋μ(d)i=1∑⌊kn⌋[d∣i]j=1∑⌊kn⌋[d∣j]
其中 i=1∑⌊kn⌋[d∣i] 和 j=1∑⌊kn⌋[d∣j] 均为区间 [i,⌊kn⌋] 内 d 的倍数,可化简为 ⌊kdn⌋ 得
k=1∑nμ(k)d=1∑⌊kn⌋μ(d)⌊kdn⌋2
这样用线性筛筛出莫比乌斯函数并用分块、前缀和优化计算即可
AC代码:
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn = 1e7 + 5;
const int mod = 998244353;
bool IsPrime[maxn];
int Tot;
int Prime[maxn];
int Mu[maxn];
ll Prefix[maxn];
void Sieve() {
for (int i = 2; i < maxn; ++i) IsPrime[i] = true;
Mu[1] = 1;
for (ll i = 2; i < maxn; ++i) {
if (IsPrime[i]) {
Mu[i] = -1;
Prime[Tot++] = i;
}
for (ll j = 0; j < Tot && i * Prime[j] < maxn; ++j) {
IsPrime[i * Prime[j]] = false;
if (i % Prime[j] == 0) {
Mu[i * Prime[j]] = 0;
break;
}
Mu[i * Prime[j]] = -Mu[i];
}
}
for (ll i = 1; i < maxn; ++i)
Prefix[i] = (Prefix[i - 1] + Mu[i]) % mod;
}
ll CalF(ll Key) {
ll Ans = 0;
for (ll l = 1, tp; l <= Key; l = tp + 1) {
tp = Key / (Key / l);
Ans = (Ans + (Prefix[tp] - Prefix[l - 1] + mod) % mod * (Key / l) % mod * (Key / l) % mod) % mod;
}
return Ans;
}
ll CalAns(ll Key) {
ll Ans = 0;
for (ll l = 1, tp; l <= Key; l = tp + 1) {
tp = Key / (Key / l);
Ans = (Ans + (Prefix[tp] - Prefix[l - 1] + mod) * CalF(Key / l) % mod) % mod;
}
return Ans;
}
ll N;
int main(int argc, char *argv[]) {
Sieve();
scanf("%lld", &N);
printf("%lld\n", CalAns(N));
return 0;
}
其实算式 k=1∑nμ(k)d=1∑⌊kn⌋μ(d)⌊kdn⌋2 还可以化简
令 kd=T 并化简可得
T=1∑n⌊Tn⌋2k∣T∑μ(k)μ(kT)
其中 k∣T∑μ(k)μ(kT) 部分即为莫比乌斯函数 μ 与莫比乌斯函数 μ 的狄利克雷卷积,我们设它为 f
f(x)=k∣x∑μ(k)μ(kT)
所以现在的任务就是计算 f 这个函数的前缀和,因为莫比乌斯函数 μ 为积性函数,所以其卷积也为积性函数
积性函数的前缀和我们可以利用杜教筛来计算
令 S(n)=i=1∑nf(i) , S 为积性函数,根据杜教筛的核心式有
S(n)=i=1∑n(f∗1)(i)−i=2∑nS(⌊in⌋)
将 f 拆开
S(n)=i=1∑n(μ∗μ∗1)(i)−i=2∑nS(⌊in⌋)
其中 μ∗1=[n=1]=ϵ ,而 μ∗ϵ=μ ,所以上式可化简为
S(n)=i=1∑nμ(i)−i=2∑nS(⌊in⌋)
这样我们就可以利用杜教筛在 O(n32) 的时间复杂度内递归计算出莫比乌斯函数 μ 与莫比乌斯函数 μ 的狄利克雷卷积 f 的前缀和
但是我们注意到式子中还有一个 i=1∑nμ(i) 的莫比乌斯函数前缀和,而这个积性函数的前缀和也需要利用杜教筛来计算,具体请看
51nod 1244 莫比乌斯函数之和(杜教筛)
这样我们就可以通过Div1版本的小清新数论了
题目链接
AC代码:
#include <bits/stdc++.h>
using namespace std;
const int maxn = 1e7 + 5;
const int mod = 998244353;
bool IsPrime[maxn];
int Tot;
int Prime[maxn];
int Mu[maxn];
long long PrefixMu[maxn];
long long PrefixF[maxn];
void Sieve() {
for (int i = 2; i < maxn; ++i) IsPrime[i] = true;
Mu[1] = 1; PrefixF[1] = 1;
for (long long i = 2; i < maxn; ++i) {
if (IsPrime[i]) {
Mu[i] = -1;
PrefixF[i] = -2;
Prime[Tot++] = i;
}
for (long long j = 0; j < Tot && i * Prime[j] < maxn; ++j) {
IsPrime[i * Prime[j]] = false;
Mu[i * Prime[j]] = -Mu[i];
if (i % Prime[j] == 0) {
Mu[i * Prime[j]] = 0;
if ((i / Prime[j]) % Prime[j] == 0) PrefixF[i * Prime[j]] = 0;
else PrefixF[i * Prime[j]] = PrefixF[i / Prime[j]];
break;
}
else PrefixF[i * Prime[j]] = -2 * PrefixF[i];
}
}
for (int i = 1; i < maxn; ++i) {
PrefixMu[i] = (PrefixMu[i - 1] + Mu[i]) % mod;
PrefixF[i] = (PrefixF[i - 1] + PrefixF[i]) % mod;
}
}
long long __mod(long long Key1, long long Key2) {
long long Ans = Key1 - Key2;
while (Ans < 0) Ans += mod;
return Ans % mod;
}
map<long long, long long> SumMu;
long long SigmaMu(long long Key) {
if (Key < maxn) return PrefixMu[Key];
if (SumMu[Key]) return SumMu[Key];
long long Ans = 1;
for (long long l = 2, tp; l <= Key; l = tp + 1) {
tp = Key / (Key / l);
Ans -= (tp - l + 1) * SigmaMu(Key / l);
}
SumMu[Key] = Ans;
return SumMu[Key];
}
map<long long, long long> SumF;
long long SigmaF(long long Key) {
if (Key < maxn) return PrefixF[Key];
if (SumF[Key]) return SumF[Key];
long long Ans = SigmaMu(Key);
for (long long l = 2, tp; l <= Key; l = tp + 1) {
tp = Key / (Key / l);
Ans = __mod(Ans, (tp - l + 1) * SigmaF(Key / l) % mod);
}
SumF[Key] = Ans;
return Ans;
}
long long Cal(long long Key) {
long long Ans = 0;
for (long long l = 1, tp; l <= Key; l = tp + 1) {
tp = Key / (Key / l);
Ans = (Ans + __mod(SigmaF(tp), SigmaF(l - 1)) * (Key / l) % mod * (Key / l) % mod) % mod;
}
return Ans;
}
long long N;
int main(int argc, char *argv[]) {
Sieve();
scanf("%lld", &N);
printf("%lld\n", Cal(N));
return 0;
}