保序回归问题
作者：互联网
保序回归问题
基本形式
\[f(x)=\sum_{i=1}^nw_i|a_i-b_i|^k \]有一些要求形如
\[b_x\le b_y \]最小化\(f(x)\)
一般解法
我们可以整体二分，对于值域区间\([l,r]\)，我们二分\(mid=\frac{l+r}{2}\)，对于当前需要考虑的变量，判断它的取值在\([l,mid]\)更优还是\((mid,r]\)更优。
有结论：
取值在\([l,mid]\)更优还是\((mid,r]\)更优等价于取\(mid\)还是\(mid+\theta\)更优，其中\(\theta\)是根据题目要求的精度而确定的极小值
所以我们现在判断一个点取\(mid\)更优还是\(mid+\theta\)更优即可
考虑把这个问题转化为最小权闭合子图问题。我们先钦定所有点选择\(mid\)，如果一个点变成了\(mid+\theta\)，那么所有有要求\(b_x\le b_y\)的点的\(y\)也必须选择\(mid+\theta\)。也就是每个点的点权为\(f(mid+\theta)-f(mid)\)，求最小权闭合子图。
根据我们网络流的结果判断每个点取值在\([l,mid]\)更优还是\((mid,r]\)更优，分治解决即可。
例题
问题的形式是
\[f(x)=\sum_{i=1}^nw_i(a_i-b_i)^2 \]我们只需要找到偏序关系即可。
对\(A,B\)分别建线性基。如果我们将一个元素\(x\in A\)从线性基中删除，而另一个元素\(y\not\in A\)可以插入，那么用\(y\)替换\(x\)可以使得集合大小不变，所以\(b_x\le b_y\)。同理，如果将一个元素\(x\in B\)从线性基中删除，而另一个元素\(y\not\in B\)可以插入，要求\(b_x\ge b_y\)。
而这样就一定能覆盖所有的情况吗？答案是肯定的。因为\(A,B\)已经是满足条件的最大的集合了，所以必须一一替换。而因为线性基可以用唯一的方案表示一个数，所以如果同时删去线性基中一个集合\(S\)，可以被插入的所有数一定需要有\(S\)里的元素来表示，而\(S\)里的所有元素能够替换的数都已经连好了边，所以可以插入的所有数一定是\(S\)里每个元素能够替换的并集。
因为取值必须为整数，\(mid\) 取$ \left\lfloor\frac{l+r}{2}\right\rfloor$， \(\theta=1\)。然后就分治即可
#include <bits/stdc++.h>
using namespace std;
const int N=1005,M=N*N*4,inf=1e9;
#define ull unsigned long long
int n,m,v[N],a[N],b[N],mx;
bool visa[N],visb[N];
ull buc[64],c[N];
inline int insert(ull x){
for(int i=63;~i;--i){
if(!((x>>i)&1))continue;
if(buc[i])x^=buc[i];
else{buc[i]=x;return i;}
}
return -1;
}
namespace Flow{
struct Edge{int to,next,flow;}e[M<<1];
int head[N],ecnt=1,cur[N];
inline void adde(int u,int v,int flow){
e[++ecnt]=(Edge){v,head[u],flow};head[u]=ecnt;
e[++ecnt]=(Edge){u,head[v],0};head[v]=ecnt;
}
struct Queue{
int q[N],head,tail;
inline void init(){head=1;tail=0;}
inline void push(int x){q[++tail]=x;}
inline int front(){return q[head];}
inline void pop(){++head;}
inline bool empty(){return head>tail;}
}q;
int s,t,dep[N];
bool bfs(){
memset(dep+1,-1,t<<2);
dep[s]=0;q.init();q.push(s);cur[s]=head[s];
while(!q.empty()){
int u=q.front();q.pop();
for(int i=head[u],v;i;i=e[i].next){
if(e[i].flow<=0)continue;v=e[i].to;
if(dep[v]==-1){dep[v]=dep[u]+1;cur[v]=head[v];q.push(v);}
}
}
return dep[t]>=0;
}
int dinic(int u,int flow){
if(u==t)return flow;
int ret=0,f;
for(int i=cur[u],v;i;i=e[i].next){
v=e[i].to;cur[u]=i;
if(e[i].flow<=0||dep[v]!=dep[u]+1)continue;
f=dinic(v,min(flow,e[i].flow));
ret+=f;flow-=f;e[i].flow-=f;e[i^1].flow+=f;
if(flow<=0)break;
}
return ret;
}
inline int Maxflow(int S,int T){
s=S;t=T;
int ans=0;
while(bfs())ans+=dinic(s,inf);
return ans;
}
inline void init(){
memset(head,0,sizeof(head));
ecnt=1;
}
}
vector<int> to[N];
#define pb push_back
int f[N],q[N];
#define ll long long
inline ll calc(int x,int y){
return 1ll*(y-v[x])*(y-v[x]);
}
int id[N],q1[N],q2[N],vis[N],inde;
void work(int l,int r,int ql,int qr){
if(l==r){
for(int i=ql;i<=qr;++i)f[q[i]]=l;
return;
}
if(ql>qr)return;
int mid=(l+r)>>1;++inde;
for(int i=ql;i<=qr;++i)id[q[i]]=i-ql+1,vis[q[i]]=inde;
int s=qr-ql+2,t=s+1;
Flow::init();
for(int i=ql,u;i<=qr;++i){
u=q[i];
for(int v:to[u])if(vis[v]==inde)
Flow::adde(id[u],id[v],inf);
ll x=calc(u,mid)-calc(u,mid+1);
if(x>0)Flow::adde(s,id[u],x);
else Flow::adde(id[u],t,-x);
}
Flow::Maxflow(s,t);
int tot1=0,tot2=0;
for(int i=ql;i<=qr;++i)
if(Flow::dep[id[q[i]]]==-1)q1[++tot1]=q[i];
else q2[++tot2]=q[i];
memcpy(q+ql,q1+1,tot1<<2);
memcpy(q+ql+tot1,q2+1,tot2<<2);
work(l,mid,ql,ql+tot1-1);
work(mid+1,r,ql+tot1,qr);
}
int main(){
scanf("%d %d",&n,&m);
for(int i=1;i<=n;++i)scanf("%llu",&c[i]);
for(int i=1;i<=n;++i)scanf("%d",&v[i]);
for(int i=1;i<=m;++i)scanf("%d",&a[i]),visa[a[i]]=1;
for(int i=1;i<=m;++i)scanf("%d",&b[i]),visb[b[i]]=1;
for(int i=1;i<=m;++i){
memset(buc,0,sizeof(buc));
for(int j=1;j<=m;++j)if(i!=j)insert(c[a[j]]);
for(int j=1,p;j<=n;++j){
if(visa[j])continue;
p=insert(c[j]);
if(p>=0)buc[p]=0,to[a[i]].pb(j);
}
}
for(int i=1;i<=m;++i){
memset(buc,0,sizeof(buc));
for(int j=1;j<=m;++j)if(i!=j)insert(c[b[j]]);
for(int j=1,p;j<=n;++j){
if(visb[j])continue;
p=insert(c[j]);
if(p>=0)buc[p]=0,to[j].pb(b[i]);
}
}
for(int i=1;i<=n;++i)q[i]=i,mx=max(mx,v[i]);
work(0,mx,1,n);
ll ans=0;
for(int i=1;i<=n;++i)ans+=calc(i,f[i]);
printf("%lld\n",ans);
return 0;
}
如果不考虑图的特殊性（基环树）的话，直接用保序回归的一般解法即可。
但是图具有特殊性，我们考虑使用更高效的DP做法替换最小权闭合子图的网络流做法。
也就是说，我们需要在线性复杂度内，确定每个点选\(mid\)还是\(mid+1\)。
- 如果图是一棵树
设\(f_{u,0/1}\)表示该点选择\(mid/mid+1\)的最小代价。
如果存在限制\(b_u\le b_v\)，\(f_{u,0}\leftarrow f_{v,0}+f_{v,1},\ \ f_{u,1}\leftarrow f_{v,1}\)
否则\(b_u\ge b_v\)，\(f_{u,1}\leftarrow f_{v,0}+f_{v,1},\ \ f_{u,0}\leftarrow f_{v,0}\)
随便选一个点为根DFS，确定了根的颜色之后染色。如果\(b_{fa}\le b_{u}\)且\(b_{fa}\)取了\(mid+1\)，那么\(u\)必须取\(mid+1\)，如果\(b_{fa}\ge b_{u}\)且\(b_{fa}\)取了\(mid\)，那么\(u\)必须取\(mid\)。其它情况根据DP值取。
- 图是一个基环树
设\(f_{u,0/1,0/1}\)表示根选择\(mid/mid+1\)，该点选择\(mid/mid+1\)的最小代价。其中根为环上任意一点
转移与上面类似。\(f_{?,0,?}\)和\(f_{?,1,?}\)分别转移
但是因为基环树多一条边的限制，假设这一条环边是\((rt,u)\)，那么需要根据\(b_{rt}\)和\(b_{u}\)的大小关系，将\(f_{u,0,1}\)或\(f_{u,1,0}\)设置为\(inf\)，即该状态不合法。对于\(rt\)，只有\(f_{rt,0,0}\)和\(f_{rt,1,1}\)合法。
染色和上面类似。
分治+树形DP即可
const int N=2e5+5;
#define ll long long
const ll inf=2e18;
inline void Min(ll &a,ll b){if(a>b)a=b;}
struct Edge{int to,next;}e[N<<1];
int head[N],ecnt=1;
inline void adde(int u,int v){e[++ecnt]=(Edge){v,head[u]};head[u]=ecnt;}
int inq[N],inde;
int n,m,k,a[N],root,dep[N];
void dfs1(int u,int pre){
for(int i=head[u],v;i&&!root;i=e[i].next){
v=e[i].to;if((i^1)==pre||inq[v]!=inde)continue;
if(!dep[v]) dep[v]=dep[u]+1,dfs1(v,i);
else if(dep[v]<dep[u])root=u;
}
}
inline int find_rt(int x){
if(m==n-1)return x;
else{root=0;dep[x]=1;dfs1(x,0);return root?root:x;}
}
inline ll update(ll x){return min(x,inf);}
ll f[N][2][2],ans[N],w[N];
inline ll Pow(ll x){return k==1?x:x*x;}
inline ll calc(int x,int y){return w[x]*Pow(abs(a[x]-y));}
bool vis[N];
void dfs(int u,int rt,int pre,int mid){
int flag=0;
f[u][0][0]=f[u][1][0]=calc(u,mid);
f[u][0][1]=f[u][1][1]=calc(u,mid+1);
vis[u]=1;
for(int i=head[u],v;i;i=e[i].next){
v=e[i].to;if(inq[v]!=inde||pre==(i^1))continue;
if(v==rt) flag=i;
else if(!vis[v]){
dfs(v,rt,i,mid);
bool cur=(i&1)^1;
for(int j=0;j<2;++j){
f[u][j][cur]=update(f[u][j][cur]+f[v][j][cur]);
f[u][j][cur^1]=update(f[u][j][cur^1]+min(f[v][j][cur],f[v][j][cur^1]));
}
}
}
if(flag){
if(flag&1)f[u][1][0]=inf;
else f[u][0][1]=inf;
}
if(u==rt) for(int i=0;i<2;++i)f[u][i][i^1]=inf;
}
bool col[N],coled[N];
void dfs2(int u,int pre,bool colu,bool c){
col[u]=c;coled[u]=1;
for(int i=head[u],v;i;i=e[i].next){
v=e[i].to;if((i^1)==pre||coled[v])continue;
if((i&1)!=c)dfs2(v,i,colu,c);
else dfs2(v,i,colu,f[v][colu][0]>f[v][colu][1]);
}
}
int q[N],q1[N],q2[N];
void work(int l,int r,int ql,int qr){
if(ql>qr)return;
if(l==r){for(int i=ql;i<=qr;++i)ans[q[i]]=l;return;}
int mid=(l+r)>>1;++inde;
for(int i=ql;i<=qr;++i)inq[q[i]]=inde,vis[q[i]]=0,coled[q[i]]=0,dep[q[i]]=0;
for(int i=ql;i<=qr;++i)if(!vis[q[i]]){
int rt=find_rt(q[i]);
dfs(rt,rt,0,mid);
if(f[rt][0][0]<f[rt][1][1])dfs2(rt,0,0,0);
else dfs2(rt,0,1,1);
}
int tot1=0,tot2=0;
for(int i=ql;i<=qr;++i)
if(!col[q[i]])q1[++tot1]=q[i];
else q2[++tot2]=q[i];
memcpy(q+ql,q1+1,tot1<<2); memcpy(q+ql+tot1,q2+1,tot2<<2);
work(l,mid,ql,ql+tot1-1); work(mid+1,r,ql+tot1,qr);
}
int main(){
n=read();m=read();k=read();
int mx=0;
for(int i=1;i<=n;++i)a[i]=read(),mx=max(mx,a[i]);
for(int i=1;i<=n;++i)w[i]=read();
for(int i=1,u,v;i<=m;++i){
u=read();v=read();bool val=readtype();
if(val)adde(u,v),adde(v,u);
else adde(v,u),adde(u,v);
}
for(int i=1;i<=n;++i)q[i]=i;
work(0,mx,1,n);
ll x=0;
for(int i=1;i<=n;++i)x+=calc(i,ans[i]);
printf("%lld\n",x);
return 0;
}
标签：ll,更优,int,回归,mid,问题,theta,保序,ql 来源： https://www.cnblogs.com/harryzhr/p/14891588.html