XProfan

WindowTitle "    LINEARE OPTIMIERUNG MITTELS SIMPLEX-ALGORITHMUS  "
'Engl. Original source in Fortran-77 by John H Matthews 1986-04
'translated to XProfan-11.2a free/German by p.specht 2021-05
'Experimental DEMO only, no warranties whatsoever! Keine Gewähr!
WindowStyle 24:font 2:randomize:set("decimals",17)
var co&=rgb(200+rnd(56),200+rnd(56),200+rnd(56))
Window 0,0-%maxx,%maxy-40:CLS co&
declare i&,j&,k&,MMAX&,NMAX&,M&,M1&,M2&,M3&,F$,F1$,mm&,kp&,bmax!
declare icase&,n&,EPS!,R!,nl1&,nl2&,ir&,q1!,iafb&,kk&,ii&,iabf&
declare test!,ip&,q!,qp!,q0!,i1&,k1&,PIV!,kh&,iis&,M12&,ente&
MMAX& = 25 : NMAX& = 25
declare A![MMAX&,NMAX&],iposv&[MMAX&],izrov&[NMAX&]
PRINT
PRINT " Anzahl Variablen der Zielfunktion: "; : INPUT N&
PRINT " (Nichtnegativitätsbedingungen werden automatisch erzeugt!)\n"
PRINT " Anzahl der <= (Kapazitätsgrenzen): "; : INPUT M1&
PRINT " Anzahl der >= (Mindestauslastung): "; : INPUT M2&
PRINT " Anzahl der ==      (Fixsetzungen): "; : INPUT M3&
PRINT
M& = M1& + M2& + M3&'Total number of constraints
F$  = " ########0.0###############;-########0.0###############;        0.0"'"%g"
F1$ = " ########0.0###############;-########0.0###############;        0.0"'"%g"

whileloop M&+2:i&=&Loop

    whileloop N&+1:j&=&Loop

        A![i&,j&]=0.0

    endwhile

endwhile

PRINT "\n Eingabe der Zielfunktion "
print " (Maximiert mit positiven Zielkoeffizienten, minimiert mit Negativen!)\n"

whileloop 2,n&+1:i&=&Loop

    PRINT " Koeffizient "; int(i&-1); ": ";:INPUT A![1,i&]

endwhile

PRINT " Vorkonstante [meist 0]: ";:INPUT A![1,1]
PRINT:print "\n Beachte die Reihenfolge der Ungleichungen: <=, >=, ==\n"

whileloop m&:i&=&Loop

    PRINT " Eingabe der Ungleichung Nr."; i&; ":"

    whileloop 2,n&+1:j&=&Loop

        PRINT " Koeffizient "; int(j& - 1); ": ";:INPUT R!
        A![i&+1,j&] = -1*R!

    endwhile

    PRINT " Zugehöriger Grenzwert: ";:INPUT A![i&+1,1]
    print

endwhile

PRINT
CLS co&
'  PRINT " Kontroll-Tableau: ":print
'  whileloop m&+1:i&=&Loop
'    whileloop n&+1:j&=&Loop
'      print tab((j&-1)*27);A![i&,j&];
'    endwhile
'    PRINT
'  endwhile
print
GOSUB "S1000"'call simplx(A,M,N,M1,M2,M3,ICASE,IZROV,IPOSV)

IF icase& = 0'result ok.

    PRINT
    PRINT " Maximum der Zielfunktion = ";
    PRINT format$(F1$, A![1,1])
    ente&=0

    whileloop n&:i&=&Loop

        whileloop m&:j&=&Loop

            IF iposv&[j&] = i&

                PRINT " Optimales X";i&;" = ";format$(F1$,A![j&+1,1])
                ente&=1:break

            ENDIF

        endwhile

        ifnot ente&

            PRINT "  X";i&;" = ";format$(F1$,0.0)
            ente&=0

        endif

    endwhile

ELSE

    PRINT " Keine Lösung gefunden (Fehlercode = "; icase&; ")."

ENDIF

waitinput
END
'  Simplex Algorithmus
S1000:
' simplx(A[1..,1.. ], m, n, m1, m2, m3, icase, izrov, iposv)
Declare l1&[MMAX&],l2&[MMAX&],l3&[MMAX&]
EPS!=0.000001'<<<<<<<< Adaptieren an verwendete Größenordnungen!

IF M& <> (M1& + M2& + M3&)

    PRINT " Bad input constraint counts in simplx."
    RETURN

ENDIF

nl1& = N&

whileloop n&:k&=&Loop

    l1&[k&] = k&
    izrov&[k&] = k&

endwhile

nl2& = M&

whileloop m&:i&=&Loop

    IF A![i&+1,1] < 0.0

        PRINT " Bad input tableau in simplx, Constants b_i must be nonnegative."
        RETURN

    ENDIF

    l2&[i&] = i&
    iposv&[i&] = N& + i&

endwhile

whileloop m2&:i&=&Loop

    l3&[i&] = 1

endwhile

ir& = 0
case (M2& + M3&) = 0 : GOTO "G30"
ir& = 1

whileloop n&+1:k&=&Loop

    q1!=0.0

    whileloop M1&+1,M&:i&=&Loop

        q1!=q1!+A![i&+1,k&]

    endwhile

    A![M&+2,k&]= -1*q1!

endwhile

G10:
iabf&=0 : mm&=M&+1
GOSUB "S2000"'simp1(a,m+1,l1,nl1,0,kp,bmax)

IF (bmax! <= EPS!) AND (A![M&+2,1] < (-1*EPS!))

    icase& = -1  :print " Auxiliary objective function is still negative and can’t be improved! No feasible solution exists."
    RETURN

ELSEIF (bmax! <= EPS!) AND (A![M&+2,1] <= EPS!)

    m12& = M1& + M2& + 1

    IF m12& <= M&'THEN

        whileloop m12&,m&:ip&=&Loop

            IF iposv&[ip&] = ip& + N&

                iabf& = 1: mm& = ip&
                GOSUB "S2000"'simp1(a,ip,l1,nl1,1,kp,bmax)
                case bmax! > EPS!:GOTO "G1"

            ENDIF

        endwhile

    ENDIF

    ir& = 0
    m12& = m12& - 1
    case (M1&+1) > m12& :GOTO "G30"

    whileloop m1&+1,M1&+M2&

        IF l3&[i& - M1&]=1

            whileloop n&+1:k&=&Loop

                A![i&+1,k&] = -1.0*A![i&+1,k&]

            endwhile

        ENDIF

    endwhile

    GOTO "G30"'Go phase two.

ENDIF

GOSUB "S2100"'call simp2(a,m,n,l2,nl2,ip,kp,q1)

IF ip& = 0

    icase& = -1 :print " Ohne Grenzen: Es existiert keine ermittelbare konkrete Lösung!"
    RETURN

ENDIF

G1:
i1&=M&+1:k1&=N&
GOSUB "S2200"'call simp3(a,m+1,n,ip,kp)

IF iposv&[ip&] >= (N&+M1&+M2&+1)

    whileloop nl1&:k&=&Loop

        case l1&[k&] = kp& : GOTO "G2"

    endwhile

    G2:
    nl1& = nl1& - 1

    whileloop k&, nl1& : iis&=&Loop

        l1&[iis&] = l1&[iis&+1]

    endwhile

ELSE

    case iposv&[ip&]<(N&+M1&+1):GOTO "G20"
    kh& = iposv&[ip&]-M1&-N&
    case l3&[kh&] = 0 : GOTO "G20"
    l3&[kh&] = 0

ENDIF

A![M&+2,kp&+1]=A![M&+2,kp&+1]+1.0

whileloop m&+2:i&=&Loop

    A![i&,kp&+1] = -1.0*A![i&,kp&+1]

endwhile

G20:
iis& = izrov&[kp&]
izrov&[kp&] = iposv&[ip&]
iposv&[ip&] = iis&
case ir& <> 0 : GOTO "G10"
G30:
iabf& = 0 : mm& = 0
GOSUB "S2000"'call simp1(a,0,l1,nl1,0,kp,bmax)

IF bmax! <= EPS!

    icase& = 0
    RETURN

ENDIF

GOSUB "S2100"'call simp2(a,m,n,l2,nl2,ip,kp,q1)

IF ip& = 0 :Print " Die Zielfunktion ist unbeschränkt!"

    icase& = 1
    RETURN

ENDIF

i1& = M& : k1& = N&
GOSUB "S2200"'call simp3(a,m,n,ip,kp)
GOTO "G20"
RETURN
S2000:
' Subroutine simp1(a, mm, l1, nl1, iabf, kp, bmax)
kp& = l1&[1]
bmax! = A![mm&+1,kp&+1]
case nl1&<2:RETURN

whileloop 2,nl1&:k&=&Loop

    IF iabf& = 0'THEN

        test! = A![mm&+1,l1&[k&]+1]-bmax!

    ELSE

        test! = ABS(A![mm&+1,l1&[k&]+1]) - ABS(bmax!)

    ENDIF

    IF test! > 0'THEN

        bmax! = A![mm&+1,l1&[k&]+1]
        kp& = l1&[k&]

    ENDIF

endwhile

RETURN
S2100:
'Subroutine simp2(a, m, n, l2, nl2, ip, kp, q1)
EPS! = 0.000001
ip& = 0
case nl2& < 1 : RETURN

whileloop nl2&:i&=&Loop

    case A![i&+1,kp&+1] < (-1*EPS!) : GOTO "G2102"

endwhile

RETURN'No possible pivots. Return with message.
G2102:
q1! = -1*A![l2&[i&]+1,1]/A![l2&[i&]+1,kp&+1]
ip& = l2&[i&]
case (i&+1) > nl2&:RETURN

WhileLoop i&+1,nl2&:i&=&Loop

    ii& = l2&[i&]

    IF A![ii&+1,kp&+1] < (-1*EPS!)

        q! = -1*A![ii&+1,1] / A![ii&+1,kp&+1]

        IF q! < q1!

            ip& = ii&
            q1! = q!

        ELSEIF q! = q1! :print " Degenerierter Fall!"

            whileloop n& : k&=1

                qp! = -1*A![ip&+1,k&+1] / A![ip&+1,kp&+1]
                q0! = -1*A![ii&+1,k&+1] / A![ii&+1,kp&+1]
                case q0! <> qp! : BREAK'GOTO "G2106"

            endwhile

            G2106:
            case q0! < qp! : ip& = ii&

        ENDIF

    ENDIF

endwhile

RETURN
S2200:
' Subroutine simp3(a, i1,k1,ip,kp)
piv! = 1 / A![ip&+1,kp&+1]

IF i1&>=0

    whileloop i1&+1:ii&=&Loop

        IF (ii&-1) <> ip&

            A![ii&,kp&+1] = A![ii&,kp&+1]*piv!

            whileloop k1&+1:kk&=&Loop

                IF (kk&-1) <> kp&

                    A![ii&,kk&]=A![ii&,kk&]-A![ip&+1,kk&]*A![ii&,kp&+1]

                ENDIF

            endwhile

        ENDIF

    endwhile

ENDIF

whileloop k1&+1:kk&=&Loop

    case (kk&-1)<>kp&: A![ip&+1,kk&]= -1*A![ip&+1,kk&]*piv!

endwhile

A![ip&+1,kp&+1] = piv!
RETURN
' end of file simp.xprf