[ZJOI2014] 星系调查(树上差分 + 数学推式子)
problem
luogu-P3340
題面寫得那么長(zhǎng),其實(shí)說(shuō)白了就是求一條直線,使得若干個(gè)點(diǎn)到這條直線的距離平方的和最小,求這個(gè)最小值。
solution
我超愛(ài)數(shù)學(xué),數(shù)學(xué)就是我的命,我一天不學(xué)數(shù)學(xué)我就難受!
假設(shè)擬合出的直線為:y=kx+by=kx+by=kx+b。
點(diǎn)到直線 Ax0+By0+C=0Ax_0+By_0+C=0Ax0?+By0?+C=0 的距離公式:d=∣Ax+By+C∣A2+B2d=\frac{|Ax+By+C|}{\sqrt{A^2+B^2}}d=A2+B2?∣Ax+By+C∣? 不會(huì)就去問(wèn)數(shù)學(xué)老師
接下來(lái)就是推式子:
∑(kxi?yi+b)2k2+1=∑k2xi2+yi2+b2?2kxiyi+2kbxi?2byik2+1\frac{\sum(kx_i-y_i+b)^2}{k^2+1}=\frac{\sum k^2x_i^2+y_i^2+b^2-2kx_iy_i+2kbx_i-2by_i}{k^2+1} k2+1∑(kxi??yi?+b)2?=k2+1∑k2xi2?+yi2?+b2?2kxi?yi?+2kbxi??2byi??
xi,yix_i,y_ixi?,yi? 都可被視作已知量,所以只有 k,bk,bk,b 是未知數(shù)。
答案要最小化這個(gè)式子的值,k,bk,bk,b 沒(méi)有任何等量關(guān)系,所以我們兩個(gè)各自最小化即可。
先最小化 bbb 的求值,把 kkk 假想成已知量。
=∑(b2+(2kxi?2yi)b+yi2?2kxiyi+k2xi2)k2+1=nb2+∑(2kxi?2yi)b+∑(yi2?2kxiyi+k2xi2)k2+1=\frac{\sum(b^2+(2kx_i-2y_i)b+y_i^2-2kx_iy_i+k^2x_i^2)}{k^2+1}=\frac{nb^2+\sum(2kx_i-2y_i)b+\sum(y_i^2-2kx_iy_i+k^2x_i^2)}{k^2+1} =k2+1∑(b2+(2kxi??2yi?)b+yi2??2kxi?yi?+k2xi2?)?=k2+1nb2+∑(2kxi??2yi?)b+∑(yi2??2kxi?yi?+k2xi2?)?
顯然這是個(gè)關(guān)于 bbb 的二次函數(shù)。根據(jù)數(shù)學(xué)知識(shí),我們知道當(dāng) b=∑(2yi?2kxi)2nb=\frac{\sum(2y_i-2kx_i)}{2n}b=2n∑(2yi??2kxi?)? 時(shí)取到最低點(diǎn)。
且 b=∑(2yi?2kxi)2n=∑yi?k∑xinb=\frac{\sum(2y_i-2kx_i)}{2n}=\frac{\sum y_i-k\sum x_i}{n}b=2n∑(2yi??2kxi?)?=n∑yi??k∑xi??,恰好發(fā)現(xiàn) ∑xin=xˉ\frac{\sum x_i}{n}=\bar{x}n∑xi??=xˉ(平均數(shù)),∑yin=yˉ\frac{\sum y_i}{n}=\bar{y}n∑yi??=yˉ?。
然后現(xiàn)在把 bbb 當(dāng)作常量 b=yˉ?kxˉb=\bar{y}-k\bar{x}b=yˉ??kxˉ,把 kkk 當(dāng)作未知數(shù)。
=n(yˉ?kxˉ)2+∑(2kxi?2yi)(yˉ?kxˉ)+∑(yi2?2kxiyi+k2xi2)k2+1=\frac{n(\bar{y}-k\bar{x})^2+\sum(2kx_i-2y_i)(\bar{y}-k\bar{x})+\sum(y_i^2-2kx_iy_i+k^2x_i^2)}{k^2+1} =k2+1n(yˉ??kxˉ)2+∑(2kxi??2yi?)(yˉ??kxˉ)+∑(yi2??2kxi?yi?+k2xi2?)?
=nyˉ2?2nkxˉyˉ+nk2xˉ2+∑(2kxiyˉ?2yiyˉ?2k2xixˉ+2yikxˉ)+∑(yi2?2kxiyi+k2xi2)k2+1=\frac{n\bar{y}^2-2nk\bar{x}\bar{y}+nk^2\bar{x}^2+\sum(2kx_i\bar{y}-2y_i\bar{y}-2k^2x_i\bar{x}+2y_ik\bar{x})+\sum(y_i^2-2kx_iy_i+k^2x_i^2)}{k^2+1} =k2+1nyˉ?2?2nkxˉyˉ?+nk2xˉ2+∑(2kxi?yˉ??2yi?yˉ??2k2xi?xˉ+2yi?kxˉ)+∑(yi2??2kxi?yi?+k2xi2?)?
=∑(xˉ2?2xixˉ+xi2)k2+∑(?2xˉyˉ+2xiyˉ+2yixˉ?2xiyi)k+∑(yˉ2?2yiyˉ+yi2)k2+1=\frac{\sum(\bar{x}^2-2x_i\bar{x}+x_i^2)k^2+\sum(-2\bar{x}\bar{y}+2x_i\bar{y}+2y_i\bar{x}-2x_iy_i)k+\sum(\bar{y}^2-2y_i\bar{y}+y_i^2)}{k^2+1} =k2+1∑(xˉ2?2xi?xˉ+xi2?)k2+∑(?2xˉyˉ?+2xi?yˉ?+2yi?xˉ?2xi?yi?)k+∑(yˉ?2?2yi?yˉ?+yi2?)?
為了式子的簡(jiǎn)潔美觀,我們記:
{A=∑xˉ2?∑2xix+∑xi2B=?2nxˉyˉ+∑2xiyˉ+∑2yixˉ?∑2xiyiC=nyˉ2?∑2yiyˉ+∑yi2\begin{cases} A=\sum\bar{x}^2-\sum 2x_ix+\sum x_i^2\\ B=-2n\bar{x}\bar{y}+\sum 2x_i\bar{y}+\sum2y_i\bar{x}-\sum2x_iy_i\\ C=n\bar{y}^2-\sum2y_i\bar{y}+\sum y_i^2 \end{cases} ??????A=∑xˉ2?∑2xi?x+∑xi2?B=?2nxˉyˉ?+∑2xi?yˉ?+∑2yi?xˉ?∑2xi?yi?C=nyˉ?2?∑2yi?yˉ?+∑yi2??
則有 Ak2+Bk+Ck2+1→ans\frac{Ak^2+Bk+C}{k^2+1}\rightarrow ansk2+1Ak2+Bk+C?→ans。
ans=Ak2+Bk+Ck2+1?ans(k2+1)=Ak2+Bk+C?(A?ans)k2+Bk+C?ans=0ans=\frac{Ak^2+Bk+C}{k^2+1}\Rightarrow ans(k^2+1)=Ak^2+Bk+C\Rightarrow (A-ans)k^2+Bk+C-ans=0 ans=k2+1Ak2+Bk+C??ans(k2+1)=Ak2+Bk+C?(A?ans)k2+Bk+C?ans=0
因?yàn)轭}目肯定是有解的,所以這個(gè)二元一次方程要有根,即 Δ≥0\Delta\ge 0Δ≥0。
KaTeX parse error: Undefined control sequence: \ at position 75: …ns+B^2-4AC\ge 0\? ?
又得到了一個(gè)關(guān)于 ansansans 的開(kāi)口向下的二次函數(shù),ansansans 最小取值就是其方程的較小根,再用一次求根公式即可。
所以我們只需要維護(hù) n,∑xi,∑yi,∑xi2,∑yi2,∑xiyin,\sum x_i,\sum y_i,\sum x_i^2,\sum y_i^2,\sum x_iy_in,∑xi?,∑yi?,∑xi2?,∑yi2?,∑xi?yi? 這幾個(gè)信息。
一條路徑上的信息維護(hù),發(fā)現(xiàn)樹(shù)上差分即可做到,記從根到該點(diǎn)一路上的信息之和。
基環(huán)樹(shù)存一下環(huán)的點(diǎn),然后把環(huán)拍成鏈,考慮有兩種走環(huán)的方式,比較最小值即可。
code
#include <bits/stdc++.h> using namespace std; #define maxn 100005 struct node {int n, sumx, sumy, sumx2, sumy2, sumxy;void insert( int x, int y ) { n ++, sumx += x, sumy += y, sumx2 += x * x, sumy2 += y * y, sumxy += x * y; }double calc() {double var_x = sumx * 1.0 / n, var_y = sumy * 1.0 / n;double A = sumx2 - 2 * var_x * sumx + n * var_x * var_x;double B = - 2 * sumxy + 2 * sumx * var_y + 2 * sumy * var_x - 2 * n * var_x * var_y;double C = sumy2 - 2 * sumy * var_y + n * var_y * var_y;double a = 4, b = - 4 * ( A + C ), c = 4 * A * C - B * B;double delta = sqrt( b * b - 4 * a * c );return ( - b - delta ) / ( 2 * a );} }f[maxn], g[maxn]; node operator + ( node HM, node NB ) { HM.n += NB.n, HM.sumx += NB.sumx, HM.sumy += NB.sumy, HM.sumx2 += NB.sumx2, HM.sumy2 += NB.sumy2, HM.sumxy += NB.sumxy; return HM; } node operator - ( node HM, node NB ) { HM.n -= NB.n, HM.sumx -= NB.sumx, HM.sumy -= NB.sumy, HM.sumx2 -= NB.sumx2, HM.sumy2 -= NB.sumy2, HM.sumxy -= NB.sumxy; return HM; } struct point { int x, y; }p[maxn]; int n, m, Q, cnt; int fa[maxn], dep[maxn], top[maxn], root[maxn], vis[maxn], siz[maxn], son[maxn], num[maxn], circle[maxn]; vector < int > G[maxn];int lca( int u, int v ) {while( top[u] ^ top[v] )if( dep[top[u]] < dep[top[v]] ) v = fa[top[v]];else u = fa[top[u]];return dep[u] < dep[v] ? u : v; }void dfs1( int u, int rt ) {root[u] = rt, siz[u] = vis[u] = 1, son[u] = 0;f[u] = f[fa[u]], f[u].insert( p[u].x, p[u].y );for( int v : G[u] ) if( ! vis[v] and v ^ fa[u] ) {fa[v] = u;dep[v] = dep[u] + 1;dfs1( v, rt );siz[u] += siz[v];if( siz[son[u]] < siz[v] ) son[u] = v;} }void dfs2( int u, int t ) {top[u] = t;if( son[u] ) dfs2( son[u], t );else return;for( int v : G[u] )if( fa[v] == u and son[u] ^ v ) dfs2( v, v ); }void work() {for( int u = 1;u <= n;u ++ )for( int v : G[u] )if( fa[v] ^ u and fa[u] ^ v ) {if( dep[u] > dep[v] ) swap( u, v );while( v ^ u ) {circle[++ cnt] = v;num[v] = cnt;g[cnt] = g[cnt - 1];g[cnt].insert( p[v].x, p[v].y );vis[v] = 1; v = fa[v];}circle[++ cnt] = u;num[u] = cnt;g[cnt] = g[cnt - 1];g[cnt].insert( p[u].x, p[u].y );vis[u] = 1;return;} }int main() {scanf( "%d %d", &n, &m );for( int i = 1;i <= n;i ++ ) scanf( "%d %d", &p[i].x, &p[i].y );for( int i = 1, u, v;i <= m;i ++ ) {scanf( "%d %d", &u, &v );G[u].push_back( v );G[v].push_back( u );}dfs1( 1, 1 );memset( vis, 0, sizeof( vis ) );if( n == m ) work();else circle[cnt = 1] = 1;memset( fa, 0, sizeof( fa ) );for( int i = 1;i <= cnt;i ++ ) dfs1( circle[i], circle[i] ), dfs2( circle[i], circle[i] );scanf( "%d", &Q );while( Q -- ) {int x, y;scanf( "%d %d", &x, &y );if( root[x] == root[y] )printf( "%.5f\n", (f[x] + f[y] - f[lca( x, y )] - f[fa[lca( x, y )]] ).calc() );else {if( num[root[x]] > num[root[y]] ) swap( x, y );node u = f[x] - f[root[x]] + f[y] - f[root[y]] + g[num[root[y]]] - g[num[root[x]] - 1];node v = f[x] - f[root[x]] + f[y] - f[root[y]] + g[num[root[x]]] + g[cnt] - g[num[root[y]] - 1];printf( "%.5f\n", min( u.calc(), v.calc() ) );}}return 0; }總結(jié)
以上是生活随笔為你收集整理的[ZJOI2014] 星系调查(树上差分 + 数学推式子)的全部?jī)?nèi)容,希望文章能夠幫你解決所遇到的問(wèn)題。
- 上一篇: 北戴河旅游景点介绍 五大必去景点
- 下一篇: 微信拍一拍怎么改字 微信拍一拍改字方法