codevs 切糕 网络流

题目描述 Description

经过千辛万苦小 A 得到了一块切糕，切糕的形状是长方体，小 A 打算拦腰将切糕切成两半分给小 B。出于美观考虑，小 A 希望切面能尽量光滑且和谐。于是她找到你，希望你能帮她找出最好的切割方案。 
出于简便考虑，我们将切糕视作一个长 P、宽 Q、高 R 的长方体点阵。我们将位于第 z层中第 x 行、第 y 列上(1≤x≤P, 1≤y≤Q, 1≤z≤R)的点称为(x,y,z)，它有一个非负的不和谐值 v(x,y,z)。一个合法的切面满足以下两个条件： 
1. 与每个纵轴(一共有 P*Q 个纵轴)有且仅有一个交点。即切面是一个函数 f(x,y)，对于所有 1≤x≤P, 1≤y≤Q,我们需指定一个切割点 f(x,y),且 1≤f(x,y)≤R。 
2. 切面需要满足一定的光滑性要求，即相邻纵轴上的切割点不能相距太远。对于所有的 1≤x,x’≤P 和 1≤y,y’ ≤Q，若|x-x’|+|y-y’|=1，则|f(x,y)-f(x’,y’)| ≤D，其中 D 是给定的一个非负整数。 
可能有许多切面f 满足上面的条件，小A 希望找出总的切割点上的不和谐值最小的那个，即 ∑v(x,y, f(x,y))最小。

输入描述 Input Description

输入文件第一行是三个正整数P,Q,R，表示切糕的长P、宽Q、高R。第二行有一个非负整数D，表示光滑性要求。接下来是R个P行Q列的矩阵，第z个矩阵的第x行第y列是v(x,y,z) (1≤x≤P, 1≤y≤Q, 1≤z≤R)。

输出描述 Output Description

输出仅包含一个整数，表示在合法基础上最小的总不和谐值。

样例输入 Sample Input

input1
2 2 2 
1 
6 1 
6 1 
2 6 
2 6 
input2 
2 2 2 
0 
5 1 
5 1 
2 5 
2 5

样例输出 Sample Output

output1
6 
output2
12

数据范围及提示 Data Size & Hint

100%的数据满足P,Q,R≤40，0≤D≤R，且给出的所有的不和谐值不超过1000。 

var
    map:array[0..70000,-5..5]of longint;
    dis,his,pre:array[0..70000]of longint;
    vh:array[0..70000]of longint;
    fx:array[-5..5]of longint;
    a,b,c,d,flow:longint;
 
procedure init;
var
    i,j,k:longint;
begin
    read(b,c,a,d);
    for i:=1 to a do
      for j:=1 to b do
        for k:=1 to c do
          begin
            read(map[(i-1)*b*c+(j-1)*c+k,1]);
            if i>d then
            begin
              if k>1 then map[(i-1)*b*c+(j-1)*c+k,2]:=100000000;
              if j>1 then map[(i-1)*b*c+(j-1)*c+k,3]:=100000000;
              if k<c then map[(i-1)*b*c+(j-1)*c+k,4]:=100000000;
              if j<b then map[(i-1)*b*c+(j-1)*c+k,5]:=100000000;
            end;
          end;
    fx[1]:=b*c;
    fx[2]:=-d*b*c-1;
    fx[3]:=-d*b*c-c;
    fx[4]:=-d*b*c+1;
    fx[5]:=-d*b*c+c;
    for i:=1 to 5 do
      fx[-i]:=-fx[i];
end;
 
function max(x,y:longint):longint;
begin
    if x>y then exit(x);
    exit(y);
end;
 
procedure work;
var
    i,j,aug,min:longint;
    flag:boolean;
begin
    vh[0]:=a*b*c+2;
    i:=0;
    aug:=maxlongint;
    while dis[0]<=a*b*c+1 do
      begin
        flag:=false;
        his[i]:=aug;
        if i=0 then
          begin
            for j:=1 to b*c do
              if dis[0]=dis[j]+1 then
              begin
                flag:=true;
                pre[j]:=-1;
                break;
              end;
            if flag then i:=j;
          end
        else
          for j:=-5 to 5 do
            if i+fx[j]>0 then
            if (map[i,j]>0)and(dis[i]=dis[i+fx[j]]+1) then
            begin
              flag:=true;
              pre[i+fx[j]]:=-j;
              if aug>map[i,j] then aug:=map[i,j];
              inc(i,fx[j]);
              if i>a*b*c then
              begin
                inc(flow,aug);
                while i<>0 do
                  begin
                    inc(map[i,pre[i]],aug);
                    dec(map[max(i+fx[pre[i]],0),-pre[i]],aug);
                    inc(i,fx[pre[i]]);
                    if i<0 then i:=0;
                  end;
                aug:=maxlongint;
              end;
              break;
            end;
        if flag then continue;
        min:=a*b*c+1;
        if i=0 then
          begin
            for j:=1 to b*c do
              if min>dis[j] then min:=dis[j];
          end
        else
          for j:=-5 to 5 do
            if i+fx[j]>0 then
            if (map[i,j]>0)and(dis[i+fx[j]]<min) then min:=dis[i+fx[j]];
        dec(vh[dis[i]]);
        if vh[dis[i]]=0 then break;
        dis[i]:=min+1;
        inc(vh[dis[i]]);
        if i<>0 then
        begin
          inc(i,fx[pre[i]]);
          if i<0 then i:=0;
          aug:=his[i];
        end;
      end;
    write(flow);
end;
 
begin
    init;
    work;
end.