type
  str60  = string[60];
  str30  = string[30];
  all_char = set of ' '..'}';
var
   d  : array[1..50,1..50] of real;  { Matrix to solve }
   b  : array[1..40,1..6]  of real;  { Beam Data }
   c  : array[1..40,1..2]  of real;  { column data }
   e  : array[1..10,1..4]  of real;  { wind loading }
   ld : array[1..50,1..4]  of real;  { loads, story, load, location, length }
   mp : array[1..20,1..2]  of integer; { Macro pointer of loadings }
   r  : array[0..200] of real;  { array of moments, deflections }

   s,p   : integer;  { # stories, # portals }
   q     : integer;  { size of matrix to solve }
   z     : integer;  { type of base joints 1 or 2 }
   l2,l3 : real;   { total height, total width feet }

   lp    : integer;  { load pointer, points to next row in load matrix }
   mc    : integer;  { macro counter, multible loads per macro }
   lc    : integer;  { load count of loads applied }
   lm    : array[1..50] of integer; { save load#'s of all loads applied }
   psl   : integer;  { plot scale member length }
   psm   : integer;  { plot scale moment kip-ft/inch}
   proj_descrip : str60;
   load_descrip : str60;
   in_file      : text;
   out_file     : text;
   result_file  : str30;
   time         : str30;
   option       : integer;
   shear        : boolean;
   response     : char;
   plot_file    : text;

const
   mf_plt  = 1016.00; { unit/inch on the plotter }
   term    = ^C;      { plotter label terminator }

procedure quest(x,y :integer; q_line:str60; q_set:all_char);
  begin
    if x*y>0 then gotoxy(x,y);
    clreol;
    write(q_line);
    repeat
      write(^G);
      read(kbd,response);
      response:=upcase(response);
    until response in q_set;
    write(response);
  end;

function parse_up(s:str60):str60;
  var i: integer;
  begin
    for i:=1 to length(s) do s[i]:=upcase(s[i]);
    parse_up:=s;
  end;

function exist(filename:str30):boolean;
var
 fil : file;
 i   : integer;
begin
  assign(fil,filename);
  {$I-} reset(fil); {$I+}
  i:=ioresult;
  if i=0 then close(fil);
  exist:=(I=0);
end;

function  when:str30;
  Type
    RegisterSet=Record Case Integer Of
                  1: (AX,BX,CX,DX,BP,SI,DI,DS,ES,Flags: Integer);
                  2: (AL,AH,BL,BH,CL,CH,DL,DH: Byte);
                End;
    Time=Record
           Hours,Min,Sec,Hundreths: Byte;
         End;
    DOW=(Sun,Mon,Tue,Wed,Thu,Fri,Sat);
    Date=Record
           Month,Day: Byte;
           Year: Integer;
           DayOfWeek: DOW;
         End;
    String2=String[2];
  Const
    DayName: Array [DOW] Of String[9]=
      ('Sun','Mon','Tue','Wed','Thu','Fri','Sat');
    MonName: Array [1..12] Of String[9]=
      ('Jan','Feb','Mar','Apr','May','Jun','Jul',
       'Aug','Sep','Oct','Nov','Dec');

  Var
    Regs: RegisterSet;
    i:string[4];

  Procedure GetTime(Var T:Time);

    Begin
      With Regs,T Do
       Begin
        AH:=$2C;
        MsDos(Regs);
        Hours:=CH;
        Min:=CL;
        Sec:=DH;
        Hundreths:=DL;
       End;
    End;

  Procedure GetDate(Var D:Date);

    Begin
      With Regs,D Do
       Begin
        AH:=$2A;
        MsDos(Regs);
        Month:=DH;
        Day:=DL;
        Year:=CX;
        DayOfWeek:=DOW(AL);
       End;
    End;

 Function V2(I:Integer):String2;
    Begin
      V2:=Chr(48+I Div 10)+Chr(48+I Mod 10);
    End;

  Var
    T1: Time;
    D1: Date;
    AP: String2;

  Begin
    GetTime(T1);
    GetDate(D1);
    AP:='am';
    If t1.Hours>11 Then AP:='pm';
    t1.Hours:=t1.Hours Mod 12;
    If t1.Hours=0 Then t1.Hours:=12;
    str(d1.year:4,i);
    when:=v2(t1.hours)+':'+v2(t1.min)+ap+'  '+
      dayname[d1.dayofweek]+', '+monname[d1.month]+' '+v2(d1.day)+', '+i;
  End;

function pw(n:real; p:integer):real;
  var i:integer;
      a:real;
  begin
    a:=1;
    for i:=1 to p do a:=a*n;
    pw:=a;
  end;

procedure readr( var r:real);
  var
    i,x,y : integer;
  begin
    x:=wherex; y:=wherey;
    repeat
      if paramcount=0 then gotoxy(x,y-1);
      {$I-} read(in_file,r); {$I+}
      i:=IOResult;
      if (I<>0) and (paramcount>0) then readln(in_file)
      else if (i<>0) and (paramcount=0) then
             begin gotoxy(x,y); clreol; end
    until i=0;
    if paramcount=0 then gotoxy(x,y);
  end;

procedure readi( var r:integer);
  var
    i,x,y : integer;
  begin
    x:=wherex; y:=wherey;
    repeat
      if paramcount=0 then gotoxy(x,y-1);
      {$I-} read(in_file,r); {$I+}
      i:=IOResult;
      if (I<>0) and (paramcount>0) then readln(in_file)
      else if (i<>0) and (paramcount=0) then
             begin gotoxy(x,y); clreol; end
    until i=0;
    if paramcount=0 then gotoxy(x,y);
  end;

procedure reads(var s:str60);
  var
    x,y:integer;
  begin
    x:=wherex; y:=wherey;
    if paramcount=0 then gotoxy(x,y-1);
    read(in_file,s);
    if paramcount=0 then gotoxy(x,y-1);
  end;

procedure signon;
  begin
    writeln(out_file,'RIGID FRAME ANALYSIS  By Gary Argraves - 8701.15');
    writeln(out_file);
    if paramcount>0 then
      begin
        writeln('Results of this run are in file: ',result_file);
        writeln;
        writeln('...Working    ...Please wait.');
      end;
    write(out_file,'Project Designation ? '); reads(proj_descrip);
    writeln(out_file,proj_descrip);
    if paramcount>0 then
      begin
        writeln(out_file,'Result file: ',result_file,'  Date: ',time);
        writeln(out_file);
      end;
  end;

 procedure fatal_err(s:str60);
   begin
     writeln(out_file,^G,s);
     close(in_file);
     close(out_file);
     halt;
   end;

procedure initial_var;
  var i,j:integer;
  begin
    lp:=1;
    mc:=0;
    lc:=0;
    for i:=1 to 50 do for j:=1 to 50 do d[i,j]:=0;
    for i:=1 to 40 do for j:=1 to 6 do b[i,j]:=0;
    for i:=1 to 40 do for j:=1 to 2 do c[i,j]:=0;
    for i:=1 to 10 do for j:=1 to 4 do e[i,j]:=0;
    for i:=1 to 50 do for j:=1 to 4 do ld[i,j]:=0;
    for i:=1 to 20 do for j:=1 to 2 do mp[i,j]:=0;
    for i:=1 to 50 do lm[i]:=0;
  end;

procedure get_size;
  begin

     writeln(out_file,'___INPUT FRAME SIZE___');
     write(out_file,'How many Stories ? ':25); readi(s);  writeln(out_file,s:6);
     write(out_file,'How many Portals ? ':25); readi(p);  writeln(out_file,p:6);
     writeln(out_file);
     q:=(s+1)*(p+1)+s;
  end;

procedure get_geometry;
  var
   i,j,z2  : integer;
   l       : real;
  begin
    writeln(out_file,'___INPUT FRAME GEOMETRY (units=feet)___');
    writeln(out_file,'   Input story heights.');
    l3:=0;
    for i:=1 to s do
      begin
        write(out_file,'Story #':12,i:2,' ? '); readr(l); writeln(out_file,l:8:2);
        l3:=l3+l;
        z2:=(p+1)*(s-i)+1;
        for j:=z2 to z2+p do c[j,2]:=l;
      end;
    writeln(out_file,'  Total height=',l3:7:2);
    writeln(out_file);
    writeln(out_file,'   Input beam lengths.');
    l2:=0;
    for i:=1 to p do
      begin
        write(out_file,'Beam #':12,i:2,' ? '); readr(l); writeln(out_file,l:8:2);
        l2:=l2+l;
        j:=i; repeat b[j,2]:=l; j:=j+p; until j>i+(s-1)*p;
      end;
    writeln(out_file,'  Total length=',l2:7:2);
    psl:=10;
    psm:=100;
    while l2/psl>9{inches} do
      begin psl:=psl+10; psm:=psm+100; end; (*** modify scale as required ***)
    writeln(out_file);
  end;

procedure get_member_sizes;  { by story }
  var
    i,j:integer;
  begin
    writeln(out_file,'___INPUT MEMBER SIZES - Moment of Inertia (I) (units=inch^4)___');
    for i:=1 to s do
      begin
        writeln(out_file);
        writeln(out_file,'STORY #':12,i:2);
        writeln(out_file,'-------':12,'--');
        for j:=1 to p do
          begin
            write(out_file,'Beam #':12,j:2,' ? ');
            readr(b[j+p*(s-i),1]); writeln(out_file,b[j+p*(s-i),1]:10:2);
          end;
        for j:=1 to p+1 do
          begin
            write(out_file,'Column #':12,j:2,' ? ');
            readr(c[j+(s-i)*(p+1),1]); writeln(out_file,c[j+(s-i)*(p+1),1]:10:2);
          end;
      end;
    writeln(out_file);
  end;

procedure get_base_joint_type;
  begin
    writeln(out_file,'___BASE JOINT TYPE___');
    writeln(out_file,'1. Fixed ':15);
    writeln(out_file,'2. Pinned':15);
    repeat
      write(out_file,'Enter 1 or 2 ? ':12); readi(z); writeln(out_file,z:3);;
    until z in [1,2];
    writeln(out_file);
  end;

procedure get_macro_loads;
  var
    i     : integer;
    w,j,k : real;    { load parameters }
    l4    : real;
  begin
    writeln(out_file,'  ___________ MACRO LOAD DEFINITION ____________');
    writeln(out_file);
    writeln(out_file,'            +-A--+---C---+');
    writeln(out_file,'         W  !    #########--W');
    writeln(out_file,'          \ ====================== -- Story 2');
    writeln(out_file,'        -+-#I      I      I      I');
    writeln(out_file,'         | #I      I      I      I');
    writeln(out_file,'      C--| #I      I      I      I');
    writeln(out_file,'         | #====================== -- Story 1');
    writeln(out_file,'        -+-#I      I      I      I');
    writeln(out_file,'      A--|  I  P1  I  P2  I  P3  I');
    writeln(out_file,'       __|__I______I______I______I___ Ground');
    writeln(out_file);
    writeln(out_file,'  S=0 :Lateral Load   S=-1 :Next Macro Definition');
    writeln(out_file,'  S>0 :Gravity Load   S=-2 :Quit Macro Definition');
    writeln(out_file);
(*    if lp<1 then lp:=1; *)
    repeat
      mc:=mc+1;
      mp[mc,1]:=lp;
      repeat
        writeln(out_file,'Macro #':12,mc:2,'  Load #',lp:2);
        writeln(out_file,'-------':12,'--','  ------','--');
        write(out_file,'Enter Story # S=? ':22); readi(i); writeln(out_file,i:5);
        if (i>s) then fatal_err('*** Error: Incorrect stories number ***');
        if i=0 then l4:=l3 else l4:=l2;
        if i>=0 then
          begin
             write(out_file,'Intensity W=? ':22); readr(w);  writeln(out_file,w:8:2);
             write(out_file,'Location  A=? ':22); readr(j);  writeln(out_file,j:8:2);
             write(out_file,'Length    C=? ':22); readr(k);  writeln(out_file,k:8:2);
             if j+k>l4 then fatal_err('*** Error: A+C> Total Length of Frame ***');
            ld[lp,1]:=i;
            ld[lp,2]:=w;
            ld[lp,3]:=j;
            ld[lp,4]:=k;
            mp[mc,2]:=lp;
            lp:=lp+1; { bump load pointer }
          end
        else if mp[mc,2]=0 then mc:=mc-1;
        writeln(out_file);
      until i<0;
    until i<-1;
  end;

procedure apply_macro_loads;
  var
    u,i,m : integer;
    w,j,k,t,t2 : real;

  procedure process_lateral_load; { data, FEM's & Shears }
    var
      h,o   : integer;
      a,b2,c2,l : real;

    begin
     for H:=1 TO S do
       begin
         O:=1+(S-H)*(P+1);                    { o is the column number }
         T:=T+C[O,2];                         { T is the height to top of column }
         T2:=T-C[O,2];                        { T2 is the height to column base }
         IF (J<T) and (J+K>T2) THEN           { Is the load in this section ? }
           begin
             A:=0;                            { Yes, A is the relative start position }
             IF J>T2 THEN A:=J-T2;            { Load starts in from segment beginning }
             B2:=C[O,2]-A;                    { B is the load length in this segment }
             IF J+K<T THEN B2:=J+K-A-T2;      { Load does not go to end of segment }
             L:=C[O,2];                       { L is the column length }
             C2:=A+B2;
             E[H,1]:=E[H,1]-(W*ABS(1/12/sqr(L)*(6*sqr(L)*(sqr(C2)-sqr(A))-8*L*(pw(C2,3)-pw(A,3))+3*(pw(C2,4)-pw(A,4)))));
             E[H,2]:=E[H,2]+(W*ABS(1/12/sqr(L)*(4*L*(pw(C2,3)-pw(A,3))-3*(pw(C2,4)-pw(A,4)))));
             E[H,3]:=E[H,3]+(W*B2*(A+B2/2));
             E[H,4]:=E[H,4]+(W*B2*(L-A-B2/2));
           end;
       end;
    end;

  procedure process_beam_load;  { FEM'S & shears }
    var
      h,o   : integer;
      a,b2,c2,l : real;
    begin
      FOR H:=1 TO P do
        begin
          O:=H+(S-I)*P;
          T:=T+B[O,2];
          T2:=T-B[O,2];
          IF (J<T) and (J+K>T2) THEN
            begin
              A:=0; IF J>T2 THEN A:=J-T2;
              B2:=B[O,2]-A; IF J+K<T THEN B2:=J+K-A-T2;
              L:=B[O,2];
              C2:=A+B2;
              B[O,3]:=B[O,3]-(W*ABS(1/12/sqr(L)*(6*sqr(L)*(sqr(C2)-sqr(A))-8*L*(pw(C2,3)-pw(A,3))+3*(pw(C2,4)-pw(A,4)))));
              B[O,4]:=B[O,4]+(W*ABS(1/12/sqr(L)*(4*L*(pw(C2,3)-pw(A,3))-3*(pw(C2,4)-pw(A,4)))));
              B[O,5]:=B[O,5]+(W*B2*(A+B2/2));
              B[O,6]:=B[O,6]+(W*B2*(L-A-B2/2));
            end;
        end;
    end;

  begin
     writeln(out_file,'___APPLY LOADS TO FRAME___');
     writeln(out_file);
     writeln(out_file,'  # > 0 :Macro Loading');
     writeln(out_file,'  # < 0 :Single Load');
     writeln(out_file,'  # = 0 :Quit loading');
     writeln(out_file);
     lc:=0;
     write(out_file,'Enter Loading Designation ? ');
     if paramcount>0 then readln(in_file);
     reads(load_descrip); load_descrip:=parse_up(load_descrip);
     writeln(out_file,load_descrip);
     repeat
       write(out_file,'Enter # ? ');  readi(m); writeln(out_file,m:3);
       if (m>0) and (m<=mc) then begin
           for u:=mp[m,1] to mp[m,2] do
             begin
               lc:=lc+1; lm[lc]:=u;  { save all loads applied for moment dia }
               i:=trunc(ld[u,1]); w:=ld[u,2]; j:=ld[u,3]; k:=ld[u,4]; t:=0; t2:=0;
               if i=0 then process_lateral_load else process_beam_load;
             end;
         end
       else if (m<0) and (m>-lp) then
              begin
                u:=abs(m);
                lc:=lc+1; lm[lc]:=u;
                i:=trunc(ld[u,1]); w:=ld[u,2]; j:=ld[u,3]; k:=ld[u,4]; t:=0; t2:=0;
                if i=0 then process_lateral_load else process_beam_load;
              end
            else if m<>0 then fatal_err('*** Error: No such Macro or Load. ***');
     until m=0;
   writeln(out_file);
  end;

procedure setup_matrix;
  var
    z2,u,i,j : integer;
    o : real;
  label
    NEXT_J;
  begin
    Z2:=S+1;
    IF Z=1 THEN Z2:=S;
    for U:=1 TO Z2 do
      for J:=S+U*(1+P)-P TO S+U*(1+P) do
        begin
          IF U<>1 then begin
            O:=C[J-S-P-1,1]/C[J-S-P-1,2];
            IF J=S+U*(1+P)-P THEN D[J,Q+1]:=E[S+2-U,1];
            D[J,J]:=4*O+D[J,J];
            D[J,J-P-1]:=2*O;
            O:=O/C[J-S-P-1,2];
            D[J,U-1]:=6*O+D[J,U-1];
            IF U=S+1 THEN GOTO NEXT_J;
            D[U,J]:=-6*O;
            D[U,J-P-1]:=-6*O;
            D[U,U]:=D[U,U]+12*O/C[J-S-P-1,2];
            D[U,U-1]:=D[U,U-1]-12*O/C[J-S-P-1,2];
            D[J,U]:=D[J,U]-6*O;
          end;
          O:=C[J-S,1]/C[J-S,2];
          IF J=S+U*(1+P)-P THEN D[J,Q+1]:=D[J,Q+1]+E[S+1-U,2];
          D[J,J]:=D[J,J]+4*O;
          IF (Z<>1) or (u<>s) THEN begin
            D[U,J+P+1]:=6*O/C[J-S,2];
            D[J,J+P+1]:=2*O;
          end;
          O:=O/C[J-S,2];
          D[U,J]:=D[U,J]+6*O;
          D[J,U]:=D[J,U]+6*O;
          D[U,U]:=D[U,U]+12*O/C[J-S,2];
          IF U<>S THEN begin
            D[J,U+1]:=D[J,U+1]-6*O;
            D[U,U+1]:=D[U,U+1]-12*O/C[J-S,2];
          end;
          IF J<>S+U*(1+P) THEN begin
            O:=B[J-S-U+1,1]/B[J-S-U+1,2];
            D[J,J+1]:=2*O;
            D[J,J]:=D[J,J]+4*O;
            D[J,Q+1]:=D[J,Q+1]+B[J-S-U+1,3];
          end;
          IF J<>S+U*(1+P)-P THEN begin
            O:=B[J-S-U,1]/B[J-S-U,2];
            D[J,J]:=D[J,J]+4*O;
            D[J,J-1]:=D[J,J-1]+2*O;
            D[J,Q+1]:=D[J,Q+1]+B[J-S-U,4];
          end;
         NEXT_J:
        end{for j};

     FOR I:=1 TO S do begin
       U:=S+1-I;
       D[I,Q+1]:=(E[U,3]+E[U,2]+E[U,1])/C[(I-1)*(P+1)+1,2];
       IF I>1 THEN D[I,Q+1]:=D[I,Q+1]+(E[U+1,4]-E[U+1,2]-E[U+1,1])/C[(I-2)*(P+1)+1,2];
     end;
  end;

procedure solve_matrix;
  var
    i,j,l : integer;
    m : real;
    matrix_err,go,restart : boolean;

begin
  matrix_err:=false;
  repeat
    restart:=false;
    i:=1;
    repeat
      IF (Z=1) and (I=Q-P) THEN exit;
      IF D[I,I]=0 THEN
        begin
          IF I=Q THEN matrix_err:=true;
          if not matrix_err then
            begin
              j:=q; go:=false; restart:=false;
              repeat
                IF (J<>I) and (D[J,I]<>0) then
                  begin
                    FOR L:=1 TO Q+1 do D[I,L]:=D[I,L]+D[J,L];
                    if j>i then go:=true
                    else begin restart:=true; i:=q; end;
                  end;
                j:=j-1;
                if (not go) and (not restart) and (j=0) then matrix_err:=true;
              until (j=0) or (go) or (restart);
            end;
          if matrix_err then fatal_err('*** Matrix solution Error ***');
        end;
      if not restart then
        begin
          M:=D[I,I];
          FOR L:=1 TO Q+1 do D[I,L]:=D[I,L]/M;
            FOR J:=1 TO Q do
              IF J<>I THEN
                begin M:=D[J,I];
                      FOR L:=1 TO Q+1 do D[J,L]:=D[J,L]-M*D[I,L];
                end;
        end;
      i:=i+1;
    until i>q;
  until not restart;
  end;

procedure plot_beam_mom(story,beam:integer; m1,v1:real);
  var
    w,a,c2,l: real;     { load, start load, length load, length beam }
    sb,eb   : real;     { start beam, end beam length }
    sh      : real;     { story height }
    last_lc : integer;  { last load count }
    i,k     : integer;  { general }
    ox,oy   : integer;
    x       : real;     { distance along member }

  procedure find_load(var w,a,c:real);
    var
      sl,el : real; { start load, end load }
    begin
      w:=0;
      while (last_lc<=lc) and (w=0) do
        begin
          if ld[lm[last_lc],1]=story then
            begin
              sl:=ld[lm[last_lc],3];
              el:=sl+ld[lm[last_lc],4];
              if (sl<sb) and (el>sb) then
                begin w:=ld[lm[last_lc],2]; a:=0; c:=el-sb; end
              else if (sl>=sb) and (sl<eb) then
                     begin w:=ld[lm[last_lc],2]; a:=sl-sb; c:=el-sl; end;
              if (w<>0) and (a+c>l) then c:=l-a;
            end;
          last_lc:=last_lc+1;
        end;
    end;

  procedure tab_moments;
    var
      x,mx : real;
    begin
      x:=0.5; k:=1;
      repeat
        mx:=0;
        if (x>a) and (x<=a+c2) then  mx:=w/2*sqr(x-a)
        else if x>a+c2 then mx:=c2*w*(x-(a+c2/2));
        r[k]:=r[k]+mx;
        x:=x+0.5; k:=k+1;
      until x>=l;
      r[k]:=r[k]+c2*w*(l-(a+c2/2));
    end;

  procedure tab_shears;
    var
      x,vx : real;
    begin
      x:=0.5; k:=1;
      repeat
        vx:=0;
        if (x>a) and (x<=a+c2) then  vx:=(x-a)*w
        else if x>a+c2 then vx:=c2*w;
        r[k]:=r[k]-vx;
        x:=x+0.5; k:=k+1;
      until x>=l;
      r[k]:=r[k]-c2*w;
    end;


  begin
    sb:=0; for i:=1 to beam-1 do sb:=sb+b[p*(s-story)+i,2];
    eb:=sb+b[p*(s-story)+beam,2];
    sh:=0; for i:=1 to story do sh:=sh+c[(p+1)*(s-i)+1,2];
    l:=eb-sb;
    last_lc:=1;
    for i:=0 to 200 do r[i]:=0;
    i:=0; x:=0;
    repeat
      if shear then r[i]:=v1 else r[i]:=m1-v1*x;
      i:=i+1; x:=x+0.5;
    until x>=l; if shear then r[i]:=v1 else r[i]:=m1-v1*l;
    repeat
      find_load(w,a,c2);
      if w<>0 then begin if shear then tab_shears else tab_moments; end;
    until last_lc>lc;
    k:=i;
    ox:=round((sb/psl+0.10*(beam)+0.50)*mf_plt);
    oy:=round((sh/psl+0.10*(story))*mf_plt);
(*    writeln(plot_file,'SP 1;');  *)
    writeln(plot_file,'PU ',ox:5,',',oy:5,';');
    writeln(plot_file,'PD ',round(l/psl*mf_plt)+ox:5,',',oy:5,';');
    writeln(plot_file,'PU ',ox:5,',',oy:5,';');
(*    writeln(plot_file,'SP 1;'); *)
    x:=0;
    for i:=0 to k-1 do
      begin
        writeln(plot_file,'PD ',round(x/psl*mf_plt)+ox:5,',',round(r[i]/psm*mf_plt)+oy:5,';');
        x:=x+0.5;
      end;
    writeln(plot_file,'PD ',round(l/psl*mf_plt)+ox:5,',',round(r[k]/psm*mf_plt)+oy:5,';');
    writeln(plot_file,'PD ',round(l/psl*mf_plt)+ox:5,',',oy:5,';');
  end{plot_mom_dia};

procedure output_beam;
  var
    i,j,h : integer;
    o,m,m1 : real;
  begin
    if not shear then begin
      writeln(out_file,'Beam':11, 'Moment-L':10, 'Moment-R':10,'Shears-L':10,'Shears-R':10);
      writeln(out_file, '(#)':11,'(Kip-ft)':10,'(Kip-ft)':10,'(Kips)':10,'(Kips)':10);
      writeln(out_file,'----':11,'--------':10,'--------':10,'------':10,'------':10);
    end;
    for i:=0 to s-1 do
      begin
        if not shear then writeln(out_file,'Story',s-i:2);
        for j:=1 to p do
          begin
             H:=I*P+J;
             O:=B[H,1]/B[H,2];
             M:=O*(4*D[S+I+H,Q+1]+2*D[S+I+H+1,Q+1])-B[H,3];
             M1:=O*(2*D[S+I+H,Q+1]+4*D[S+I+H+1,Q+1])-B[H,4];
             D[H,1]:=(M+M1+B[H,6])/B[H,2];
             D[H,2]:=(B[H,5]-M-M1)/B[H,2];
             if not shear then writeln(out_file,j:11,m:10:3,m1:10:3,d[h,1]:10:3,d[h,2]:10:3);
             plot_beam_mom(s-i,j,m,d[h,1]);
          end;
      end;
    if not shear then writeln(out_file);
  end;

procedure plot_column_mom(story,column:integer; m1,v1:real);
  var
    w,a,c2,l: real;     { load, start load, length load, length beam }
    sc,ec   : real;     { start colum, end column length }
    bw      : real;     { beam width }
    last_lc : integer;  { last load count }
    i,k     : integer;  { general }
    ox,oy   : integer;
    x       : real;

  procedure find_load(var w,a,c:real);
    var
      sl,el : real; { start load, end load }
    begin
      w:=0;
      while (last_lc<=lc) and (w=0) do
        begin
          if (ld[lm[last_lc],1]=0) and (column=1) then
            begin
              sl:=ld[lm[last_lc],3];
              el:=sl+ld[lm[last_lc],4];
              if (sl<sc) and (el>sc) then
                begin w:=ld[lm[last_lc],2]; a:=0; c:=el-sc; end
              else if (sl>=sc) and (sl<ec) then
                     begin w:=ld[lm[last_lc],2]; a:=sl-sc; c:=el-sl; end;
              if (w<>0) and (a+c>l) then c:=l-a;
            end;
          last_lc:=last_lc+1;
        end;
    end;

  procedure tab_moments;
    var
      x,mx : real;
    begin
      x:=0.5; k:=1;
      repeat
        mx:=0;
        if (x>a) and (x<=a+c2) then  mx:=w/2*sqr(x-a)
        else if x>a+c2 then mx:=c2*w*(x-(a+c2/2));
        r[k]:=r[k]-mx;
        x:=x+0.5; k:=k+1;
      until x>=l;
      r[k]:=r[k]-c2*w*(l-(a+c2/2));
    end;

  procedure tab_shears;
    var
      x,vx : real;
    begin
      x:=0.5; k:=1;
      repeat
        vx:=0;
        if (x>a) and (x<=a+c2) then vx:=w*(x-a)
        else if x>a+c2 then vx:=w*c2;
        r[k]:=r[k]+vx;
        x:=x+0.5; k:=k+1;
      until x>=l;
      r[k]:=r[k]+c2*w;
    end;


  begin
    sc:=0; for i:=1 to story-1 do sc:=sc+c[(p+1)*(s-i)+1,2];
    ec:=sc+c[(p+1)*(s-story)+1,2];
    bw:=0; for i:=1 to column-1 do bw:=bw+b[p*(s-story)+i,2];
    l:=ec-sc;
    last_lc:=1;
    for i:=0 to 200 do begin r[i]:=0; end;
    i:=0; x:=0;
    repeat
      if shear then r[i]:=v1 else r[i]:=-(m1+v1*x);
      i:=i+1; x:=x+0.5;
    until x>=l; if shear then r[i]:=v1 else r[i]:=-(m1+v1*l);
    repeat
      find_load(w,a,c2);
      if w<>0 then begin if shear then tab_shears else tab_moments; end;
    until last_lc>lc;
    k:=i;
    ox:=round((bw/psl+0.05*(2*column-1)+0.50)*mf_plt);
    oy:=round((sc/psl+0.05*(2*story-1))*mf_plt);
(*    writeln(plot_file,'SP 1;'); *)
    writeln(plot_file,'PU ',ox:5,',',oy:5,';');
    writeln(plot_file,'PD ',ox:5,',',round(l/psl*mf_plt)+oy:5,';');
    writeln(plot_file,'PU ',ox:5,',',oy:5,';');
(*    writeln(plot_file,'SP 1;'); *)
    x:=0;
    for i:=0 to k-1 do
      begin
        writeln(plot_file,'PD ',round(r[i]/psm*mf_plt)+ox:5,',',round(x/psl*mf_plt)+oy:5,';');
        x:=x+0.5;
      end;
    writeln(plot_file,'PD ',round(r[k]/psm*mf_plt)+ox:5,',',round(l/psl*mf_plt)+oy:5,';');
    writeln(plot_file,'PD ',ox:5,',',round(l/psl*mf_plt)+oy:5,';');
    if (story=1) and (column=p+1) then
      begin
        writeln(plot_file,'PU ',mf_plt*1.97:5:0,',',7.28*mf_plt:5:0,';');
        writeln(plot_file,'SI 0.28,0.35;'); { set character height }
        if shear then writeln(plot_file,'LBRIGID FRAME ANALYSIS SHEAR DIAGRAM',term)
        else writeln(plot_file,'LBRIGID FRAME ANALYSIS MOMENT DIAGRAM',term);
        writeln(plot_file,'SI 0.20,0.27;');
        writeln(plot_file,'PU ',1.18*mf_plt:5:0,',',7.09*mf_plt:5:0,';');
        writeln(plot_file,'LBProject Name: ',proj_descrip,term);
        writeln(plot_file,'PU ',1.18*mf_plt:5:0,',',6.93*mf_plt:5:0,';');
        writeln(plot_file,'LBResult File: ',result_file,'  Date: ',time,term);
        writeln(plot_file,'PU ',1.18*mf_plt:5:0,',',6.79*mf_plt:5:0,';');
        writeln(plot_file,'LBLoad Designation: ',load_descrip,term);
        writeln(plot_file,'PU ',1.18*mf_plt:5:0,',',6.63*mf_plt:5:0,';');
        if shear then
          writeln(plot_file,'LBScales: X-',psl:3,' Feet/Inch  Y-',psm:4,' Kips/Inch',term)
        else writeln(plot_file,'LBScales: X-',psl:3,' Feet/Inch  Y-',psm:4,' Kip-Feet/Inch',term);
        writeln(plot_file,'PU 0,0;');
      end;
  end{plot_mom_dia};

procedure output_column;
  var
    i,j,h : integer;
    o,l,m,m1 : real;
  begin
    if not shear then begin
      writeln(out_file,'Column':11,'Moment':10,'Shears':10,'Axial':10);
      writeln(out_file, '(#)':11,'(Kip-ft)':10,'(Kips)':10,'(Kips)':10);
      writeln(out_file,'----':11,'--------':10,'------':10,'------':10);
    end;
    for i:=0 to s-1 do
      begin
        if not shear then writeln(out_file,'Story',s-i:2);
        for j:=1 to p+1 do
          begin
            H:=I*(P+1)+J;
            O:=C[H,1]/C[H,2];
            L:=6/C[H,2]*D[I+1,Q+1];
            M:=O*(4*D[S+H,Q+1]+2*D[S+H+P+1,Q+1]+L);
            M1:=O*(2*D[S+H,Q+1]+4*D[S+H+P+1,Q+1]+L);
            IF I<>S-1 THEN begin
              L:=O/C[H,2]*6*D[I+2,Q+1];
              M1:=M1-L;
              M:=M-L;
            end;
            IF J<=1 THEN begin
              M:=M-E[S-I,2];
              M1:=M1-E[S-I,1];
            end;
            O:=0;
            IF J<P+1 THEN O:=D[I*P+J,1];
            IF J>1 THEN O:=O+D[I*P+J-1,2];
            IF I>0 THEN O:=O+D[J,3];
            D[J,3]:=O;
            L:=0; IF J=1 THEN L:=E[S-I,3];
            if not shear then writeln(out_file,j:11,m:10:3,(m+m1-l)/c[h,2]:10:3,o:10:3);
            L:=0; IF J=1 THEN L:=E[S-I,4];
            if not shear then writeln(out_file,m1:21:3,(m+m1+l)/-c[h,2]:10:3);
            plot_column_mom(s-i,j,m1,((m+m1+l)/-c[h,2]));
          end;
      end;
    if not shear then writeln(out_file);
  end;

procedure output_rotation;
  var
    i,j,l : integer;

  begin
    writeln(out_file,'Joint':11,'Rotation':15,'Displacement':15);
    writeln(out_file, '(#)':11,'(Radians)':15,'(Feet)':15);
    writeln(out_file,'----':11,'---------':15,'------':15);
    L:=0;
    for i:=s+1 to q do
      begin
        J:=I-S;
        IF (J-1)/(P+1)=int((J-1)/(P+1)) THEN
          begin
            writeln(out_file,'Story',trunc(S-(J-1)/(P+1)):2);
            L:=L+1;
            IF L<=S THEN
              writeln(out_file,j:11,d[i,q+1]*0.0048:15:9,d[l,q+1]*0.0048:15:8)
            else writeln(out_file,j:11,d[i,q+1]*0.0048:15:9);
          end
        else writeln(out_file,j:11,d[i,q+1]*0.0048:15:9);
      end;
    writeln(out_file);
  end;

procedure reset_matrix;
  var
    i,j : integer;
  begin
    for i:=1 to q+1 do for j:=1 to q+1 do d[i,j]:=0;
    for i:=1 to 40  do
      for j:=3 to 6 do
        begin
          b[i,j]:=0;
          if i<11 then e[i,j-2]:=0;
        end;
  end;

procedure set_IO_devices;
  begin
    result_file:='****';
    if paramcount>0 then
      begin
        if not exist(paramstr(1)) then
          begin writeln('Can not find file: ',paramstr(1)); halt; end
        else if paramcount=2 then result_file:=paramstr(2)
             else if pos('.',paramstr(1))=0 then result_file:=paramstr(1)+'.OUT'
                  else result_file:=copy(paramstr(1),1,pos('.',paramstr(1))-1)+'.OUT';
      end;
    result_file:=parse_up(result_file);
    case paramcount of
      0:begin assign(in_file,'con'); assign(out_file,'con'); end;
      1:begin assign(in_file,paramstr(1)); assign(out_file,result_file); end;
      2:begin assign(in_file,paramstr(1)); assign(out_file,paramstr(2)); end;
    end{case};
    reset(in_file);
    rewrite(out_file);
  end;

procedure user_option;
  begin
    writeln(out_file,'============================================');
    writeln(out_file,'1. Define More Macros.':30);
    writeln(out_file,'2. Reload the Frame.  ':30);
    writeln(out_file,'3. Exit Program.      ':30);
    writeln(out_file);
    repeat
      option:=0;
      write(out_file,'Enter Option ? '); readi(option); writeln(out_file,option:3);
    until option in [1,2,3];
    writeln(out_file);
  end;

  procedure change_plot_page;
    begin
      writeln;
      quest(0,0,'Put New Paper in Plotter.  ...Press <SPACE BAR>',[' ']);
      writeln;
      writeln(plot_file,'IN; SP 1;');
    end;

begin (****  M A I N   L I N E  ****)
  assign(plot_file,'AUX');
  rewrite(plot_file);
  time:=when;
  set_IO_devices;
  signon;

  initial_var;
  get_size;
  get_geometry;
  get_member_sizes;
  get_base_joint_type;
  option:=1;
  repeat
    if option=1 then get_macro_loads;
    apply_macro_loads;
    setup_matrix;
    solve_matrix;
    change_plot_page; shear:=false; output_beam; output_column;
    change_plot_page; shear:=true;  output_beam; output_column;
    output_rotation;
    reset_matrix;
    user_option;
  until option=3;
  writeln(out_file,'*****  E N D   O F   R U N  *****');
  close(in_file);
  close(out_file);
  close(plot_file);
end.