problem

luogu-P3340

題面寫得那么長，其實說白了就是求一條直線，使得若干個點到這條直線的距離平方的和最小，求這個最小值。

solution

我超愛數學，數學就是我的命，我一天不學數學我就難受！

假設擬合出的直線為：$y=kx+b$。

點到直線 $A{x}_0+B{y}_0+C=0$ 的距離公式：$d=\frac{|Ax+By+C|}{\sqrt{A^2+B^2}}$ 不會就去問數學老師

接下來就是推式子：
$$\frac{\sum (k{x}_i?y_i+b{)}^2}{k^2+1}=\frac{\sum k^2{x}_i^2+y_i^2+b^2?2k{x}_i{y}_i+2kb{x}_i?2b{y}_i}{k^2+1}$$
$x_i,y_i$ 都可被視作已知量，所以只有 $k,b$ 是未知數。

答案要最小化這個式子的值，$k,b$ 沒有任何等量關系，所以我們兩個各自最小化即可。

先最小化 $b$ 的求值，把 $k$ 假想成已知量。
$$=\frac{\sum (b^2+(2k{x}_i?2y_i)b+y_i^2?2k{x}_i{y}_i+k^2{x}_i^2)}{k^2+1}=\frac{n{b}^2+\sum (2k{x}_i?2y_i)b+\sum (y_i^2?2k{x}_i{y}_i+k^2{x}_i^2)}{k^2+1}$$
顯然這是個關于 $b$ 的二次函數。根據數學知識，我們知道當 $b=\frac{\sum (2y_i?2k{x}_i)}{2n}$ 時取到最低點。

且 $b=\frac{\sum (2y_i?2k{x}_i)}{2n}=\frac{\sum y_i?k\sum x_i}{n}$，恰好發現 $\frac{\sum x_i}{n}=\bar{x}$（平均數），$\frac{\sum y_i}{n}=\bar{y}$。

然后現在把 $b$ 當作常量 $b=\bar{y}?k\bar{x}$，把 $k$ 當作未知數。
$$=\frac{n(\bar{y}?k\bar{x}{)}^2+\sum (2k{x}_i?2y_i)(\bar{y}?k\bar{x})+\sum (y_i^2?2k{x}_i{y}_i+k^2{x}_i^2)}{k^2+1}$$

$$=\frac{n{\bar{y}}^2?2nk\bar{x}\bar{y}+n{k}^2{\bar{x}}^2+\sum (2k{x}_i\bar{y}?2y_i\bar{y}?2k^2{x}_i\bar{x}+2y_{i}k\bar{x})+\sum (y_i^2?2k{x}_i{y}_i+k^2{x}_i^2)}{k^2+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_i{y}_i)k+\sum ({\bar{y}}^2?2y_i\bar{y}+y_i^2)}{k^2+1}$$

為了式子的簡潔美觀，我們記：
$$?\begin{array}{l}A=\sum {\bar{x}}^2?\sum 2x_{i}x+\sum x_i^2\\ B=?2n\bar{x}\bar{y}+\sum 2x_i\bar{y}+\sum 2y_i\bar{x}?\sum 2x_i{y}_i\\ C=n{\bar{y}}^2?\sum 2y_i\bar{y}+\sum y_i^2\end{array}$$
則有 $\frac{A{k}^2+Bk+C}{k^2+1}\to ans$。
$$ans=\frac{A{k}^2+Bk+C}{k^2+1}?ans(k^2+1)=A{k}^2+Bk+C?(A?ans)k^2+Bk+C?ans=0$$
因為題目肯定是有解的，所以這個二元一次方程要有根，即 $\Delta \ge 0$。
KaTeX parse error: Undefined control sequence: \ at position 75: …ns+B^2-4AC\ge 0\? ?
又得到了一個關于 $ans$ 的開口向下的二次函數，$ans$ 最小取值就是其方程的較小根，再用一次求根公式即可。

所以我們只需要維護 $n,\sum x_i,\sum y_i,\sum x_i^2,\sum y_i^2,\sum x_i{y}_i$ 這幾個信息。

一條路徑上的信息維護，發現樹上差分即可做到，記從根到該點一路上的信息之和。

基環樹存一下環的點，然后把環拍成鏈，考慮有兩種走環的方式，比較最小值即可。

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; }

總結

以上是生活随笔為你收集整理的[ZJOI2014] 星系调查（树上差分 + 数学推式子）的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。