BZOJ5020: [THUWC 2017]在美妙的数学王国中畅游

2024-02-04 11:49:04

其实题面好像有点不全，建议去洛谷：

P4546 [THUWC2017]在美妙的数学王国中畅游

这里还是$BZOJ$的题面。

Description

数字和数学规律主宰着这个世界。机器的运转，生命的消长，宇宙的进程，这些神秘而又美妙的过程无不可以用数学的语言展现出来。这印证了一句古老的名言： “学好数理化，走遍天下都不怕。” 学渣小R被大学的数学课程虐得生活不能自理，微积分的成绩曾是他在教室里上的课的最低分。然而他的某位陈姓室友却能轻松地在数学考试中得到满分。为了提升自己的数学课成绩，有一天晚上（在他睡觉的时候），他来到了数学王国。数学王国中，每个人的智商可以用一个属于 [0,1]的实数表示。数学王国中有 n 个城市，编号从 0 到 n−1 ，这些城市由若干座魔法桥连接。每个城市的中心都有一个魔法球，每个魔法球中藏有一道数学题。每个人在做完这道数学题之后都会得到一个在 [0,1] 区间内的分数。一道题可以用一个从 [0,1] 映射到 [0,1]的函数 f(x) 表示。若一个人的智商为 x ，则他做完这道数学题之后会得到 f(x)分。函数 f有三种形式：正弦函数 sin(ax+b) (a∈[0,1],b∈[0,π],a+b∈[0,π]) 指数函数 e^(ax+b) (a∈[−1,1],b∈[−2,0],a+b∈[−2,0]) 一次函数 ax+b (a∈[−1,1],b∈[0,1],a+b∈[0,1] 数学王国中的魔法桥会发生变化，有时会有一座魔法桥消失，有时会有一座魔法桥出现。但在任意时刻，只存在至多一条连接任意两个城市的简单路径（即所有城市形成一个森林）。在初始情况下，数学王国中不存在任何的魔法桥。数学王国的国王拉格朗日很乐意传授小R数学知识，但前提是小R要先回答国王的问题。这些问题具有相同的形式，即一个智商为 x 的人从城市 u 旅行到城市 v（即经过 u 到 v 这条路径上的所有城市，包括 u和 v ）且做了所有城市内的数学题后，他所有得分的总和是多少。

Input

第一行两个正整数 n,m 和一个字符串 type 。表示数学王国*有 n 座城市，发生了 m 个事件，该数据的类型为 type 。 typet 字符串是为了能让大家更方便地获得部分分，你可能不需要用到这个输入。其具体含义在【数据范围与提示】中有解释。接下来 n 行，第 i 行表示初始情况下编号为 i 的城市的魔法球中的函数。一个魔法用一个整数 f表示函数的类型，两个实数 a,b 表示函数的参数，若 f=1,则函数为 f(x)=sin(ax+b)(a∈[0,1],b∈[0,π],a+b∈[0,π]) f=2,则函数为 f(x)=e^(ax+b)(a∈[−1,1],b∈[−2,0],a+b∈[−2,0]) f=3,则函数为 f(x)=ax+b(a∈[−1,1],b∈[0,1],a+b∈[0,1]) 接下来 m行，每行描述一个事件，事件分为四类。 appear u v 表示数学王国中出现了一条连接 u 和 v 这两座城市的魔法桥 (0≤u,v<n,u≠v) ，保证连接前 u和 v 这两座城市不能互相到达。 disappear u v 表示数学王国中连接 u 和 v 这两座城市的魔法桥消失了，保证这座魔法桥是存在的。 magic c f a b 表示城市 c 的魔法球中的魔法变成了类型为 f ，参数为 a,b 的函数 travel u v x 表示询问一个智商为 x 的人从城市 u 旅行到城市 v （即经过 u到 v 这条路径上的所有城市，包括 u 和 v ）后，他得分的总和是多少。若无法从 u 到达 v ，则输出一行一个字符串 unreachable。 1≤n≤100000,1≤m≤200000

Output

对于每个询问，输出一行实数，表示得分的总和。

Sample Input

3 7 C1
1 1 0
3 0.5 0.5
3 -0.5 0.7
appear 0 1
travel 0 1 0.3
appear 0 2
travel 1 2 0.5
disappear 0 1
appear 1 2
travel 1 2 0.5

Sample Output

9.45520207e-001
1.67942554e+000
1.20000000e+000

题解Here!

数字和数学规律主宰着这个世界。。。

学好数理化，走遍天下都不怕。。。

这个题要是没有询问就是一个沙茶$LCT$。但是这个询问很难搞。。。现在的$LCT$怎么都和一些奇奇怪怪的东西搞在一起。。。

其实原题面有一个提示，我们称它为——泰勒展开。

如下：（选自洛谷）

若函数$f(x)$的$n$阶导数在$[a,b]$区间内连续，则对$f(x)$在$x_0(x_0\in[a,b])$处使用$n$次拉格朗日中值定理可以得到带拉格朗日余项的泰勒展开式：

$$f(x)=f(x_0)+\frac{f'(x_0)(x-x_0)}{1!}+\frac{f''(x_0)(x-x_0)^2}{2!}+ \cdots +\frac{f^{(n-1)}(x_0)(x-x_0)^{n-1}}{(n-1)!}+\frac{f^{(n)}(\xi)(x-x_0)^n}{n!},x\in[a,b]$$

其中，当$x>x_0$时，$\xi\in[x_0,x]$；当$x<x_0$时，$\xi\in[x,x_0]$。

$f^{(n)}(x)$表示函数$f(x)$的$n$阶导数。

其实别看这么多乱七八糟的名词，关键就在于那个式子。

我们注意到，等号的左边是$f(x)$！

也就是说，我们可以不用搞其它的东西，只要把$f(x)$的$n$阶导数搞出来就好。

什么？你不知道导数？

回去学《数学选修2-2》去吧。。。

这里给出此题要用到的导数计算：

基本公式：

$$f(x)=x\quad\Rightarrow\quad f'(x)=1$$

$$f(x)=\sin(x)\quad\Rightarrow\quad f'(x)=\cos(x)$$

$$f(x)=\cos(x)\quad\Rightarrow\quad f'(x)=-\sin(x)$$

$$f(x)=e^x\quad\Rightarrow\quad f'(x)=e^x$$

计算公式：

$$[f(x)+g(x)]'=f'(x)+g'(x)$$

$$[f(x)g(x)]'=f'(x)g(x)+f(x)g'(x)$$

$$\text{链式法则：}[f(g(x))]'=f'(g(x))\times g'(x)$$

于是我们可以愉快地推式子了！

以$f(x)=sin(x)$为例：

泰勒展开式中选取了一个$x_0$，我们为了方便，直接$x_0=0$。

设$f^n(x)$为$f(x)$的$n$阶导数。

那么：$$f(x)=f(0)+\frac{f^1(0)(x-0)}{1!}+\frac{f^2(0)(x-0)^2}{2!}+ \cdots +\frac{f^{n-1}(0)(x-0)^{n-1}}{(n-1)!}+\frac{f^{n}(\xi)(x-0)^n}{n!}$$

化简一下：$$f(x)=f(0)+\frac{f^1(0)}{1!}x+\frac{f^2(0)}{2!}x^2+ \cdots +\frac{f^{n-1}(0)}{(n-1)!}x^{n-1}+\frac{f^{n}(\xi)}{n!}x^n$$

但是这个$n$到底是多大呢？

其实，我们发现$n!$在$n$比较小时就是暴增函数。

也就是说，越往后的项对答案贡献越小，而精度也有$SPJ$，所以我们并不需要末尾那么多项。

这里我选取了$n=16$这个值，这是可以过的。

于是我们只要对每个点维护$16$项多项式就好。

但是又有人问了：题目中不是求$f(x)=sin(ax+b)$吗？

其实根据上面那个链式法则我们可以得出：$$f(x)=\sin(ax+b)\quad\Rightarrow\quad f'(x)=a\cos(ax+b)$$

同理对于另外两个式子有：

$$f(x)=ax+b\quad\Rightarrow\quad f'(x)=a$$

$$f(x)=e^{ax+b}\quad\Rightarrow\quad f'(x)=ae^{ax+b}$$

但是这只是一阶导数啊，高阶导数怎么办？

对于$f(x)=ax+b$，只有一阶导数$f'(x)=a$，高阶导数全部为$0$。

对于$f(x)=e^{ax+b}$，一阶导数$f'(x)=ae^{ax+b}$，然后每高一阶多乘一个$a$，如二阶导数为$f''(x)=a^2e^{ax+b}$。

对于$f(x)=\sin(ax+b)$，一阶导数$f'(x)=a\cos(ax+b)$，然后每高一阶多乘一个$a$，$\text{阶数}\mod 2==0$的导数为$\sin$，否则为$\cos$，并且$\text{阶数}\mod 4\leq1$的导数为正，否则为负。

由于精度会出锅，再加上每次计算$\sin(x),\cos(x),e^x$会比较慢，于是每次计算之前预处理一下即可。

至于阶乘。。。我一开始写了个逆元，发现模数不知道。。。

我这个大制杖好久才反应过来$16!=2.0922789888\times 10^{13}\leq 10^{18}$。。。

所以直接开$double$预处理逆元就好。

$LCT$应该不用多说，维护多项式的每一项单独的和。

答案就是多项式每一项之和。

至于那个什么科学计数法。。。大多数人直接$\%.8lf$。。。

其实我们可以直接用$C++$的自带版本——$\%.8e$！

自动将小数转成科学计数法！

还有啥细节那就看代码吧。。。

附代码：

#include<iostream>
#include<algorithm>
#include<cstdio>
#include<cmath>
#include<cstring>
#define MAXN 100010
#define MAXM 20
using namespace std;
int n,m;
double fact[MAXM],val[MAXN][MAXM];
int top,stack[MAXN];
struct Link_Cut_Tree{
	int f,flag,son[2];
	double v[MAXM];
}a[MAXN];
inline int read(){
	int date=0,w=1;char c=0;
	while(c<'0'||c>'9'){if(c=='-')w=-1;c=getchar();}
	while(c>='0'&&c<='9'){date=date*10+c-'0';c=getchar();}
	return date*w;
}
inline bool isroot(int rt){
	return a[a[rt].f].son[0]!=rt&&a[a[rt].f].son[1]!=rt;
}
inline void pushup(int rt){
	if(!rt)return;
	for(int i=0;i<=15;i++)a[rt].v[i]=a[a[rt].son[0]].v[i]+a[a[rt].son[1]].v[i]+val[rt][i];
}
inline void pushdown(int rt){
	if(!rt||!a[rt].flag)return;
	a[a[rt].son[0]].flag^=1;a[a[rt].son[1]].flag^=1;a[rt].flag^=1;
	swap(a[rt].son[0],a[rt].son[1]);
}
inline void turn(int rt){
	int x=a[rt].f,y=a[x].f,k=a[x].son[0]==rt?1:0;
	if(!isroot(x)){
		if(a[y].son[0]==x)a[y].son[0]=rt;
		else a[y].son[1]=rt;
	}
	a[rt].f=y;a[x].f=rt;a[a[rt].son[k]].f=x;
	a[x].son[k^1]=a[rt].son[k];a[rt].son[k]=x;
	pushup(x);pushup(rt);
}
void splay(int rt){
	top=0;
	stack[++top]=rt;
	for(int i=rt;!isroot(i);i=a[i].f)stack[++top]=a[i].f;
	while(top)pushdown(stack[top--]);
	while(!isroot(rt)){
		int x=a[rt].f,y=a[x].f;
		if(!isroot(x)){
			if((a[y].son[0]==x)^(a[x].son[0]==rt))turn(rt);
			else turn(x);
		}
		turn(rt);
	}
}
void access(int rt){
	for(int i=0;rt;i=rt,rt=a[rt].f){
		splay(rt);
		a[rt].son[1]=i;
		pushup(rt);
	}
}
inline void makeroot(int rt){access(rt);splay(rt);a[rt].flag^=1;}
int findroot(int rt){
	access(rt);splay(rt);
	while(a[rt].son[0])rt=a[rt].son[0];
	return rt;
}
inline void split(int x,int y){makeroot(x);access(y);splay(y);}
inline void link(int x,int y){makeroot(x);a[x].f=y;}
inline void cut(int x,int y){
	split(x,y);
	if(a[y].son[0]==x&&a[x].f==y&&!a[x].son[1])a[x].f=a[y].son[0]=0;
}
void change(int x,int f,double a,double b){
	memset(val[x],0,sizeof(val[x]));
	if(f==1){
		double t=1,Sin=sin(b),Cos=cos(b);
		for(int i=0;i<=15;i++){
			val[x][i]=((i%2)?Cos:Sin)*((i%4<=1)?1:-1)*t/fact[i];
			t*=a;
		}
	}
	else if(f==2){
		double t=1,Exp=exp(b);
		for(int i=0;i<=15;i++){
			val[x][i]=Exp*t/fact[i];
			t*=a;
		}
	}
	else{
		val[x][0]=b;
		val[x][1]=a;
	}
}
double query(int x,int y,double k){
	split(x,y);
	double s=0,t=1;
	for(int i=0;i<=15;i++){
		s+=a[y].v[i]*t;
		t*=k;
	}
	return s;
}
void work(){
	int u,v;
	double x,y,k;
	char ch[20];
	while(m--){
		scanf("%s",ch);u=read();v=read();
		switch(ch[0]){
			case 'a':{
				u++;v++;
				link(u,v);
				break;
			}
			case 'd':{
				u++;v++;
				cut(u,v);
				break;
			}
			case 'm':{
				u++;
				access(u);splay(u);
				scanf("%lf%lf",&x,&y);
				change(u,v,x,y);
				pushup(u);
				break;
			}
			case 't':{
				u++;v++;scanf("%lf",&k);
				if(findroot(u)!=findroot(v))printf("unreachable\n");
				else printf("%.8e\n",query(u,v,k));
				break;
			}
		}
	}
}
void init(){
	int f;
	double x,y;
	char ch[5];
	n=read();m=read();scanf("%s",ch);
	fact[0]=1;
	for(int i=1;i<=16;i++)fact[i]=fact[i-1]*i;
	for(int i=1;i<=n;i++){
		f=read();
		scanf("%lf%lf",&x,&y);
		change(i,f,x,y);
	}
}
int main(){
	init();
	work();
    return 0;
}

码农公寓