[HDOJ3359 Kind of a Blur]高斯消元、例题

【题目大意】

给定一个矩阵,A[i][j]=满足abs(i-ii)+abs(j-jj)<=d的T[ii][jj]的平均数。

让你还原成T[i,j]

【算法分析】

列方程,高斯消元。

【其它】

WA很多次,需要注意的是:消元时,那个rate必须是主元的那行做分母,否则有些情况下不行。。。

为什么会这样暂时不是很清楚,等下再YY一下。。。

【CODE】

#include #include #include #include #include using namespace std;

int m,n,d,cnt[122],pos[12][12];
double q[12][12],a[122][122],x[122];

bool input(){
scanf("%d%d%d",&n,&m,&d);
if (!n && !m && !d) return false;
for (int i=1;i<=m;i++)
for (int j=1;j<=n;j++)
scanf("%lf",&q[i][j]);
return true;
}

void init(){
memset(a,0,sizeof(a));
memset(cnt,0,sizeof(cnt));
memset(x,0,sizeof(x));
int i,j,ii,jj,tmp=1;
for (i=1;i<=m;i++)
for (j=1;j<=n;j++,tmp++)
pos[i][j]=tmp;
for (i=1;i<=m;i++)
for (j=1;j<=n;j++)
for (ii=1;ii<=m;ii++)
for (jj=1;jj<=n;jj++)
if (abs(ii-i)+abs(jj-j)<=d){
cnt[pos[i][j]]++;
a[pos[i][j]][pos[ii][jj]]=1.0;
}
for (i=1;i<=m;i++)
for (j=1;j<=n;j++)
a[pos[i][j]][m*n+1]=q[i][j]*cnt[pos[i][j]];
}

void guass(){
int i,j,k; double sum,rate;
for (k=1;k for (i=k,j=k;i<=n;i++) j=fabs(a[i][k])>fabs(a[j][k])?i:j;
for (i=k;i<=n+1;i++) swap(a[j][i],a[k][i]);
for (i=k+1;i<=n;i++)
for (rate=a[i][k]/a[k][k],j=k;j<=n+1;j++)
a[i][j]-=a[k][j]*rate;
/* if (a[i][k]!=0)
for (rate=a[k][k]/a[i][k],j=k;j<=n+1;j++)
a[i][j]=a[i][j]*rate-a[k][j];
这样写是错的!!!!!!!!!!!!
*/
}
for (i=n;i>=1;i–){
for (sum=0,j=i+1;j<=n;sum+=a[i][j]*x[j],j++);
x[i]=(a[i][n+1]-sum)/a[i][i];
}
}

void output(){
for (int i=1;i<=m;i++,printf("n"))
for (int j=1;j<=n;printf("%8.2lf",x[pos[i][j]]),j++);
}

int main(){
freopen("input.txt","r",stdin);
freopen("output.txt","w",stdout);
int Tc=0;
while (input()){
if (Tc) printf("n");
Tc++;
init();
n*=m;
guass();
n/=m;
output();
}
}

[APIO2008 roads]并差集、生成树、贪心**

【题目】http://www.apio.olympiad.org/

【算法分析】

使用构造法。

先假设所有的的水泥路都已经连好,然后利用鹅卵路进行合并块,做生成树。

然后用到的这些鹅卵路看作是必需的。

然后清空这个图,并只连上必需的鹅卵路。

一条一条的连鹅卵路使得到达K条。并且整个图没有环。

然后在这个图是用水泥路做生成树就可以了。

【CODE】

#include #include #include const int maxn=21111;
const int maxm=101111;
int n,m,k,haveuse;
int lx[maxm],ly[maxm],lt[maxm],fa[maxn];
bool use[maxm];//必须用的边

void input(){
     scanf("%d%d%d",&n,&m,&k);
     for (int i=1;i<=m;i++) scanf("%d%d%d",&lx[i],&ly[i],&lt[i]);
}

int gf(int k){
    if (fa[k]==k) return k;
    fa[k]=gf(fa[k]);
    return fa[k];
}

void make_need(){
     memset(use,false,sizeof(use));
     for (int i=1;i<=n;i++) fa[i]=i;
     for (int i=1;i<=m;i++)
       if (lt[i] && gf(lx[i])!=gf(ly[i])) fa[fa[lx[i]]]=fa[ly[i]];
     for (int i=1;i<=m;i++)
       if (!lt[i] && gf(lx[i])!=gf(ly[i])){
         fa[fa[lx[i]]]=fa[ly[i]];
         use[i]=true;
         haveuse++;
       }
}

void add_egg_edge(){
     for (int i=1;i<=n;i++) fa[i]=i;
     for (int i=1;i<=m;i++)
       if (use[i]) fa[gf(lx[i])]=gf(ly[i]);
     for (int i=1;i<=m;i++)
       if (haveuse         use[i]=true;
         fa[fa[lx[i]]]=fa[ly[i]];
         haveuse++;
       }
}

void add_water_muddy_edge(){
     for (int i=1;i<=m;i++)
       if (lt[i] && gf(lx[i])!=gf(ly[i])){
         fa[fa[lx[i]]]=fa[ly[i]];
         use[i]=true;
       }
}

void output(){
     for (int i=1;i<=m;i++)
       if (use[i]) printf("%d %d %dn",lx[i],ly[i],lt[i]);
}

int main(){
    freopen("roads.in","r",stdin);
    freopen("roads.out","w",stdout);
    input();
    make_need();
    add_egg_edge();
    if (haveuse!=k) puts("no solutionn");
    else{
         add_water_muddy_edge();
         output();
    }
}

[APIO2008 dna]动态规划、统计方案数

【题目】http://www.apio.olympiad.org/

【算法分析】

F[i,j,k]表示i..m这个后缀串中,s[i]==j,且有k个逆序的位置的方案数是多少。

然后Sum[i,j,k]表示F[i,j,0]+F[i,j,1]+….+F[i,j,k]

转移很简单,就不说了。

然后从1开始搞下来,利用Sum来从前面开始判定下来处理一下就好。

【CODE】

#include #include #include #include using namespace std;
typedef long long big;
big n,K,R;
int s[51111],ans[51111];
big F[51111][5][11],Sum[51111][5][11];

void input(){
     cin >> n >> K >> R;
     char c;
     for (int i=1;i<=n;i++){
         for (c=getchar();c==’n’ || c==’ ‘;c=getchar());
         switch (c){
                case ‘A’:
                     s[i]=1;
                     break;
                case ‘C’:
                     s[i]=2;
                     break;
                case ‘G’:
                     s[i]=3;
                     break;
                case ‘T’:
                     s[i]=4;
                     break;
                case ‘N’:
                     s[i]=0;
                     break;
         }
     }
     if (!s[n])
       for (int j=0;j<5;j++)
         F[n][j][0]=1;
     else
       F[n][s[n]][0]=1;
}

void dp(){
     int i,j,k,l;
     for (i=n-1;i>=1;i–)
       if (!s[i])
         for (j=1;j<5;j++)
           for (k=0;k             for (l=1;l<5;l++)
               if (l                 F[i][j][k+1]+=F[i+1][l][k];
               else
                 F[i][j][k]+=F[i+1][l][k];
       else{
         j=s[i];
           for (k=0;k             for (l=1;l<5;l++)
               if (l                 F[i][j][k+1]+=F[i+1][l][k];
               else
                 F[i][j][k]+=F[i+1][l][k];
       }
       for (i=1;i<=n;i++)
         for (j=1;j<5;j++)
           for (Sum[i][j][0]=F[i][j][0],k=1;k             Sum[i][j][k]=Sum[i][j][k-1]+F[i][j][k];
}

void output(){
     int i,j,k,kn=K;
     big ns,ss=0;
     for (i=1;i<=n;i++){
         kn–;
         for (j=1;j<5;j++){
             if (j==ans[i-1]) kn++;
             if (ss+Sum[i][j][kn]>=R) break;
             else ss+=Sum[i][j][kn];
         }
         ans[i]=j;
     }
//========================================================================
     char c;
     for (i=1;i<=n;i++){
         switch (ans[i]){
                case 1:c=’A’; break;
                case 2:c=’C’; break;
                case 3:c=’G’; break;
                case 4:c=’T’; break;
         }
         putchar(c);
     }
}

int main(){
    freopen("dna.in","r",stdin);
    freopen("dna.out","w",stdout);
    input();   
    dp();
    output();
}

[APIO2007 zoo]状态压缩动态规划

【题目】http://www.oi.bbjy.com/n38c10.aspx

【算法分析】

这个比较好想,就是F[i][j]表示到达第i个点,状态为j的时候最多让多少个boy满意。

j是一个4位的2进制数,表示前面4个位置是搬走还是留下。在我的程序中,1表示搬走,0表示留下。

然后由于一个圈回来的时候连接的部分有些麻烦,所以直接枚举开始的4个位,然后再转移。

另外转移的部分比较纠结,我一开始用5个指针一起线性遍历,最慢的那个点变成3s+…..

然后后来围观了野牛的,发现他的快得很变态…..

呃,原来他用得位运算判重….

但是悲剧的是,我这个改了以后还是慢他的一倍

【其他】时间比较惨剧,在我学校的电脑里,不开任何东西,最慢那个点1.6s,开了-O999,变成0.8s

我也不知道为什么可以开-O999这东西。。。但是开-O999确实比开-O2快= =,那位神牛来解释一下这是啥?

【CODE】

#include #include const int N=51111;
const int INF=0x7FFFFFFF;
struct gtp{int y,next;}g[1001111];
int n,c,e,ans,hs,ls,times;
int p[N],love[N],hate[N],Hls[N],Lls[N],list[N],v[N];
int F[N][16];

inline int C(int x)          {if (x<0) x+=n; if (x>n) x-=n; return x;}
inline void addH(int x,int y) {e++; g[e].y=y; g[e].next=Hls[x]; Hls[x]=e;}
inline void addL(int x,int y) {e++; g[e].y=y; g[e].next=Lls[x]; Lls[x]=e;}
inline void max(int &x,int y) {if (y>x) x=y;}

void input(){
     scanf("%d%d",&n,&c);
     for (int i=1;i<=c;i++){
         int t1,t2,x;
         scanf("%d%d%d",&p[i],&t1,&t2);
         p[i]=C(p[i]+4);
         for (int j=1;j<=t1;j++){
             scanf("%d",&x);
             addH(x,i);
             hate[i]|=1<         }
         for (int j=1;j<=t2;j++){
             scanf("%d",&x);
             addL(x,i);
             love[i]|=1<         }
     }
}

int Cal(int i,int ns){
    int t,res;
    res=0;
    ns&=31;
    if (ns&1) t=Hls[i];
         else t=Lls[i];
    ns>>=1;
    for (;t;t=g[t].next){
        hs=hate[g[t].y]>>(C(p[g[t].y]-i)+1);
        ls=love[g[t].y]>>(C(p[g[t].y]-i)+1);
        if (hs&ns || ls&(INF^ns)) continue;
        res++;
    }
    return res;
}

void dp(int st){
     int i,ns,now,tmp;
     for (i=1;i<=n;i++)    
       for (ns=0;ns<16;ns++)
         F[i][ns]=-INF;
     F[5][(st<<1)&15]=Cal(5,st<<1);
     F[5][((st<<1)+1)&15]=Cal(5,(st<<1)+1);
     for (i=5;i       for (ns=0;ns<16;ns++)
         if (F[i][ns]!=-INF){
           max(F[i+1][(ns<<1)&15],F[i][ns]+Cal(i+1,ns<<1));
           max(F[i+1][((ns<<1)+1)&15],F[i][ns]+Cal(i+1,(ns<<1)+1));
         }
     for (ns=0;ns<16;ns++)
       if (F[n][ns]!=-INF){
         now=F[n][ns];
         now+=Cal(1,(ns<<1)+(st>>3));
         now+=Cal(2,(ns<<2)+(st>>2));
         now+=Cal(3,(ns<<3)+(st>>1));
         now+=Cal(4,(ns<<4)+(st>>0));
         ans>?=now;
       }
}

void count(){
     int now=0,i,res=0,t;
     times++;
     for (i=1;i<=n;i++)
       if (list[i]==0)
         for (t=Lls[i];t;t=g[t].next) v[g[t].y]=times;
       else
         for (t=Hls[i];t;t=g[t].next) v[g[t].y]=times;
     for (i=1;i<=c;i++)
       if (v[i]==times) res++;
     ans>?=res;
}

void dfs(int i){
     if (i>n){
       count();
       return;
     }
     list[i]=0;
     dfs(i+1);
     list[i]=1;
     dfs(i+1);
}

void work(){
     if (n>=5)
       for (int st=0;st<16;st++) dp(st);
     else dfs(1);
     printf("%dn",ans);
}

int main(){
    freopen("zoo.in","r",stdin);
    freopen("zoo.out","w",stdout);
    input();
    work();
}

[APIO2007 backup]用贪心实现特殊图的费用流增广**

【题目】http://www.oi.bbjy.com/n38c10.aspx

【算法分析】

O(nk)的dp。测得70分….

然后我YY到一个费用流,发现复杂度还是O(nk),而且常数上必然比dp更慢。。。

于是翻看Sinya神牛的题解。。。Orz,原来是用贪心去做费用流。。。彻底膜拜。

http://sinyalee.com/blog/?p=427

拜一下先….

然后我YY到的费用流就和野牛的建图一样.

下面我就说一下我对那个贪心的理解.

首先我们想象我们是一只小蝌蚪,在流中寻找增广路….

(借一下野牛的图。。。)

然后可以YY到,我们找的增广路必然是从与S相连的一条未增广的边进去,然后向左或者向右若干次“上下波动”(当然也可以不上下动,直接进去T那里),然后跑到T点去。

如在a2->a3->a4这里波动。

我们如果增广路包含这一路径的话,就可以看成把a3给取消了,然后选择了a2和a4.
这样的话,连续的(不取,取,不取,取,不取….)就可以当成一个组。

注意最左边和最右边都必须为“不取”。

分成那个组以后,组的中间必然是不能做增广的,因为与S相连的边已经饱和~

所以就可以看成必须在两端搞进去。。。搞了以后,无论在左边开始搞还是右边开始搞,结果都是一样的,都是多了一个点,都是向两边迈进了一步。(画图YY一下即可)

然后就再把相邻的两组合并到一起。因为你和她已经弄成新的OXOXOXOXOXOXO….结构了

现在大体我们已经明了了,但是边沿情况怎么处理呢?

注意到,如果a1或者an已经选了,是不可能再取消的。那就别放进决策队列。

证明:

你YY一下费用流的增广过程,

会发现那里要么和S连的只属于一个点的边已经满流,要么和T连的只属于一个店的边已经满流,你是没有办法取消它了。。。

因为费用流是对的,所以这个是对的,因为这只是在模拟增广费用流。

然后综上一下,发现堆可以达到我们需要的每次找最短增广路径的目的。

增广增加的费用其实就是SUM[L[i],R[i]]-2*w[i]。

其中sum[i,j]表示a[i]到a[j]的和,w[i]表示这个区间已经选了的点的权值和。

根据这个做堆的关键字就可以了。

【其他】膜拜Sinya

【CODE】

#include #include #include #include using namespace std;
const int N=101111;
int n,K;
int a[N],inheap[N],fa[N],L[N],R[N];
int sum[N],ans;

inline int S(int i,int j){
       if (!i) return sum[j];
       return sum[j]-sum[i-1];
}

struct HeapType{
       struct pt{int pos;int w;}h[N];
       int tot;
      
       bool cmp(int i,int j){
            return (S(L[h[i].pos],R[h[i].pos])-h[i].w*2)
                  <(S(L[h[j].pos],R[h[j].pos])-h[j].w*2);
       }
      
       void swap(int i,int j){
            int tmp;
            tmp=h[i].pos; h[i].pos=h[j].pos; h[j].pos=tmp;
            tmp=h[i].w; h[i].w=h[j].w; h[j].w=tmp;
            inheap[h[i].pos]=i; inheap[h[j].pos]=j;
       }
      
       void up(int k){
            while (k>1 && cmp(k,k/2)){
                  swap(k/2,k);
                  k/=2;
            }
       }
      
       void down(int k){
            int p;
            while (k*2<=tot){
                  p=k*2;
                  if (p+1<=tot && cmp(p+1,p)) p++;
                  if (cmp(p,k)){
                    swap(p,k);
                    k=p;
                  }else break;
            }
       }
      
       void ins(int pos,int w){
            tot++;
            h[tot].pos=pos;
            h[tot].w=w;
            inheap[pos]=tot;
            up(tot);    
       }
      
       void del(int inh){
            if (!inh) return;
            int tmp=h[inh].pos;
            swap(tot,inh);
            inheap[tmp]=0;
            tot–;
            up(inh);
            down(inh);
       }
}heap;

void input(){
     cin >> n >> K;
     int x;
     for (int i=1;i<=n;i++){
         scanf("%d",&x);
         a[i]=x;
     }
     for (int i=1;i     n–;
     sum[0]=0;
     for (int i=1;i<=n+1;i++){
         L[i]=i; R[i]=i; fa[i]=i;
         sum[i]=sum[i-1]+a[i];
     }
     L[0]=0; R[0]=0; fa[0]=0;
     heap.tot=0;
     for (int i=1;i<=n;i++) heap.ins(i,0);
}

int gf(int k){
    if (fa[k]==k) return k;
    fa[k]=gf(fa[k]);
    return fa[k];
}

void work(){    
     int i,k,ans=0,Lc,Rc,p,nw;
     for (k=1;k<=K;k++){
         p=gf(heap.h[1].pos);
         ans+=S(L[p],R[p])-heap.h[1].w*2;
         nw=S(L[p],R[p])-heap.h[1].w;
         heap.del(1);
        
         Lc=gf(L[p]-1);
         nw+=heap.h[inheap[Lc]].w;
         L[p]=L[Lc];
         fa[Lc]=p;
         heap.del(inheap[Lc]);
        
         Rc=gf(R[p]+1);
         nw+=heap.h[inheap[Rc]].w;
         R[p]=R[Rc];
         fa[Rc]=p;
         heap.del(inheap[Rc]);
        
         if (L[p]>0 && R[p]<=n) heap.ins(p,nw);
     }
     cout << ans << endl;
}

int main(){
    freopen("backup.in","r",stdin);
    freopen("backup.out","w",stdout);
    input();
    work();
    return 0;
}