stepu.pas

运筹学演示程序之运输问题

unit StepU;

interface
uses
  Windows,Messages,SysUtils,Classes,StdCtrls,Controls,Forms,
  Grids,Dialogs,ExtCtrls,Graphics,BasicU;
const
   size=30;
var
   stack:array[1..size,1..4] of integer;
   test : matrix;

function step_init(pc,sc:byte;var amount:matrix):boolean;
                                          //返回值为是否出现退化
function step_init2(pc,sc:byte;var amount:matrix;cs:matrix):boolean;
                                          //返回值为是否出现退化
function prj_step(pc,sc:byte;produce,sale:array of real;
      cost:matrix;var amount:matrix):real;
function draw_step(pc,sc:byte;produce,sale:array of real;
      cost:matrix;var amount:matrix;
      var grid1:TStringGrid;var grid2:tstringgrid;ptime:dword;
      var edt:Tedit):boolean;

implementation

type
  stack_of_point=array[1..size,1..2] of integer;
var
  ar1,ar2,st:stack_of_point;

FUNCTION step_init2(pc,sc:byte;var amount:matrix;cs:matrix):boolean;
var tuihua:boolean;
    c,i,j:integer;
    mi,mj:integer;
  procedure stack_init(var stk:stack_of_point);
  var i:integer;
  begin
    for i:=1 to size do
      stk[i,1]:=-1;
  end;
  function push_stack(var stk:stack_of_point;x,y:byte):boolean;
  var i:integer;
  begin
    i:=0;
    repeat
      i:=i+1;
    until (i>=size)or(stk[i,1]=-1);
    if (stk[i,1]<>-1) then result:=false
    else
    begin
      stk[i,1]:=x;
      stk[i,2]:=y;
      result:=true;
    end;
  end;
  function in_stack(stk:stack_of_point;x,y:byte):boolean;
  var i:integer;
  begin
    i:=0;
    repeat
      i:=i+1;
    until (i>=size)or((stk[i,1]=x)and(stk[i,2]=y));
    if ((stk[i,1]=x)and(stk[i,2]=y)) then result:=true
    else result:=false;
  end;
  function searchone(var stk:stack_of_point):boolean;
  var i,j,x,y:integer;
      u,v:array of boolean;
      change:boolean;
  begin
    x:=-1;  y:=-1;
    for i:=0 to pc-1 do
      for j:=0 to sc-1 do
        if (amount[i,j]>=0) then
        begin
          if (not in_stack(ar1,i,j))and(not in_stack(ar2,i,j)) then
          begin
            x:=i;
            y:=j;
          end;
        end;
    if (x>-1) then
    begin
      //begin search
      setlength(u,pc);
      setlength(v,sc);
      for i:=0 to pc-1 do u[i]:=false;
      for j:=0 to sc-1 do v[j]:=false;
      u[x]:=true;
      v[y]:=true;
      repeat     //寻找一个分部。
        change:=false;
        for i:=0 to pc-1 do
          if u[i] then
          begin
            for j:=0 to sc-1 do
              if ((amount[i,j]>=0)and(not v[j])) then
              begin
                v[j]:=true;
                change:=true;
              end;
          end;
        for j:=0 to sc-1 do
          if v[j] then
          begin
            for i:=0 to pc-1 do
              if ((amount[i,j]>=0)and(not u[i])) then
              begin
                u[i]:=true;
                change:=true;
              end;
          end;
      until (not change);
      for i:=0 to pc-1 do     //找到的点压入堆栈。
        for j:=0 to sc-1 do
        if (u[i])and(v[j])and(amount[i,j]>=0) then
          push_stack(stk,i,j);
      result:=true;
    end
    else result:=false;
  end;
  procedure getpoints;
  var u1,v1,u2,v2:array of boolean;
      i,j:integer;
  begin
    setlength(u1,pc);
    setlength(v1,sc);
    setlength(u2,pc);
    setlength(v2,sc);
    for i:=0 to pc-1 do u1[i]:=false;
    for j:=0 to sc-1 do v1[j]:=false;
    for i:=0 to pc-1 do u2[i]:=false;
    for j:=0 to sc-1 do v2[j]:=false;
    i:=0;
    repeat
      i:=i+1;
      if (ar1[i,1]<>-1) then
      begin
        u1[ar1[i,1]]:=true;
        v1[ar1[i,2]]:=true;
      end;
    until ((i>=size)or(ar1[i,1]=-1));
    i:=0;
    repeat
      i:=i+1;
      if (ar2[i,1]<>-1) then
      begin
        u2[ar2[i,1]]:=true;
        v2[ar2[i,2]]:=true;
      end;
    until ((i>=size)or(ar2[i,1]=-1));
    for i:=0 to pc-1 do     //找到的点压入堆栈。
      if (u1[i]) then
      begin
        for j:=0 to sc-1 do
        if (v2[j])and(amount[i,j]<0) then
          push_stack(st,i,j);
      end;
    for i:=0 to pc-1 do
      if (u2[i]) then
      begin
        for j:=0 to sc-1 do
        if (v1[j])and(amount[i,j]<0) then
          push_stack(st,i,j);
      end;
  end;
  function find_min(var x:integer;var y:integer):boolean;
  var i,j,m:integer;
      min:real;
      sign:array of array of boolean;
      zero:boolean;
  begin
    // 用sign 来标志矩阵中哪些点是虚拟点所对应的运费点
    // true 表示为实变量
    setlength(sign,pc,sc);
    for i:=0 to pc-1 do
      for j:=0 to sc-1 do
        sign[i,j]:=true;
    for i:=0 to pc-1 do
    begin
      zero:=true;
      for j:=0 to sc-1 do if (amount[i,j]<>0) then zero:=false;
      if zero then
        for j:=0 to sc-1 do sign[i,j]:=false;
    end;
    for j:=0 to sc-1 do
    begin
      zero:=true;
      for i:=0 to pc-1 do if (amount[i,j]<>0) then zero:=false;
      if zero then
        for i:=0 to c-1 do sign[i,j]:=false;
    end;
    // 标志完毕
    // 优先选择实变量
    min:=-1;
    i:=-1;
    j:=-1;
    m:=0;
    repeat
      m:=m+1;
      if (st[m,1]<>-1) then
        if (sign[st[m,1],st[m,2]])and((min=-1)or(cs[st[m,1],st[m,2]]<min)) then
        begin
          min:=cs[st[m,1],st[m,2]];
          i:=st[m,1];
          j:=st[m,2];
        end;
    until ((m>=size)or(st[m,1]=-1));
    if (min<>-1) then
    begin
      x:=i;
      y:=j;
      result:=true;
    end
    else           // 如果没有实变量可供选择...那只好拿虚拟点来凑数了。
    begin          // 重新再找一遍，不考虑实与需
      m:=0;
      repeat
        m:=m+1;
        if (st[m,1]<>-1) then
          if ((min=-1)or(cs[st[m,1],st[m,2]]<min)) then
          begin
            min:=cs[st[m,1],st[m,2]];
            i:=st[m,1];
            j:=st[m,2];
          end;
      until ((m>=size)or(st[m,1]=-1));
      if (min<>-1) then
      begin
        x:=i;
        y:=j;
        result:=true;
      end
      else
        result:=false;
    end;
  end;
BEGIN
  c:=0;
  for i:=0 to pc-1 do
    for j:=0 to sc-1 do
      begin
        if (amount[i,j]=0) then amount[i,j]:=-1 //设置-1以作为区别
        else c:=c+1;
      end;
  tuihua:=false;
  if (c<(pc+sc-1)) then tuihua:=true;
  if (tuihua) then
  begin
    messagedlg('初始可行解中出现了退化问题。'+chr(13)
             +'先对初始可行解进行处理，补充基变量。',mtinformation,
             [mbok],0);
    repeat
      stack_init(ar1);
      stack_init(ar2);
      stack_init(st);
      searchone(ar1);
      searchone(ar2);
      getpoints;
      if find_min(mi,mj) then
      begin
        amount[mi,mj]:=0;
        c:=c+1;
      end;
    until (c>=(pc+sc-1));
  end;
  result:=tuihua;
END;

FUNCTION step_init(pc,sc:byte;var amount:matrix):boolean;
var tuihua,getpoint:boolean;
    i,j,c,x,y:byte;
    u,v : array of real;
    ub,vb : array of boolean;
  procedure check_u_v;
  var i,j:byte;
      notgetone:boolean;
  begin
    for i:=0 to pc-1 do ub[i]:=false;
    for j:=0 to sc-1 do vb[j]:=false;
    for i:=0 to pc-1 do u[i]:=-1;
    for j:=0 to sc-1 do v[j]:=-1;
    if (random(2)=1) then
    begin
      i:=random(pc);
      ub[i]:=true;
      u[i]:=0;
    end
    else
    begin
      j:=random(sc);
      vb[j]:=true;
      v[j]:=0;
    end;
    repeat
      notgetone:=true;
      for j:=0 to sc-1 do
        if (vb[j]) then
          for i:=0 to pc-1 do
            if ((not ub[i])and(amount[i,j]>=0)) then
            begin
              u[i]:=amount[i,j]-v[j];
              ub[i]:=true;
              notgetone:=false;
            end;
      for i:=0 to pc-1 do
        if (ub[i]) then
          for j:=0 to sc-1 do
            if ((not vb[j])and(amount[i,j]>=0)) then
            begin
              v[j]:=amount[i,j]-u[i];
              vb[j]:=true;
              notgetone:=false;
            end;
    until notgetone;
  end;
BEGIN
  setlength(u,pc);
  setlength(v,sc);
  setlength(ub,pc);
  setlength(vb,sc);
  c:=0;
  for i:=0 to pc-1 do
    for j:=0 to sc-1 do
      begin
        if (amount[i,j]=0) then amount[i,j]:=-1 //设置-1以作为区别
        else c:=c+1;
      end;
  tuihua:=false;
  if (c<(pc+sc-1)) then tuihua:=true;
  if (tuihua) then
  begin
   messagedlg('初始可行解中出现了退化问题。'+chr(13)
             +'先对初始可行解进行处理，补充基变量。',mtinformation,
             [mbok],0);
   repeat
    check_u_v;
    getpoint:=false;
    i:=0;  j:=0;  x:=0;  y:=0;
    while ((i<=pc-1)and(ub[i])) do i:=i+1;
    if (ub[i]) then
      while ((j<=sc-1)and(vb[j])) do j:=j+1;
    if (not ub[i]) then
    begin
      for j:=0 to sc-1 do
        if (vb[j]) then
        begin
          x:=i; y:=j;
          getpoint:=true;
        end;
    end
    else if (not vb[j]) then
    begin
      for i:=0 to pc-1 do
        if (ub[i]) then
        begin
          x:=i; y:=j;
          getpoint:=true;
        end;
    end;
    if (not getpoint) then
    begin
      repeat
        i:=random(pc);
        j:=random(sc);
      until (((ub[i])and(not vb[j]))or((not ub[i])and(vb[j])));
      x:=i;   y:=j;
    end;
    if (amount[x,y]<0) then
    begin
      amount[x,y]:=0;
      c:=c+1;
    end
    else
    begin
      messagedlg('位势法初始化失败，请重试。',mterror,[mbok],0);
      result:=false;
      exit;
    end;
   until (c>=(pc+sc-1));
  end;
  result:=tuihua;
END;

FUNCTION prj_step(pc,sc:byte;produce,sale:array of real;
      cost:matrix;var amount:matrix):real;
var
  u,v : array of real;
  ub,vb : array of boolean;
  x,y:byte;
  do_fail:boolean;
  all_done:boolean;
  rc:real;
  function u_v_counted:boolean;
  var d:boolean;i:byte;
  begin
    d:=true;
    for i:=0 to pc-1 do if (not ub[i]) then d:=false;
    for i:=0 to sc-1 do if (not vb[i]) then d:=false;
    result:=d;
  end;
  procedure check_u_v;
  var i,j:byte;
  begin
    for i:=0 to pc-1 do ub[i]:=false;
    for j:=0 to sc-1 do vb[j]:=false;
    if (random(2)=1) then
    begin
      i:=random(pc);
      ub[i]:=true;
      u[i]:=0;
    end
    else
    begin
      i:=random(sc);
      vb[i]:=true;
      v[i]:=0;
    end;
    repeat
      for j:=0 to sc-1 do
        if (vb[j]) then
          for i:=0 to pc-1 do
            if ((not ub[i])and(amount[i,j]>=0)) then
            begin
              u[i]:=cost[i,j]-v[j];
              ub[i]:=true;
            end;
      for i:=0 to pc-1 do
        if (ub[i]) then
          for j:=0 to sc-1 do
            if ((not vb[j])and(amount[i,j]>=0)) then
            begin
              v[j]:=cost[i,j]-u[i];
              vb[j]:=true;
            end;
    until u_v_counted;
  end;
  procedure do_test;
  var i,j:byte;
  begin
    for i:=0 to pc-1 do
      for j:=0 to sc-1 do
        test[i,j]:=u[i]+v[j];
  end;
  function find_x_y(var x:byte;var y:byte):boolean;
  var i,j:byte; min,k:real; g:boolean;
  begin
    min:=0; g:=false;
    for i:=0 to pc-1 do
      for j:=0 to sc-1 do
      begin
        k:=cost[i,j]-test[i,j];
        if (k<min) then
        begin
          min:=k;
          x:=i;
          y:=j;
          g:=true;
        end;
      end;
    result:=g;
  end;
  function north(var x:byte;var y:byte):boolean;
  var i:byte;  get:boolean;
  begin
    if (x<=0) then
    begin
      result:=false;
      exit;
    end;
    i:=x-1;
    get:=false;
    while ((i>0)and(amount[i,y]<0)) do
      i:=i-1;
    if(amount[i,y]>=0) then
    begin
      x:=i;
      get:=true;
    end;
    result:=get;
  end;
  function south(var x:byte;var y:byte):boolean;
  var i:byte;  get:boolean;
  begin