傳說中的斜率最佳化

來源:互聯網
上載者:User

對於形如f[i]由其之前若干f[j]決定的dp來說,最佳化大致有單調隊列,線段樹,四邊形不等式和斜率最佳化。
單調隊列一般是累加關係,可以提前直接預先計算後面的dp值,從而決策。而四邊形不等式則多用於第i個目標支配k到j這一段。至於乘除關係的則可以使用斜率最佳化。

當我們想知道兩個決策點哪一個更優時,可以將兩個決策點所產生的決策運算式寫出來。設為f[j1]和f[j2]。令j1<j2,則如果f[j1]<f[j2],那麼可以得到關於j1,j2和i的運算式。

由於對於f[i]的最優解為k時,所有j<k滿足slope[j,k]<=s[i],又因為s[i+1]為遞增的,所以slope[j,k]<s[i+1]。因此f[i+1]的決策k1>=k。我們可以看出這類dp的決策序列一般是遞增的。

那麼我們看維護的決策序列。f[i]在決策之時選取了j1,則slope[j1,j2]>s[i].之後對f[i+1]進行決策時,最優決策點向後移動,那麼就有slope[j1,j2]<=s[i+1].但是對於f[i+1]的最優取值集合來說,slope[j2,j3]>s[i+1]。顯然,slope[j2,j3]>slope[j1,j2],因此從這個角度也可以說明應該維護斜率遞增序列。

可用於斜率最佳化的運算式整理出後左邊必須形如(g1[j1]-g1[j2])/(g2[j1]-g2[j2])的形式。而且左邊必須可以不含帶i的項(否則無法維護隊列)。由於式子本身具有單調性,分母一般為負。如果為正,則後面維護隊列的時候取值標準也應該剛好相反。

由於我們想要找到(g1[j1]-g1[j2])/(g2[j1]-g2[j2])>S(i)的項,故維護隊列的時候滿足斜率遞增的順序。前面的項如果(g1[j1]-g1[j2])/(g2[j1]-g2[j2])<=S(i)則離開隊列,因為這表示j2比j1優。而之所以維護斜率上升,是為了找到一個突變點,在此點處,(g1[j1]-g1[j2])/(g2[j1]-g2[j2])>S(i)。說明j1優於j2。這時候j1為拋物線的最低點,也就是最優值的決策點。

斜率最佳化用的時候要自己推一下運算式,比傻瓜式地直接套用資料結構最佳化要好些。呵呵。

下面有兩題,一個是HNOI2008的題,這個是裸的,比較基礎。另一個是pku3709,北大月賽的題目,維護隊列的時候有個小技巧。

因為在斜率還沒有超過S(i)的時候,最優決策值一定是取最後面的一個。而pku3709在維護隊列的時候,不能算完了直接加進去。因為對於i+1來說,i是不能拿來更新的,而i的加入則可能會將本來能用來決策的j給去掉。

所以說,先不要放進去,需要它的時候再放。什麼時候需要呢?就是在k個位置之後了。所以說計算i之前先把i-k放進隊列,再進行相關計算和維護。

/*HNOI2008Toy*/

#include<iostream>
#include<cstdio>
#include<cstring>
using namespace std;
const int maxn = 50010;
int n;
long long l;
long long f[maxn],c[maxn];
int q[maxn];
inline bool slide(int x,int y,int z)
{
     long long sx = f[x]+(x+c[x])*(x+c[x]);
     long long sy = f[y]+(y+c[y])*(y+c[y]);
     long long sz = f[z]+(z+c[z])*(z+c[z]);
     if((sx-sy)*(y+c[y]-z-c[z])<(sy-sz)*(x+c[x]-y-c[y]))return true;
     return false;
}
int main()
{
   int i,j;
  // freopen("test2.txt","r",stdin);
  // freopen("out1.txt","w",stdout);
   while(scanf("%d%I64d",&n,&l)!=EOF)
   {
      c[0] = 0;
      for(i = 1;i<=n;i++)
        scanf("%I64d",&c[i]),c[i]+=c[i-1];
      int head = 0,tail = 0;
      q[tail++] = 0;      //0這個決策一定要放進去!!
      f[0] = 0;
      f[1] = (c[1]-l)*(c[1]-l);
      q[tail++] = 1;
      for(i = 2;i<=n;i++)
      {
         long long T = i+c[i]-l-1;
         while(1)
         {
            if(head+1>=tail)break;
            long long s = f[q[head]]+(q[head]+c[q[head]])*(q[head]+c[q[head]]);
            long long t = f[q[head+1]]+(q[head+1]+c[q[head+1]])*(q[head+1]+c[q[head+1]]);
            long long x = q[head]+c[q[head]]-q[head+1]-c[q[head+1]];
            if(s-t<2*T*x)break;
            head++;
         }
         if(head+1<tail)
         {
            f[i] = f[q[head]]+(i-q[head]-1+c[i]-c[q[head]]-l)*(i-q[head]-1+c[i]-c[q[head]]-l);
            while(head+1<tail&&!slide(q[tail-2],q[tail-1],i))tail--;
            q[tail++] = i;
         }
         else
         {
            if(f[q[head]]+(i-q[head]-1+c[i]-c[q[head]]-l)*(i-q[head]-1+c[i]-c[q[head]]-l)<=(c[i]+i-1-l)*(c[i]+i-1-l))
              f[i] = f[q[head]]+(i-q[head]-1+c[i]-c[q[head]]-l)*(i-q[head]-1+c[i]-c[q[head]]-l);
            else
            {
               head++;
               f[i] = (c[i]+i-1-l)*(c[i]+i-1-l);
            }
            q[tail++] = i;
         }
      }
      printf("%I64d\n",f[n]);
   }
   return 0;
}

/*pku3709*/

#include<iostream>
#include<cstdio>
#include<cstring>
using namespace std;
const int maxn = 500010;
int T,n,k;
long long f[maxn],a[maxn],sum[maxn];
int q[maxn];
inline bool slide(int x,int y,int z)
{
   long long xv = f[x]-sum[x]+x*a[x+1];
   long long yv = f[y]-sum[y]+y*a[y+1];
   long long zv = f[z]-sum[z]+z*a[z+1];
   if((xv-yv)*(a[y+1]-a[z+1])<(yv-zv)*(a[x+1]-a[y+1]))return true;
   return false;
}
int main()
{
    int i,j;
  //  freopen("test2.txt","r",stdin);
   // freopen("out1.txt","w",stdout);
    scanf("%d",&T);
    while(T--)
    {
        scanf("%d%d",&n,&k);
        sum[0] = 0;
        for(i = 1;i<=n;i++)
         scanf("%I64d",&a[i]),sum[i] = sum[i-1]+a[i];
        int head = 0,tail = 0;
        f[0] = 0;
        q[tail++] = 0;
        for(i = k;i<2*k;i++)
          f[i] =  sum[i]-i*a[1];
        for(i = 2*k;i<=n;i++)
        {
           while(head+1<tail&&!slide(q[tail-2],q[tail-1],i-k))tail--;
           q[tail++] = i-k;
        //   printf("i=%d tail = %d q[tail] = %d\n",i,tail-1,q[tail-1]);
           while(1)
           {
              if(head+1>=tail)break;
              long long x = f[q[head]]-sum[q[head]]+q[head]*a[q[head]+1];
              long long y = f[q[head+1]]-sum[q[head+1]]+q[head+1]*a[q[head+1]+1];
              long long z = a[q[head]+1]-a[q[head+1]+1];
              if(x-y<i*z)break;
              head++;
           }
       //    printf("q[head]=%d\n",q[head+1]);
           if(head+1<tail)
           {
              f[i] = f[q[head]]+sum[i]-sum[q[head]]-(i-q[head])*a[q[head]+1];
           }
           else
           {
              if(f[q[head]]+sum[i]-sum[q[head]]-(i-q[head])*a[q[head]+1]<sum[i]-i*a[1])
                f[i] = f[q[head]]+sum[i]-sum[q[head]]-(i-q[head])*a[q[head]+1];
              else
              {
                head++;
                f[i] = sum[i]-i*a[1];
              }
           }
        }
   /*     for(i = 1;i<=n;i++)
         printf("f[%d] = %lld\n",i,f[i]);*/
        printf("%I64d\n",f[n]);
    }
    return 0;
}

轉載自:http://blog.163.com/myq_952/blog/static/863906320112711750378/

聯繫我們

該頁面正文內容均來源於網絡整理,並不代表阿里雲官方的觀點,該頁面所提到的產品和服務也與阿里云無關,如果該頁面內容對您造成了困擾,歡迎寫郵件給我們,收到郵件我們將在5個工作日內處理。

如果您發現本社區中有涉嫌抄襲的內容,歡迎發送郵件至: info-contact@alibabacloud.com 進行舉報並提供相關證據,工作人員會在 5 個工作天內聯絡您,一經查實,本站將立刻刪除涉嫌侵權內容。

A Free Trial That Lets You Build Big!

Start building with 50+ products and up to 12 months usage for Elastic Compute Service

  • Sales Support

    1 on 1 presale consultation

  • After-Sales Support

    24/7 Technical Support 6 Free Tickets per Quarter Faster Response

  • Alibaba Cloud offers highly flexible support services tailored to meet your exact needs.