/* Rational simplex 21.08.2020 */ /* global bigrat */ /* Determine initial variables */ var EPSILON = bigRat(1, 1E10), FACTOR = 25, INFINITY = nines(), TAB = 7, ZERO = bigRat(), a = new Array(), b = new Array(), e = new Array(), n1 = new Array(), n2 = new Array(), s = new Array(), bv = new Array(), mbv = new Array(), nbv = new Array(), tab = new Array(), i, j, k, l, m, n, beq, bet, ble, cb, cn, ms, ns, prog = 1, rule = 1, op = 1, opt = 1, max, min, few = true, less, nor1 = false, nor2 = false, output = true, perts, red = false, scal = false, solve, steps, status, str, tmp; function form() { /* Build and write form */ var form = '
'; form += ""; form += ""; form += ''; form += ''; form += ''; form += '
Bitte lineares Programm eingeben (Trennzeichen Leerzeichen,
letzte Zeile Zielfunktion, erste Spalte rechte Seiten wie im Beispiel):
Ergebnis (es kommt nicht auf die Anzahl der korrekten Dezimalstellen an, sondern auf die Zeitkomplexität und die Anzahl der Schritte des Simplexverfahrens):
'; form += ''; form += ''; if (output) form += ''; else form += ''; if (prog == 1) form += ''; else form += ''; if (prog == 2) form += ''; else form += ''; if (rule == 1) form += ''; else form += ''; if (rule == 2) form += ''; else form += ''; if (rule == 3) form += ''; else form += ''; if (opt == 1) form += ''; else form += ''; if (opt == 2) form += ''; else form += ''; if (op == 1) form += ''; else form += ''; if (op == 2) form += '
'; else form += '
'; form += ''; if (red) form += ''; else form += ''; if (few) form += ''; else form += ''; if (nor1) form += ''; else form += ''; if (scal) form += ''; else form += ''; if (nor2) form += '
'; else form += ''; document.writeln(form); } function decimal(val) { /* Format single value to decimal */ return val.toDecimal(16).replace(/\./, ",").toUpperCase(); } function format(val) { /* Format single value */ var arg = (isFinite(val)) ? String(val).replace(/\./, ",").toUpperCase() : (val > 0) ? "∞" : "-∞"; return (arg.length > 2 && arg.indexOf("/1") == arg.length - 2) ? arg.substr(0, arg.length - 2) : arg; } function space(val) { /* Determine space */ return " ".substr(0, Math.min(val, 128)); } function finite() { /* Set small values to zero and check if there are infinite values */ for (i = 0; i < m; i++) for (j = 0; j < n; j++) if (!isFinite(a[i][j].numerator) || !isFinite(a[i][j].denominator)) return solve = false; return solve = true; } function nines() { /* Determine bigRat for infinity */ var big = "99999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999"; big += big; return bigRat(big + big); } function norm_scale() { /* Normalise and scale if demanded */ var min, obj, sca; if (nor1 || nor2) { if (output) str += "Die 1. Normierung wird durchgeführt.\n"; for (i = 1; i < m; i++) { sca = ZERO; for (j = 1; j < n; j++) sca = sca.add(a[i][j].times(a[i][j])); sca = bigRat(Math.sqrt(sca)); for (j = 0; j < n; j++) a[i][j] = a[i][j].divide(sca); n1[i] = sca.reciprocate(); } } if (scal) { if (output) str += "Die Skalierung wird durchgeführt.\n"; for (j = 1; j < n; j++) { min = INFINITY; for (i = 1; i < m; i++) { obj = a[i][j].abs(); if (obj.gt(ZERO) && obj.lt(min)) min = obj; } min = min.reciprocate(); for (i = 0; i < m; i++) a[i][j] = a[i][j].times(min); e[tab[j]] = min; } if (nor1 && nor2) { /* Normalise a second time */ if (output) str += "Die 2. Normierung wird durchgeführt.\n"; for (i = 1; i < m; i++) { sca = ZERO; for (j = 1; j < n; j++) sca = sca.add(a[i][j].times(a[i][j])); sca = bigRat(Math.sqrt(sca)); for (j = 0; j < n; j++) a[i][j] = a[i][j].divide(sca); n2[i] = sca.reciprocate(); } } } } function redundant() { /* Check whether constraint is <= 0 and increase counter by one if it is new */ var test; for (i = 1; i < m; i++) { test = bv[i]; if (mbv[test]) { for (j = 1; j < n && a[i][j].leq(ZERO); j++); if (j == n) { mbv[test] = false; if (test < ns) cn++; else cb++; } } } } function trans(row, col) { /* State that the boundary of computation was transcended */ k = row; l = col; solve = false; tableau("", 2); str += "Die Rechengrenze wurde überschritten.\n"; } function check() { /* Check wrong characters */ var error = "", input = document.ratplex.data.value + " ", rest = input.replace(/([\d\-\/,eE\s])/g,""), sign = "", len, rex, rows, signs, erg = input.replace(/((\s*(\-?\d{1,16},\d{1,16}((e|E)\-?\d+)?|\-?\d+((e|E)\d+)?(\/0*[1-9]\d*((e|E)\d+)?)?)\s+)+)/g, "1"); if (rest != "") { signs = rest.match(/(.)/g); signs.sort(); for (i = 0; i < signs.length; i++) { if (signs[i] != sign) { sign = signs[i]; error += sign + ", "; } } if (sign == signs[0]) alert("Das Zeichen " + sign + " ist nicht erlaubt."); else alert("Die folgenden Zeichen sind nicht erlaubt: " + error.substr(0, error.length - 2) + "."); return false; } /* Check numbers */ if (erg == "1") { input = input.replace(/(,)/g,".").replace(/([\n\r])/g,"$").replace(/((\s*\$\s*)+)/g,"$").toUpperCase(); if (input.charAt(0) == "$") input = input.substr(1); if (input.charAt(input.length - 1) == "$") input = input.substr(0, input.length - 1); rex = /(\-?\d{1,16}\.\d{1,16}(E\-?\d+)?|\-?\d+(E\d+)?(\/\d+(E\d+)?)?)/g; rows = input.split("$"); a[1] = rows[0].match(rex); m = rows.length; n = a[1].length; for (k = 2; k < m; k++) { a[k] = rows[k - 1].match(rex); /* Check matrix */ len = a[k].length; if (n != len) { if (len == 1) alert("Die Matrix hat in Zeile " + k + " leider eine Spalte statt " + n + "."); else alert("Die Matrix hat in Zeile " + k + " leider " + len + " Spalten statt " + n + "."); return false; } } if (m > 1) { a[0] = rows[k - 1].match(rex); len = a[0].length; if (n != len) { if (len == 1) alert("Die Matrix hat in Zeile " + k + " leider eine Spalte statt " + n + "."); else alert("Die Matrix hat in Zeile " + k + " leider " + len + " Spalten statt " + n + "."); return false; } } if (m < 2 || n < 2) { alert("Es sind mindestens zwei Zeilen und mindestens zwei Spalten erforderlich."); return false; } /* Transform matrix to big rationals and check its finiteness */ for (i = 0; i < m; i++) for (j = 0; j < n; j++) a[i][j] = bigRat(a[i][j]); if (finite()) { return true; } else { for (l = 1; l < n; l++) nbv[l] = l; for (k = 1; k < m; k++) bv[k] = n + k - 1; str = "Das lineare Programm hat " + m + " Zeilen und " + n + " Spalten.\n"; trans(i, j); return false; } } else { alert("Die Eingabe ist entweder nicht ausreichend oder so nicht zulässig, weil z. B. der Zahlenaufbau nicht stimmt."); return false; } } function random() { /* Create random matrix */ var len = new Array(n), maxs = new Array(n), max, mem, x; m = document.ratplex.m.selectedIndex + 2; n = document.ratplex.n.selectedIndex + 2; document.ratplex.x.selectedIndex = 0; for (i = 0; i < m; i++) { a[i] = new Array(n); for (j = 0; j < n; j++) { x = Math.random(); a[i][j] = (x > 0.5) ? Math.round(FACTOR * x) : -Math.round(FACTOR * x); } } a[0][0] = 0; str = ""; /* Format table and compute linear programme */ for (j = 0; j < n; j++) { len[j] = new Array(m); max = 0; for (i = 0; i < m; i++) { mem = (isFinite(a[i][j])) ? a[i][j].toString().length : (a[i][j] > 0) ? 1 : 2; if (mem > max) max = mem; len[j][i] = mem; } maxs[j] = ++max; } for (i = 1; i < m; i++) { for (j = 0; j < n; j++) str += space(maxs[j] - len[j][i]) + format(a[i][j]); str += "\n"; } for (j = 0; j < n; j++) str += space(maxs[j] - len[j][0]) + format(a[0][j]); str += "\n"; document.ratplex.data.value = str; compute(); } function pivot(inc, mem, v) { var ele = a[0][v], dif, pro, quo, sca, sum; switch (rule) { case 1 : /* Compute pivot by steepest edge and maybe by best value */ sum = bigRat(1); for (i = 1; i < m; i++) sum = sum.add(a[i][v].times(a[i][v])); pro = ele.times(ele); quo = pro.divide(sum); if (quo.gt(max)) { k = mem; l = v; max = quo; } else { if (quo.eq(max)) { pro = ele.times(inc); if (pro.gt(bet)) { k = mem; l = v; bet = pro; } } } break; case 2 : /* Compute pivot by best value and maybe by steepest edge */ pro = ele.times(inc); if (pro.gt(max)) { k = mem; l = v; max = pro; } else { if (pro.eq(max)) { sum = bigRat(1); for (i = 1; i < m; i++) sum = sum.add(a[i][v].times(a[i][v])); pro = ele.times(ele); quo = pro.divide(sum); if (quo.gt(bet)) { k = mem; l = v; bet = quo; } } } break; case 3 : /* Compute pivot by smallest angle and maybe by best value */ quo = ele.divide(a[mem][v]); sum = quo; sca = quo.times(quo); for (j = 1; j < n; j++) { pro = a[mem][j].times(quo); dif = a[0][j].subtract(pro); sum = sum.add(dif); pro = dif.times(dif); sca = sca.add(pro); } sum = sum.times(sum.abs()); quo = sum.divide(sca); if (quo.lt(min)) { k = mem; l = v; min = quo; } else { if (quo.eq(min)) { pro = ele.times(inc); if (pro.gt(max)) { k = mem; l = v; max = pro; } } } break; } } function transform(val) { /* Compute right-hand side */ var rowel, pivot = a[k][l].reciprocate(), pline = new Array(), i, j, r; a[k][l] = pivot; rowel = pivot.times(a[k][0]).negate(); if (val == 0) { for (i = 0; i < k; i++) a[i][0] = a[i][0].add(rowel.times(a[i][l])); for (i = k + 1; i < m; i++) a[i][0] = a[i][0].add(rowel.times(a[i][l])); } else { for (i = 0; i < val; i++) { r = s[i]; a[r][0] = a[r][0].add(rowel.times(a[r][l])); } } /* Compute pivot row */ a[k][0] = rowel.negate(); for (j = 1; j < l; j++) { pline[j] = a[k][j]; a[k][j] = pivot.times(a[k][j]); } for (j = l + 1; j < n; j++) { pline[j] = a[k][j]; a[k][j] = pivot.times(a[k][j]); } pivot = pivot.negate(); /* Compute pivot column and apply rectangle rule */ for (i = 0; i < k; i++) { a[i][l] = pivot.times(a[i][l]); rowel = a[i][l]; for (j = 1; j < l; j++) a[i][j] = a[i][j].add(rowel.times(pline[j])); for (j = l + 1; j < n; j++) a[i][j] = a[i][j].add(rowel.times(pline[j])); } for (i = k + 1; i < m; i++) { a[i][l] = pivot.times(a[i][l]); rowel = a[i][l]; for (j = 1; j < l; j++) a[i][j] = a[i][j].add(rowel.times(pline[j])); for (j = l + 1; j < n; j++) a[i][j] = a[i][j].add(rowel.times(pline[j])); } /* Check redundancy */ if (red || few) redundant(); } function tableau(string, piv) { /* Format tableau */ var i, j, len = new Array(), maxs = new Array(), nbvl = new Array(), max = TAB, mem; if (opt == 1) a[0][0] = a[0][0].negate(); str += string; len[n] = new Array(m); nbvl[n] = max; /* Determine length of entries */ for (i = 1; i < m; i++) { mem = String(bv[i]).length; if (mem > max) max = mem; len[n][i] = mem; } maxs[n] = ++max; if (l > 0) { len[0] = new Array(m); max = TAB; nbvl[0] = max; for (i = 0; i < m; i++) { tmp = a[i][0]; mem = format(tmp).length; if (mem > max) max = mem; len[0][i] = mem; } maxs[0] = ++max; } for (j = 1; j < l; j++) { len[j] = new Array(m); max = String(nbv[j]).length; nbvl[j] = max; for (i = 0; i < m; i++) { tmp = a[i][j]; mem = format(tmp).length; if (mem > max) max = mem; len[j][i] = mem; } maxs[j] = ++max; } len[l] = new Array(m); max = (l > 0) ? String(nbv[l]).length : TAB; nbvl[l] = max; for (i = 0 ; i < k; i++) { tmp = a[i][l]; mem = format(tmp).length; if (mem > max) max = mem; len[l][i] = mem; } tmp = a[k][l]; mem = format(tmp).length + piv; if (mem > max) max = mem; len[l][k] = mem; for (i = k + 1; i < m; i++) { tmp = a[i][l]; mem = format(tmp).length; if (mem > max) max = mem; len[l][i] = mem; } maxs[l] = ++max; for (j = l + 1; j < n; j++) { len[j] = new Array(m); max = String(nbv[j]).length; nbvl[j] = max; for (i = 0; i < m; i++) { tmp = a[i][j]; mem = format(tmp).length; if (mem > max) max = mem; len[j][i] = mem; } maxs[j] = ++max; } /* Output entries */ str += space(maxs[n] - TAB) + " BV:"; str += space(maxs[0] - TAB) + "RS\\NBV:"; for (j = 1; j < n; j++) str += space(maxs[j] - nbvl[j]) + nbv[j]; str += "\n"; for (i = 1; i < k; i++) { str += space(maxs[n] - len[n][i]) + bv[i]; for (j = 0; j < n; j++) str += space(maxs[j] - len[j][i]) + format(a[i][j]); str += "\n"; } if (k > 0) { str += space(maxs[n] - len[n][k]) + bv[k]; for (j = 0; j < l; j++) str += space(maxs[j] - len[j][k]) + format(a[k][j]); str += (piv == 2) ? space(maxs[l] - len[l][k]) + "»" + format(a[k][l]) + "«" : space(maxs[l] - len[l][k]) + format(a[k][l]); for (j = l + 1; j < n; j++) str += space(maxs[j] - len[j][k]) + format(a[k][j]); str += "\n"; } for (i = k + 1; i < m; i++) { str += space(maxs[n] - len[n][i]) + bv[i]; for (j = 0; j < n; j++) str += space(maxs[j] - len[j][i]) + format(a[i][j]); str += "\n"; } if (prog == 1) str += (opt == 1 ) ? space(maxs[n] - TAB) + "PP Max:" : space(maxs[n] - TAB) + "PP Min:"; else str += (opt == 1 ) ? space(maxs[n] - TAB) + "DP Min:" : space(maxs[n] - TAB) + "DP Max:"; for (j = 0; j < l; j++) str += space(maxs[j] - len[j][0]) + format(a[0][j]); str += (piv == 2 && k == 0) ? space(maxs[l] - len[l][0]) + "»" + format(a[0][l]) + "«" : space(maxs[l] - len[l][0]) + format(a[0][l]); for (j = l + 1; j < n; j++) str += space(maxs[j] - len[j][0]) + format(a[0][j]); str += "\n"; if (opt == 1) a[0][0] = a[0][0].negate(); } function loop(last) { var t = new Array(), u = new Array(), o, p, q, r, v, ele, inc, man, mem, obj, per, pro, quo, sca; l = 0; nul = 0; do { k = 0; q = 0; /* Determine the coefficients of the objective function that are greater than zero and their indices */ for (j = 1; j < n; j++) if (a[0][j].gt(ZERO)) u[q++] = j; if (q > 0) { o = 0; bet = ZERO; max = ZERO; pro = ZERO; min = INFINITY; /* Separate right-hand sides that are zero from the other ones */ if (nul == 0) { p = 0; for (i = 1; i < m; i++) if (a[i][0].eq(ZERO)) s[o++] = i; else t[p++] = i; } else { for (i = 1; i < m; i++) if (a[i][0].eq(ZERO)) o++; } if (o > 0 || nul > 0) { if (nul == 0) nul = o; /* Compute scalar products if needed */ if (o > 0) { for (i = 0; i < nul && solve; i++) { sca = ZERO; r = s[i]; for (j = 1; j < n; j++) sca = sca.add(a[r][j].times(a[r][j])); a[r][0] = a[r][0].add(bigRat(Math.sqrt(sca))); if (!isFinite(a[r][0])) { trans(r, 0); return; } } } /* Determine smallest step width if some right-hand sides are zero */ for (j = 0; j < q; j++) { v = u[j]; for (i = 0; i < nul && a[s[i]][v].leq(ZERO); i++); if (i == nul) { inc = INFINITY; for (i = 0; i < p; i++) { r = t[i]; ele = a[r][v]; if (ele.gt(ZERO)) { quo = a[r][0].divide(ele); if (quo.lt(inc)) { inc = quo; mem = r; } } } if (inc.lt(INFINITY)) pivot(inc, mem, v); } } /* Undo perturbations */ if (k > 0) { for (i = 0; i < nul; i++) a[s[i]][0] = ZERO; nul = 0; } } /* Determine smallest step width if all right-hand sides are not zero */ else { for (j = 0; j < q; j++) { v = u[j]; inc = INFINITY; for (i = 1; i < m; i++) { ele = a[i][v]; if (ele.gt(ZERO)) { quo = a[i][0].divide(ele); if (quo.lt(inc)) { inc = quo; mem = i; } } } if (inc.lt(INFINITY)) pivot(inc, mem, v); } } /* Determine smallest step width if the objective function could not be increased */ if (nul > 0) { bet = ZERO; max = ZERO; min = INFINITY; for (j = 0; j < q; j++) { v = u[j]; inc = INFINITY; for (i = 0; i < nul; i++) { r = s[i]; ele = a[r][v]; if (ele.gt(ZERO)) { quo = a[r][0].divide(ele); if (quo.lt(inc)) { inc = quo; mem = r; } } } if (inc.lt(INFINITY)) pivot(inc, mem, v); } per = (k > 0) ? ++perts : false; } else { per = false; } /* Transform tableau and build details if demanded */ if (k > 0) { mem = bv[k]; bv[k] = nbv[l]; nbv[l] = mem; steps++; if (output) { tableau("", 2); str += "Das Pivotelement ist a[" + k + "][" + (l + 1) + "] = " + format(a[k][l]) + ".\n"; if (per) str += "Eine Perturbation war erforderlich.\n"; } man = a[0][0]; transform(nul); if (last && a[0][0].neq(man)) return; if (!solve) { trans(i, j); k = 0; } } } } while (k > 0); /* Undo perturbations */ if (nul > 0 && solve) { if (output) tableau("", 0); for (i = 0; i < nul; i++) a[s[i]][0] = ZERO; if (output) str += "Tableau ohne Perturbationen:\n"; } } function numbers(string) { /* Prepare pivot steps and perturbations for output */ if (steps == 1) str += "Es wurden" + string + " ein Pivotschritt und "; else str += "Es wurden" + string + " " + steps + " Pivotschritte und "; if (perts == 0) str += "keine Perturbationen benötigt.\n"; else str += (perts == 1) ? "eine Perturbation benötigt.\n" : perts + " Perturbationen benötigt.\n"; } function compute() { function results() { /* Collect results */ var br = new Array(ms - 1), ng = new Array(); numbers(""); if (prog == 1) str += (opt == 1) ? "Das Maximum liegt bei " + format(a[0][0].negate()) + " (" + decimal(a[0][0].negate()) + ")." : "Das Minimum liegt bei " + format(a[0][0]) + " (" + decimal(a[0][0]) + ")."; else str += (opt == 1) ? "Das Minimum liegt bei " + format(a[0][0].negate()) + " (" + decimal(a[0][0].negate()) + ")." : "Das Maximum liegt bei " + format(a[0][0]) + " (" + decimal(a[0][0]) + ")."; for (j = 1; j < ns; j++) x[j] = ZERO; for (i = 1; i < ms; i++) { y[i] = ZERO; if (bv[i] < ns) x[bv[i]] = a[i][0]; } if (scal) for (j = 1; j < n; j++) x[tab[j]] = x[tab[j]].times(e[tab[j]]); /* Determine first primal solution for minimised number of constraints */ if (few && ble) { for (i = m; i < ms; i++) { bv[i] -= ns - 1; br[i] = b[bv[i]][0]; for (j = 1; j < ns; j++) br[i] = br[i].subtract(b[bv[i]][j].times(x[j])); ng[i] = br[i].lt(ZERO); } for (j = n; j < ns; j++) { min = ZERO; for (i = m; i < ms; i++) { if (ng[i] && b[bv[i]][nbv[j]].lt(ZERO)) { tmp = br[i].divide(b[bv[i]][nbv[j]]); if (tmp.gt(min)) { min = tmp; ng[i] = false; } } } x[nbv[j]] = min; } for (i = m; i < ms; i++) bv[i] += ns - 1; } /* Output renormalised solutions */ str += (prog == 1) ? "\n1. primale rationale Lösung: " : "\n1. duale rationale Lösung: "; for (j = 1; j < ns; j++) { str += format(x[j]) + " "; if (nbv[j] >= ns) y[nbv[j] - ns + 1] = a[0][j].negate(); } str += (prog == 1) ? "\n1. primale Fließkommalösung: " : "\n1. duale Fließkommalösung: "; for (j = 1; j < ns; j++) str += decimal(x[j]) + " "; if (nor1 || nor2) for (i = 1; i < m; i++) y[i] = y[i].times(n1[i]); if (nor1 && scal && nor2) for (i = 1; i < m; i++) y[i] = y[i].times(n2[i]); if (ms == 1) { str += (prog == 1) ? "\n1. duale rationale Lösung: Existiert nicht.\n1. duale Fließkommalösung: Existiert nicht." : "\n1. primale rationale Lösung: Existiert nicht.\n1. primale Fließkommalösung: Existiert nicht."; } else { str += (prog == 1) ? "\n1. duale rationale Lösung: " : "\n1. primale rationale Lösung: "; for (i = 1; i < ms; i++) str += format(y[i]) + " "; str += (prog == 1) ? "\n1. duale Fließkommalösung: " : "\n1. primale Fließkommalösung: "; for (i = 1; i < ms; i++) str += decimal(y[i]) + " "; } if (m < ms || less) str += "\nDie Restriktionsanzahl wurde verringert.\n"; else str += "\n"; } /* Insert examples */ var index = document.ratplex.x.selectedIndex; if (index > 0) { switch (index) { case 1 : document.ratplex.data.value = "15 17 -12 -2 15 -6 -4 21\n16 15 -13 -4 -1 -7 -8 0\n-8 -12 -3 -5 16 -6 22 -4\n17 -3 -12 -1 -1 21 14 0\n-7 -9 14 -11 23 -12 21 21\n21 -11 19 -2 -9 -2 16 14\n 0 -2 -3 22 16 16 -3 -11\n14 -10 17 -6 23 23 -12 24\n 0 23 23 24 -6 18 -5 -12"; break; case 2 : document.ratplex.data.value = "30 1 3 2\n25 4 2 1\n45 3 4 4\n50 2 3 5\n 0 5 8 4"; break; case 3 : document.ratplex.data.value = " 1 1 -8 -2 0 0\n 28 2 8 -2 0 0\n 42 1 -6 1 0 0\n 1 0 1 0 0 0\n0,6 0 0 1 -1 0\n8,5 -1 8 2 0 -1\n 1 -1 8 1 0 0\n0,1 -1 -9 2 0 0\n 0 1 0 0 0 0"; break; case 4 : document.ratplex.data.value = " 1 0 1 -8 0 -2 0 0 0 0\n 28 0 2 8 0 -2 0 0 0 0\n 42 0 1 -6 0 1 0 0 0 0\n 1 0 0 1 0 0 0 0 0 0\n 0,6 0 0 0 0 1 0 -1 0 0\n 8,5 0 -1 8 0 2 0 0 0 -1\n 1 0 -1 8 0 1 0 0 0 0\n 0,1 0 -1 -9 0 2 0 0 0 0\n 1 1 -8 0 -2 0 0 0 0 0\n 700 2 8 0 -2 0 0 0 0 0\n1050 1 -6 0 1 0 0 0 0 0\n 0,6 0 0 0 1 0 -1 0 0 0\n 8,5 -1 8 0 2 0 0 0 -1 0\n 1 -1 8 0 1 0 0 0 0 0\n 0,1 -1 -9 0 2 0 0 0 0 0\n 0 0 1 0 0 0 0 0 0 0"; break; case 5 : document.ratplex.data.value = "-2 -2 1 0 0 -1\n 1 -1 -2 0 0 -2\n-5 0 0 1 -1 -1\n 1 0 0 -1 1 -2\n 0 -4 -6 -2 -2 -4"; break; case 6 : document.ratplex.data.value = "0 -2 -9 1 9\n0 1/3 1 -1/3 -2\n1 1 1 1 1\n0 2 3 -1 -12"; break; case 7 : document.ratplex.data.value = "-4 -2 -1 0 0\n-6 1 -2 0 0\n-2 0 0 1 -1\n-2 0 0 -1 1\n-4 -1 -2 -1 -2\n 0 -2 1 -5 1"; break; case 8 : document.ratplex.data.value = " 4 1 1 -1\n -4 -1 -1 -1\n-0,5 -1 -1 0\n 3 -1 0 0\n 3 0 -1 0\n 3 0 0 -1\n 0 1 2 1"; break; case 9 : document.ratplex.data.value = " 5 1 0 0 0 0 0 0 0\n 25 4 1 0 0 0 0 0 0\n 125 8 4 1 0 0 0 0 0\n 625 16 8 4 1 0 0 0 0\n 3125 32 16 8 4 1 0 0 0\n 15625 64 32 16 8 4 1 0 0\n 78125 128 64 32 16 8 4 1 0\n390625 256 128 64 32 16 8 4 1\n 0 128 64 32 16 8 4 2 1"; break; case 10 : document.ratplex.data.value = " 1 1 0 0 0 0 0 0 0\n 100 20 1 0 0 0 0 0 0\n 10000 200 20 1 0 0 0 0 0\n 1000000 2000 200 20 1 0 0 0 0\n 100000000 20000 2000 200 20 1 0 0 0\n 10000000000 200000 20000 2000 200 20 1 0 0\n 1000000000000 2000000 200000 20000 2000 200 20 1 0\n100000000000000 20000000 2000000 200000 20000 2000 200 20 1\n 0 10000000 1000000 100000 10000 1000 100 10 1"; break; } } str = ""; /* Check input */ if (check()) { /* Build tableau, dualise if requested and index (non-) base variables */ k = 0; l = 0; if (Number(document.ratplex.opt.value) > 0) opt = Number(document.ratplex.opt.value); if (Number(document.ratplex.op.value) > 0) op = Number(document.ratplex.op.value); if (Number(document.ratplex.prog.value) > 0) prog = Number(document.ratplex.prog.value); if (opt == 1) a[0][0] = a[0][0].negate(); else for (j = 1; j < n; j++) a[0][j] = a[0][j].negate(); if (op == 2) for (i = 1; i < m; i++) for (j = 0; j < n; j++) a[i][j] = a[i][j].negate(); if (prog == 2) { /* Transpose matrix */ opt = (opt == 1) ? 2 : 1; min = Math.min(m, n); for (i = 0; i < min; i++) { for (j = 0; j <= i; j++) { tmp = a[i][j]; a[i][j] = a[j][i].negate(); a[j][i] = tmp.negate(); } } if (m != n) { if (m < n) { for (j = m; j < n; j++) { a[j] = new Array(m); for (i = 0; i < m; i++) { a[j][i] = a[i][j].negate(); } } } else { for (i = n; i < m; i++) { for (j = 0; j < n; j++) { a[j][i] = a[i][j].negate(); } } } tmp = m; m = n; n = tmp; } } for (j = 1; j < n; j++) nbv[j] = j; for (i = 1; i < m; i++) bv[i] = n + i - 1; str += "Das lineare Programm hat " + m + " Zeilen und " + n + " Spalten.\n"; /* Check whether the linear programme can be solved for simple constraints */ for (i = 1; i < m && solve; i++) { for (j = 1; j < n && a[i][j].geq(ZERO); j++); if (j == n && a[i][0].lt(ZERO)) { k = i; solve = false; tableau("", 2); str += "Das lineare Programm ist nicht lösbar.\n"; } } /* Initialise variables, put normalised rows into new matrix and determine trivial redundant constraints if any */ var abs, ele, f = new Array(), g = new Array(), ar = new Array(), cj = 0, eq = new Array(), le = new Array(), h, hco, hsi, ind, leng = new Array(), maxs = new Array(), lco, llim, lsi, mem, none = true, one, phase1 = false, q = 0, r, ri = 0, sca, too, ulim; red = document.ratplex.red.checked; few = document.ratplex.few.checked; less = false; if (red || few) { b[0] = new Array(); for (j = 0; j < n; j++) b[0][j] = a[0][j]; for (i = 1; i < m; i++) { b[i] = new Array(); b[i][0] = a[i][0]; max = ZERO; sca = ZERO; for (j = 1; j < n; j++) { b[i][j] = a[i][j]; abs = b[i][j].abs(); if (abs.gt(max)) max = abs; } max = (max.eq(ZERO)) ? bigRat(1) : max.reciprocate(); b[i][0] = b[i][0].times(max); for (j = 1; j < n; j++) { b[i][j] = b[i][j].times(max); sca = sca.add(b[i][j].times(b[i][j])); } sca = (sca.gt(ZERO)) ? bigRat(1 / Math.sqrt(sca)) : bigRat(1); for (j = 0; j < n; j++) b[i][j] = b[i][j].times(sca); f[i] = b[i][0].geq(ZERO); if (f[i]) { for (j = 1; j < n && b[i][j].leq(ZERO); j++); f[i] = j == n; } } /* Compare rows and compute upper and lower limits to remove further redundant rows if possible */ if (m > 2) { for (i = 1; i < m; i++) { hsi = b[i][0]; hgi = hsi.gt(ZERO); hli = hsi.lt(ZERO); hge = hsi.geq(ZERO); hle = hsi.leq(ZERO); heq = hsi.eq(ZERO); for (k = i + 1; !f[i] && k < m; k++) { if (!f[k]) { lsi = b[k][0]; if (hle && lsi.geq(ZERO)) l = i; else if (hge && lsi.leq(ZERO)) l = k; else l = 0; h = (l == i) ? k : i; if (l != 0) { llim = ZERO; ulim = INFINITY; } else { l = k; if (hgi) { llim = ZERO; ulim = hsi.divide(lsi); } else { llim = hsi.divide(lsi); ulim = INFINITY; } } one = true; for (j = 1; one && j < n; j++) { hco = b[h][j]; lco = b[l][j]; hgo = hco.gt(ZERO); hlo = hco.lt(ZERO); lgo = lco.gt(ZERO); llo = lco.lt(ZERO); if (hgo && lgo || hlo && llo) { hco = hco.divide(lco); if (lgo && llim.lt(hco)) llim = hco; else if (llo && ulim.gt(hco)) ulim = hco; if (llim.gt(ulim.add(EPSILON))) one = false; } else { if (hgo && lco.leq(ZERO)) one = false; else if (hco.geq(ZERO) && llo) one = false; } } f[h] = one; /* Alter probable redundant row if necessary and compute again */ if (!one && (hgi && lsi.gt(ZERO) || hli && lsi.lt(ZERO) || heq && lsi.eq(ZERO))) { llim = (hli) ? lsi.divide(hsi) : ZERO; ulim = (hgi) ? lsi.divide(hsi) : INFINITY; one = true; for (j = 1; one && j < n; j++) { hco = b[h][j]; lco = b[l][j]; hgo = hco.gt(ZERO); hlo = hco.lt(ZERO); lgo = lco.gt(ZERO); llo = lco.lt(ZERO); if (hgo && lgo || hlo && llo) { lco = lco.divide(hco); if (hgo && llim.lt(lco)) llim = lco; else if (hlo && ulim.gt(lco)) ulim = lco; if (llim.gt(ulim.add(EPSILON))) one = false; } else { if (hlo && lco.leq(ZERO)) one = false; else if (hco.leq(ZERO) && lgo) one = false; } } f[l] = one; } } } } } /* Output redundant rows if any and remove them from the initial matrix */ for (i = 1; i < m; i++) if (f[i]) g[q++] = i; if (q > 0) { str += (q == 1) ? "Folgende Zeile ist" : "Folgende Zeilen sind"; str += " überflüssig:\n"; for (j = 0; j < n; j++) { leng[j] = new Array(q); max = 0; for (i = 0; i < q; i++) { tmp = a[g[i]][j]; mem = (isFinite(tmp)) ? tmp.toString().length : (tmp > 0) ? 1 : 2; if (mem > max) max = mem; leng[j][i] = mem; } maxs[j] = ++max; } for (i = 0; i < q; i++) { ind = g[i]; for (j = 0; j < n; j++) str += space(maxs[j] - leng[j][i]) + decimal(a[ind][j]); str += "\n"; } q = 1; for (i = 1; i < m; i++) { if (!f[i]) { for (j = 0; j < n; j++) a[q][j] = a[i][j]; q++; } } m = q; less = true; } cb = 0; cn = 0; for (i = 0; i < m + n - 1; i++) mbv[i] = true; } /* Determine non-positive columns */ ms = m; ns = n; beq = false; ble = false; for (i = 1; i < m; i++) ar[i] = false; for (j = 1; j < n; j++) { eq[j] = a[0][j].eq(ZERO); if (eq[j]) { for (i = 1; i < m && a[i][j].eq(ZERO); i++); eq[j] = i == m; if (eq[j] && !beq) beq = cj = j; else { for (i = 1; i < m && a[i][j].leq(ZERO); i++); le[j] = i == m; if (le[j]) { if (!ble) ble = cj = j; for (i = 1; i < m; i++) if (a[i][j].lt(ZERO)) ar[i] = true; } } } } if (few) { /* Put non-positive columns to the right-hand side of the matrix b */ if (ble || beq) { q = 1; r = n; for (j = 1; j < n; j++) { if (le[j] || eq[j]) { r--; for (i = 0; i < m; i++) b[i][r] = a[i][j]; nbv[r] = j; } else { for (i = 0; i < m; i++) a[i][q] = a[i][j]; tab[q] = j; nbv[q++] = j; } } for (j = r; j < n; j++) for (i = 0; i < m; i++) a[i][j] = b[i][j]; n = q; nbv[0] = 0; /* Put the rows to compute on top of matrix a by using matrix b */ if (ble) { q = 1; r = ms; for (i = 1; i < ms; i++) { if (ar[i]) { r--; for (j = 0; j < ns; j++) b[r][j] = a[i][j]; bv[r] = i + ns - 1; } else { for (j = 0; j < ns; j++) a[q][j] = a[i][j]; bv[q++] = i + ns - 1; } } for (i = r; i < ms; i++) for (j = 0; j < ns; j++) a[i][j] = b[i][j]; m = q; } q = ns - 1; for (i = 1; i < ms; i++) for (j = 0; j < ns; j++) b[bv[i] - q][nbv[j]] = bigRat(a[i][j]); } else { for (i = 1; i < m; i++) for (j = 0; j < n; j++) b[i][j] = a[i][j]; for (j = 1; j < n; j++) tab[j] = j; } } /* Determine index for second solution */ else { for (j = 1; j < n; j++) tab[j] = j; } for (i = 1; i < m && ri == 0; i++) { if (a[i][0].eq(ZERO)) { for (j = 1; j < n && a[i][j].geq(ZERO); j++); if (j == n) ri = i; } } if (solve) { /* Check whether linear programme is limited */ var allp, alls, c = new Array(n), index, limit = true, p = n, u = new Array(n), v = new Array(m), w = new Array(n), min = ZERO, obj, x = new Array(ns - 1), y = new Array(ms - 1); mem = a[0][0]; for (j = 1; j < n && limit; j++) { if (a[0][j].gt(ZERO)) { for (i = 1; i < m && a[i][j].lt(ZERO); i++); limit = i < m; } } if (limit && m > 1) { /* Initialise variables and determine first negative value of the right-hand side if any */ if (Number(document.ratplex.rule.value) > 0) rule = Number(document.ratplex.rule.value); output = document.ratplex.output.checked; nor1 = document.ratplex.nor1.checked; scal = document.ratplex.scal.checked; nor2 = document.ratplex.nor2.checked; perts = 0; steps = 0; for (i = 1; i < m && a[i][0].geq(ZERO); i++); if (i < m) { /* Phase I: Add slack variable in the last column and determine minimal objective coefficient */ phase1 = true; k = 0; l = 0; if (output) tableau("", 0); norm_scale(); min = ZERO; for (i = 1; i < m; i++) { a[i][n] = bigRat(-1); ele = a[i][0]; if (ele.lt(min)) { k = i; min = ele; } } /* Save objective function and make tableau valid by pivot step */ for (j = 1; j < n; j++) { c[j] = a[0][j]; a[0][j] = ZERO; none &= c[j].eq(ZERO); } nbv[n] = bv[k]; bv[k] = 0; a[0][0] = ZERO; c[0] = ZERO; a[0][n] = bigRat(-1); l = n++; steps++; if (output) { tableau("Phase I:\n", 2); str += "Das Pivotelement ist a[" + k + "][" + (l + 1) + "] = -1/1.\n"; } transform(0); if (solve) { /* Loop, check feasibility and collect results */ index = k; loop(false); if (solve) { if (a[0][0].abs().gt(EPSILON)) { solve = false; // infeasible k = 0; l = 0; tableau("", 2); str += "Das lineare Programm ist nicht lösbar.\n"; numbers(""); } else { k = index; if (output) tableau("", 0); if (nbv[l] != 0) { /* Move slack variable from the base variables to the last column of the non-base variables and drop last column */ for (l = n - 1; l > 0 && a[k][l].eq(ZERO); l--); bv[k] = nbv[l]; nbv[l] = 0; steps++; transform(0); if (output) tableau("Das Pivotelement ist a[" + k + "][" + (l + 1) + "] = " + format(a[k][l].reciprocate()) + ".\n", 2); } if (solve) { phase1 = false; n--; if (nbv[l] == 0 && l < n) { /* Replace column in tableau with slack variable by last column and drop the latter */ for (i = 1; i < m; i++) a[i][l] = a[i][n]; nbv[l] = nbv[n]; } a[0][0] = mem; /* Finish phase I */ if (none) { l = 0; tableau("Ende Phase I:\n", 0); results(); } else { /* Translate objective function into new variables and check redundancy */ none = false; for (j = 1; j < n; j++) { l = nbv[j]; a[0][j] = (l < ns) ? c[l] : ZERO; } for (i = 1; i < m; i++) { l = bv[i]; if (l < ns) { obj = c[l].negate(); for (j = 0; j < n; j++) a[0][j] = a[0][j].add(obj.times(a[i][j])); } } for (j = 0; j < n; j++) if (a[0][j].abs().lt(EPSILON)) a[0][j] = ZERO; /* Check redundancy */ if (red || few) redundant(); if (output) str += "Ende Phase I:\n"; } } else { str += "Die Rechengrenze wurde überschritten.\n"; numbers(""); } } } else { numbers(""); } } else { str += "Die Rechengrenze wurde überschritten.\n"; numbers(""); } } else { /* Normalise and scale if demanded */ none = false; k = 0; l = 0; if (output && (nor1 || nor2 || scal)) tableau("", 0); norm_scale(); } /* Phase II: Loop, check whether the linear programme is limited and collect further results */ if (solve && !none) { phase1 = false; for (j = 1; j < n && a[0][j].leq(ZERO); j++); if (j < n) { if (min.lt(ZERO)) { if (output) { numbers(""); str += "Phase II:\n"; } allp = perts; alls = steps; perts = 0; steps = 0; } loop(false); if (solve) { /* Check whether linear programme is limited */ for (j = 1; j < n && a[0][j].leq(ZERO); j++); if (j < n) { solve = false; l = j; tableau("", 2); str += "Das lineare Programm ist unbeschränkt in Spalte " + ++j + ".\n"; numbers(""); } else { tableau("", 0); results(); } } else { numbers(""); } if (min.lt(ZERO)) { perts += allp; steps += alls; numbers(" insgesamt"); } } else { k = 0; l = 0; tableau("", 0); results(); } } } else { /* Output result if there are constraints */ if (m > 1 && !limit) { solve = false; k = 0; l = j - 1; tableau("", 2); str += "Das lineare Programm ist unbeschränkt in Spalte " + j + ".\n"; } if (prog == 2) { opt = (opt == 1) ? 2 : 1; } } } if (red || few) { /* Output result if there are no constraints */ if (m == 1) { for (j = 1; j < n && a[0][j].eq(ZERO); j++); if (j < n && a[0][j].gt(ZERO)) { solve = false; k = 0; l = j; tableau("", 2); str += "Das lineare Programm ist unbeschränkt in Spalte " + j + ".\n"; } else { k = 0; l = 0; perts = 0; steps = 0; tableau("", 0); results(); } } /* Output redundant rows */ m = ms; n = ns; if (cb > 0) { too = (q > 0) ? " ebenfalls" : ""; str += (cb == 1) ? "Folgende Zeile ist" : "Folgende Zeilen sind"; str += too + " überflüssig:\n"; for (j = 0; j < n; j++) { leng[j] = new Array(m); max = 0; for (i = 1; i < m; i++) { tmp = b[i][j]; mem = (isFinite(tmp)) ? tmp.toString().length : (tmp > 0) ? 1 : 2; if (mem > max) max = mem; leng[j][i] = mem; } maxs[j] = ++max; } for (i = 1; i < m; i++) { if (!mbv[i + n - 1]) { for (j = 0; j < n; j++) str += space(maxs[j] - leng[j][i]) + decimal(b[i][j]); str += "\n"; } } } if (!mbv[0]) cn--; if (cn > 0) { /* Output nonnegativity conditions */ str += (cn == 1) ? "Folgende Variable bedarf" : "Folgende Variablen bedürfen"; str += " keiner Nichtnegativitätsbedingung:\n"; for (j = 1; j < n; j++) if (!mbv[j]) str += j + " "; str += "\n"; } } if (solve) { /* Check if there are other solutions */ q = 0; mem = a[0][0]; for (i = 1; i < m; i++) { if (a[i][0].eq(ZERO)) { v[++q] = bv[i]; for (j = 1; j < n; j++) b[q][j] = a[i][j]; } } l = 0; for (j = 1; j < n; j++) { u[j] = a[0][j]; w[j] = nbv[j]; if (u[j].eq(ZERO)) { if (j > ++l) { nbv[l] = nbv[j]; for (i = 1; i < m; i++) a[i][l] = a[i][j]; } a[0][l] = bigRat(1); } } /* Compute second solution for x */ if (ble || beq) { none = false; l = 0; x[cj] = x[cj].add(bigRat(1)); } else { none = true; } if (l > 0) { n = ++l; none = red; output = false; red = false; loop(true); red = none; none = mem.eq(a[0][0]); for (i = 1; i < ns; i++) x[i] = ZERO; for (i = 1; i < ms; i++) if (bv[i] < ns) x[bv[i]] = a[i][0]; if (scal) for (i = 1; i < ns; i++) x[tab[i]] = x[tab[i]].times(e[tab[i]]); } /* Output x */ if (none) { str += (prog == 1) ? "2. primale rationale Lösung: Existiert nicht.\n2. primale Fließkommalösung: Existiert nicht." : "2. duale rationale Lösung: Existiert nicht.\n2. duale Fließkommalösung: Existiert nicht."; } else { str += (prog == 1) ? "2. primale rationale Lösung: " : "2. duale rationale Lösung: "; for (j = 1; j < ns; j++) str += format(x[j]) + " "; str += (prog == 1) ? "\n2. primale Fließkommalösung: " : "\n2. duale Fließkommalösung: "; for (j = 1; j < ns; j++) str += decimal(x[j]) + " "; } /* Compute second solution for y */ if (ri > 0) { none = false; q = 0; y[ri] = y[ri].add(bigRat(1)); } else { none = true; } if (q > 0) { m = p; n = ++q; a = new Array(m); a[0] = new Array(n); a[0][0] = mem; for (i = 1; i < m; i++) { a[i] = new Array(n); a[i][0] = u[i].negate(); bv[i] = w[i]; for (j = 1; j < n; j++) a[i][j] = b[j][i].negate(); } for (j = 1; j < n; j++) { nbv[j] = v[j]; a[0][j] = bigRat(1); } none = red; output = false; red = false; loop(true); red = none; none = mem.eq(a[0][0]); for (i = 1; i < ms; i++) y[i] = ZERO; for (i = 1; i < m; i++) if (bv[i] >= ns) y[bv[i] - ns + 1] = a[i][0]; if (nor1 || nor2) for (i = 1; i < ms; i++) y[i] = y[i].times(n1[i]); if (nor1 && scal && nor2) for (i = 1; i < ms; i++) y[i] = y[i].times(n2[i]); } /* Output y */ if (none) { str += (prog == 1) ? "\n2. duale rationale Lösung: Existiert nicht.\n2. duale Fließkommalösung: Existiert nicht." : "\n2. primale rationale Lösung: Existiert nicht.\n2. primale Fließkommalösung: Existiert nicht."; } else { str += (prog == 1) ? "\n2. duale rationale Lösung: " : "\n2. primale rationale Lösung: "; for (i = 1; i < ms; i++) str += format(y[i]) + " "; str += (prog == 1) ? "\n2. duale Fließkommalösung: " : "\n2. primale Fließkommalösung: "; for (i = 1; i < ms; i++) str += decimal(y[i]) + " "; } } } document.ratplex.result.value = str; } /* Main programme: Build form with example as input and compute linear programme */ form(); document.ratplex.data.value = "15 17 -12 -2 15 -6 -4 21\n16 15 -13 -4 -1 -7 -8 0\n-8 -12 -3 -5 16 -6 22 -4\n17 -3 -12 -1 -1 21 14 0\n-7 -9 14 -11 23 -12 21 21\n21 -11 19 -2 -9 -2 16 14\n 0 -2 -3 22 16 16 -3 -11\n14 -10 17 -6 23 23 -12 24\n 0 23 23 24 -6 18 -5 -12"; compute();