solver.dart 18.9 KB
Newer Older
1
// Copyright 2016 The Chromium Authors. All rights reserved.
2 3 4
// Use of this source code is governed by a BSD-style license that can be
// found in the LICENSE file.

5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146
import 'dart:collection';

import 'constraint.dart';
import 'expression.dart';
import 'priority.dart';
import 'result.dart';
import 'param.dart';
import 'term.dart';

enum _SymbolType { invalid, external, slack, error, dummy, }

class _Symbol {
  const _Symbol(this.type, this.tick);
  final _SymbolType type;
  final int tick;

  @override
  String toString() {
    String typeString = 'unknown';
    switch (type) {
      case _SymbolType.invalid:
        typeString = 'i';
        break;
      case _SymbolType.external:
        typeString = 'v';
        break;
      case _SymbolType.slack:
        typeString = 's';
        break;
      case _SymbolType.error:
        typeString = 'e';
        break;
      case _SymbolType.dummy:
        typeString = 'd';
        break;
    }
    return '$typeString$tick';
  }
}

class _Tag {
  _Tag(this.marker, this.other);
  _Tag.fromTag(_Tag tag)
      : this.marker = tag.marker,
        this.other = tag.other;
  _Symbol marker;
  _Symbol other;
}

class _EditInfo {
  _Tag tag;
  Constraint constraint;
  double constant;
}

bool _isValidNonRequiredPriority(double priority) {
  return (priority >= 0.0 && priority < Priority.required);
}

typedef Result _SolverBulkUpdate(dynamic item);

bool _nearZero(double value) {
  const double epsilon = 1.0e-8;
  return value < 0.0 ? -value < epsilon : value < epsilon;
}

class _Row {
  _Row(this.constant) : this.cells = new Map<_Symbol, double>();

  _Row.fromRow(_Row row)
      : this.cells = new Map<_Symbol, double>.from(row.cells),
        this.constant = row.constant;

  final Map<_Symbol, double> cells;

  double constant = 0.0;

  double add(double value) => constant += value;

  void insertSymbol(_Symbol symbol, [double coefficient = 1.0]) {
    double val = cells[symbol] ?? 0.0;

    if (_nearZero(val + coefficient)) {
      cells.remove(symbol);
    } else {
      cells[symbol] = val + coefficient;
    }
  }

  void insertRow(_Row other, [double coefficient = 1.0]) {
    constant += other.constant * coefficient;
    other.cells.forEach((_Symbol s, double v) => insertSymbol(s, v * coefficient));
  }

  void removeSymbol(_Symbol symbol) {
    cells.remove(symbol);
  }

  void reverseSign() {
    constant = -constant;
    cells.forEach((_Symbol s, double v) => cells[s] = -v);
  }

  void solveForSymbol(_Symbol symbol) {
    assert(cells.containsKey(symbol));
    double coefficient = -1.0 / cells[symbol];
    cells.remove(symbol);
    constant *= coefficient;
    cells.forEach((_Symbol s, double v) => cells[s] = v * coefficient);
  }

  void solveForSymbols(_Symbol lhs, _Symbol rhs) {
    insertSymbol(lhs, -1.0);
    solveForSymbol(rhs);
  }

  double coefficientForSymbol(_Symbol symbol) => cells[symbol] ?? 0.0;

  void substitute(_Symbol symbol, _Row row) {
    double coefficient = cells[symbol];

    if (coefficient == null) {
      return;
    }

    cells.remove(symbol);
    insertRow(row, coefficient);
  }

  @override
  String toString() {
    StringBuffer buffer = new StringBuffer();

    buffer.write(constant);

    cells.forEach((_Symbol symbol, double value) {
      buffer.write(" + " + value.toString() + " * " + symbol.toString());
    });

    return buffer.toString();
  }
}
147 148

class Solver {
149 150 151 152 153 154 155
  final Map<Constraint, _Tag> _constraints = new Map<Constraint, _Tag>();
  final Map<_Symbol, _Row> _rows = new Map<_Symbol, _Row>();
  final Map<Variable, _Symbol> _vars = new Map<Variable, _Symbol>();
  final Map<Variable, _EditInfo> _edits = new Map<Variable, _EditInfo>();
  final List<_Symbol> _infeasibleRows = new List<_Symbol>();
  final _Row _objective = new _Row(0.0);
  _Row _artificial = new _Row(0.0);
156
  int tick = 1;
157

158 159 160 161
  /// Attempts to add the constraints in the list to the solver. If it cannot
  /// add any for some reason, a cleanup is attempted so that either all
  /// constraints will be added or none.
  Result addConstraints(List<Constraint> constraints) {
162 163
    _SolverBulkUpdate applier = (Constraint c) => addConstraint(c);
    _SolverBulkUpdate undoer = (Constraint c) => removeConstraint(c);
164

165
    return _bulkEdit(constraints, applier, undoer);
166 167
  }

168
  Result addConstraint(Constraint constraint) {
169
    if (_constraints.containsKey(constraint))
170 171
      return Result.duplicateConstraint;

172 173
    _Tag tag = new _Tag(new _Symbol(_SymbolType.invalid, 0),
        new _Symbol(_SymbolType.invalid, 0));
174

175
    _Row row = _createRow(constraint, tag);
176

177
    _Symbol subject = _chooseSubjectForRow(row, tag);
178

179
    if (subject.type == _SymbolType.invalid && _allDummiesInRow(row)) {
180 181 182 183 184 185 186
      if (!_nearZero(row.constant)) {
        return Result.unsatisfiableConstraint;
      } else {
        subject = tag.marker;
      }
    }

187
    if (subject.type == _SymbolType.invalid) {
188
      if (!_addWithArtificialVariableOnRow(row))
189 190 191 192 193 194 195 196 197
        return Result.unsatisfiableConstraint;
    } else {
      row.solveForSymbol(subject);
      _substitute(subject, row);
      _rows[subject] = row;
    }

    _constraints[constraint] = tag;

198
    return _optimizeObjectiveRow(_objective);
199 200
  }

201 202 203 204 205 206 207
  Result removeConstraints(List<Constraint> constraints) {
    _SolverBulkUpdate applier = (Constraint c) => removeConstraint(c);
    _SolverBulkUpdate undoer = (Constraint c) => addConstraint(c);

    return _bulkEdit(constraints, applier, undoer);
  }

208
  Result removeConstraint(Constraint constraint) {
209
    _Tag tag = _constraints[constraint];
210
    if (tag == null)
211 212
      return Result.unknownConstraint;

213
    tag = new _Tag.fromTag(tag);
214 215 216 217
    _constraints.remove(constraint);

    _removeConstraintEffects(constraint, tag);

218
    _Row row = _rows[tag.marker];
219 220 221
    if (row != null) {
      _rows.remove(tag.marker);
    } else {
222 223
      _Symbol leaving = _leavingSymbolForMarkerSymbol(tag.marker);
      assert(leaving != null);
224

225 226
      row = _rows.remove(leaving);
      assert(row != null);
227 228 229 230 231
      row.solveForSymbols(leaving, tag.marker);
      _substitute(tag.marker, row);
    }

    return _optimizeObjectiveRow(_objective);
232 233
  }

234 235
  bool hasConstraint(Constraint constraint) {
    return _constraints.containsKey(constraint);
236 237
  }

238 239 240 241 242 243 244
  Result addEditVariables(List<Variable> variables, double priority) {
    _SolverBulkUpdate applier = (Variable v) => addEditVariable(v, priority);
    _SolverBulkUpdate undoer = (Variable v) => removeEditVariable(v);

    return _bulkEdit(variables, applier, undoer);
  }

245
  Result addEditVariable(Variable variable, double priority) {
246
    if (_edits.containsKey(variable))
247 248
      return Result.duplicateEditVariable;

249
    if (!_isValidNonRequiredPriority(priority))
250 251 252
      return Result.badRequiredStrength;

    Constraint constraint = new Constraint(
253 254 255
      new Expression([new Term(variable, 1.0)], 0.0),
      Relation.equalTo
    );
256
    constraint.priority = priority;
257

258
    assert(addConstraint(constraint) == Result.success);
259

260
    _EditInfo info = new _EditInfo();
261 262 263 264 265 266 267
    info.tag = _constraints[constraint];
    info.constraint = constraint;
    info.constant = 0.0;

    _edits[variable] = info;

    return Result.success;
268 269
  }

270 271 272 273 274 275 276 277
  Result removeEditVariables(List<Variable> variables) {
    _SolverBulkUpdate applier = (Variable v) => removeEditVariable(v);
    _SolverBulkUpdate undoer = (Variable v) =>
        addEditVariable(v, _edits[v].constraint.priority);

    return _bulkEdit(variables, applier, undoer);
  }

278
  Result removeEditVariable(Variable variable) {
279
    _EditInfo info = _edits[variable];
280
    if (info == null)
281 282
      return Result.unknownEditVariable;

283
    assert(removeConstraint(info.constraint) == Result.success);
284 285 286

    _edits.remove(variable);
    return Result.success;
287 288
  }

289 290
  bool hasEditVariable(Variable variable) {
    return _edits.containsKey(variable);
291 292
  }

293
  Result suggestValueForVariable(Variable variable, double value) {
294
    if (!_edits.containsKey(variable))
295 296 297 298 299
      return Result.unknownEditVariable;

    _suggestValueForEditInfoWithoutDualOptimization(_edits[variable], value);

    return _dualOptimize();
300 301
  }

302 303
  Set<dynamic> flushUpdates() {
    Set<dynamic> updates = new HashSet<dynamic>();
304

305
    for (Variable variable in _vars.keys) {
306 307
      _Symbol symbol = _vars[variable];
      _Row row = _rows[symbol];
308 309 310

      double updatedValue = row == null ? 0.0 : row.constant;

311 312
      if (variable.applyUpdate(updatedValue) && variable.owner != null) {
        dynamic context = variable.owner.context;
313
        if (context != null)
314
          updates.add(context);
315 316
      }
    }
317

318
    return updates;
319
  }
320

321 322 323 324 325 326
  Result _bulkEdit(
    Iterable<dynamic> items,
    _SolverBulkUpdate applier,
    _SolverBulkUpdate undoer
  ) {
    List<dynamic> applied = <dynamic>[];
327 328 329 330 331 332 333 334 335 336 337 338 339 340 341
    bool needsCleanup = false;

    Result result = Result.success;

    for (dynamic item in items) {
      result = applier(item);
      if (result == Result.success) {
        applied.add(item);
      } else {
        needsCleanup = true;
        break;
      }
    }

    if (needsCleanup) {
342
      for (dynamic item in applied.reversed)
343 344 345 346 347 348
        undoer(item);
    }

    return result;
  }

349
  _Symbol _symbolForVariable(Variable variable) {
350
    _Symbol symbol = _vars[variable];
351

352
    if (symbol != null)
353 354
      return symbol;

355
    symbol = new _Symbol(_SymbolType.external, tick++);
356 357 358 359 360
    _vars[variable] = symbol;

    return symbol;
  }

361
  _Row _createRow(Constraint constraint, _Tag tag) {
362
    Expression expr = new Expression.fromExpression(constraint.expression);
363
    _Row row = new _Row(expr.constant);
364

365
    expr.terms.forEach((Term term) {
366
      if (!_nearZero(term.coefficient)) {
367
        _Symbol symbol = _symbolForVariable(term.variable);
368

369
        _Row foundRow = _rows[symbol];
370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385

        if (foundRow != null) {
          row.insertRow(foundRow, term.coefficient);
        } else {
          row.insertSymbol(symbol, term.coefficient);
        }
      }
    });

    switch (constraint.relation) {
      case Relation.lessThanOrEqualTo:
      case Relation.greaterThanOrEqualTo:
        {
          double coefficient =
              constraint.relation == Relation.lessThanOrEqualTo ? 1.0 : -1.0;

386
          _Symbol slack = new _Symbol(_SymbolType.slack, tick++);
387 388 389
          tag.marker = slack;
          row.insertSymbol(slack, coefficient);

390
          if (constraint.priority < Priority.required) {
391
            _Symbol error = new _Symbol(_SymbolType.error, tick++);
392 393 394 395 396 397 398
            tag.other = error;
            row.insertSymbol(error, -coefficient);
            _objective.insertSymbol(error, constraint.priority);
          }
        }
        break;
      case Relation.equalTo:
399
        if (constraint.priority < Priority.required) {
400 401
          _Symbol errPlus = new _Symbol(_SymbolType.error, tick++);
          _Symbol errMinus = new _Symbol(_SymbolType.error, tick++);
402 403 404 405 406 407 408
          tag.marker = errPlus;
          tag.other = errMinus;
          row.insertSymbol(errPlus, -1.0);
          row.insertSymbol(errMinus, 1.0);
          _objective.insertSymbol(errPlus, constraint.priority);
          _objective.insertSymbol(errMinus, constraint.priority);
        } else {
409
          _Symbol dummy = new _Symbol(_SymbolType.dummy, tick++);
410 411 412 413 414 415 416 417 418 419 420 421 422
          tag.marker = dummy;
          row.insertSymbol(dummy);
        }
        break;
    }

    if (row.constant < 0.0) {
      row.reverseSign();
    }

    return row;
  }

423 424
  _Symbol _chooseSubjectForRow(_Row row, _Tag tag) {
    for (_Symbol symbol in row.cells.keys) {
425
      if (symbol.type == _SymbolType.external)
426 427 428
        return symbol;
    }

429 430
    if (tag.marker.type == _SymbolType.slack ||
        tag.marker.type == _SymbolType.error) {
431
      if (row.coefficientForSymbol(tag.marker) < 0.0)
432 433 434
        return tag.marker;
    }

435 436
    if (tag.other.type == _SymbolType.slack ||
        tag.other.type == _SymbolType.error) {
437
      if (row.coefficientForSymbol(tag.other) < 0.0)
438 439 440
        return tag.other;
    }

441
    return new _Symbol(_SymbolType.invalid, 0);
442 443
  }

444 445
  bool _allDummiesInRow(_Row row) {
    for (_Symbol symbol in row.cells.keys) {
446
      if (symbol.type != _SymbolType.dummy)
447 448 449 450 451
        return false;
    }
    return true;
  }

452
  bool _addWithArtificialVariableOnRow(_Row row) {
453
    _Symbol artificial = new _Symbol(_SymbolType.slack, tick++);
454 455
    _rows[artificial] = new _Row.fromRow(row);
    _artificial = new _Row.fromRow(row);
456 457 458 459 460 461 462 463 464

    Result result = _optimizeObjectiveRow(_artificial);

    if (result.error) {
      // FIXME(csg): Propagate this up!
      return false;
    }

    bool success = _nearZero(_artificial.constant);
465
    _artificial = new _Row(0.0);
466

467
    _Row foundRow = _rows[artificial];
468 469
    if (foundRow != null) {
      _rows.remove(artificial);
470
      if (foundRow.cells.isEmpty)
471 472
        return success;

473
      _Symbol entering = _anyPivotableSymbol(foundRow);
474
      if (entering.type == _SymbolType.invalid)
475 476 477 478 479 480 481
        return false;

      foundRow.solveForSymbols(artificial, entering);
      _substitute(entering, foundRow);
      _rows[entering] = foundRow;
    }

482
    for (_Row row in _rows.values)
483 484 485 486 487
      row.removeSymbol(artificial);
    _objective.removeSymbol(artificial);
    return success;
  }

488
  Result _optimizeObjectiveRow(_Row objective) {
489
    while (true) {
490
      _Symbol entering = _enteringSymbolForObjectiveRow(objective);
491
      if (entering.type == _SymbolType.invalid)
492 493
        return Result.success;

494 495
      _Symbol leaving = _leavingSymbolForEnteringSymbol(entering);
      assert(leaving != null);
496

497
      _Row row = _rows.remove(leaving);
498 499 500 501 502 503
      row.solveForSymbols(leaving, entering);
      _substitute(entering, row);
      _rows[entering] = row;
    }
  }

504
  _Symbol _enteringSymbolForObjectiveRow(_Row objective) {
505
    Map<_Symbol, double> cells = objective.cells;
506

507
    for (_Symbol symbol in cells.keys) {
508
      if (symbol.type != _SymbolType.dummy && cells[symbol] < 0.0)
509 510 511
        return symbol;
    }

512
    return new _Symbol(_SymbolType.invalid, 0);
513 514
  }

515
  _Symbol _leavingSymbolForEnteringSymbol(_Symbol entering) {
516
    double ratio = double.MAX_FINITE;
517
    _Symbol result;
518
    _rows.forEach((_Symbol symbol, _Row row) {
519
      if (symbol.type != _SymbolType.external) {
520 521
        double temp = row.coefficientForSymbol(entering);
        if (temp < 0.0) {
Hixie's avatar
Hixie committed
522 523 524
          double tempRatio = -row.constant / temp;
          if (tempRatio < ratio) {
            ratio = tempRatio;
525
            result = symbol;
526 527 528 529 530 531 532
          }
        }
      }
    });
    return result;
  }

533
  void _substitute(_Symbol symbol, _Row row) {
534
    _rows.forEach((_Symbol first, _Row second) {
535
      second.substitute(symbol, row);
536
      if (first.type != _SymbolType.external && second.constant < 0.0) {
537 538 539 540
        _infeasibleRows.add(first);
      }
    });
    _objective.substitute(symbol, row);
541
    if (_artificial != null)
542 543 544
      _artificial.substitute(symbol, row);
  }

545 546
  _Symbol _anyPivotableSymbol(_Row row) {
    for (_Symbol symbol in row.cells.keys) {
547 548
      if (symbol.type == _SymbolType.slack ||
          symbol.type == _SymbolType.error) {
549 550 551
        return symbol;
      }
    }
552
    return new _Symbol(_SymbolType.invalid, 0);
553
  }
554

555
  void _removeConstraintEffects(Constraint cn, _Tag tag) {
556
    if (tag.marker.type == _SymbolType.error) {
557 558
      _removeMarkerEffects(tag.marker, cn.priority);
    }
559
    if (tag.other.type == _SymbolType.error) {
560 561 562 563
      _removeMarkerEffects(tag.other, cn.priority);
    }
  }

564 565
  void _removeMarkerEffects(_Symbol marker, double strength) {
    _Row row = _rows[marker];
566 567 568 569 570 571 572
    if (row != null) {
      _objective.insertRow(row, -strength);
    } else {
      _objective.insertSymbol(marker, -strength);
    }
  }

573
  _Symbol _leavingSymbolForMarkerSymbol(_Symbol marker) {
574 575 576
    double r1 = double.MAX_FINITE;
    double r2 = double.MAX_FINITE;

577
    _Symbol first, second, third;
578

579
    _rows.forEach((_Symbol symbol, _Row row) {
580
      double c = row.coefficientForSymbol(marker);
581
      if (c == 0.0)
582
        return;
583
      if (symbol.type == _SymbolType.external) {
584
        third = symbol;
585 586 587 588
      } else if (c < 0.0) {
        double r = -row.constant / c;
        if (r < r1) {
          r1 = r;
589
          first = symbol;
590 591 592 593 594
        }
      } else {
        double r = row.constant / c;
        if (r < r2) {
          r2 = r;
595
          second = symbol;
596 597 598 599
        }
      }
    });

600
    return first ?? second ?? third;
601
  }
602 603

  void _suggestValueForEditInfoWithoutDualOptimization(
604
      _EditInfo info, double value) {
605 606 607 608
    double delta = value - info.constant;
    info.constant = value;

    {
609 610
      _Symbol symbol = info.tag.marker;
      _Row row = _rows[info.tag.marker];
611 612 613 614 615 616 617 618 619 620 621 622 623 624 625 626 627 628 629

      if (row != null) {
        if (row.add(-delta) < 0.0) {
          _infeasibleRows.add(symbol);
        }
        return;
      }

      symbol = info.tag.other;
      row = _rows[info.tag.other];

      if (row != null) {
        if (row.add(delta) < 0.0) {
          _infeasibleRows.add(symbol);
        }
        return;
      }
    }

630 631
    for (_Symbol symbol in _rows.keys) {
      _Row row = _rows[symbol];
632 633 634
      double coeff = row.coefficientForSymbol(info.tag.marker);
      if (coeff != 0.0 &&
          row.add(delta * coeff) < 0.0 &&
635
          symbol.type != _SymbolType.external) {
636 637 638 639 640 641 642
        _infeasibleRows.add(symbol);
      }
    }
  }

  Result _dualOptimize() {
    while (_infeasibleRows.length != 0) {
643 644
      _Symbol leaving = _infeasibleRows.removeLast();
      _Row row = _rows[leaving];
645 646

      if (row != null && row.constant < 0.0) {
647
        _Symbol entering = _dualEnteringSymbolForRow(row);
648

649
        assert(entering.type != _SymbolType.invalid);
650 651 652 653 654 655 656 657 658 659 660

        _rows.remove(leaving);

        row.solveForSymbols(leaving, entering);
        _substitute(entering, row);
        _rows[entering] = row;
      }
    }
    return Result.success;
  }

661
  _Symbol _dualEnteringSymbolForRow(_Row row) {
662
    _Symbol entering;
663 664 665

    double ratio = double.MAX_FINITE;

666
    Map<_Symbol, double> rowCells = row.cells;
667

668
    for (_Symbol symbol in rowCells.keys) {
669 670
      double value = rowCells[symbol];

671
      if (value > 0.0 && symbol.type != _SymbolType.dummy) {
672 673 674 675 676 677 678 679 680
        double coeff = _objective.coefficientForSymbol(symbol);
        double r = coeff / value;
        if (r < ratio) {
          ratio = r;
          entering = symbol;
        }
      }
    }

681
    return entering ?? new _Symbol(_SymbolType.invalid, 0);
682
  }
683

684
  @override
685 686 687 688 689 690 691 692 693 694
  String toString() {
    StringBuffer buffer = new StringBuffer();
    String separator = "\n~~~~~~~~~";

    // Objective
    buffer.writeln(separator + " Objective");
    buffer.writeln(_objective.toString());

    // Tableau
    buffer.writeln(separator + " Tableau");
695 696
    _rows.forEach((_Symbol symbol, _Row row) {
      buffer.writeln('$symbol | $row');
697 698 699 700
    });

    // Infeasible
    buffer.writeln(separator + " Infeasible");
701 702 703
    _infeasibleRows.forEach((_Symbol symbol) {
      buffer.writeln(symbol);
    });
704 705 706

    // Variables
    buffer.writeln(separator + " Variables");
707 708 709
    _vars.forEach((Variable variable, _Symbol symbol) {
      buffer.writeln('$variable = $symbol');
    });
710 711 712

    // Edit Variables
    buffer.writeln(separator + " Edit Variables");
713 714 715
    _edits.forEach((Variable variable, _EditInfo editinfo) {
      buffer.writeln(variable);
    });
716 717 718

    // Constraints
    buffer.writeln(separator + " Constraints");
719 720 721
    _constraints.forEach((Constraint constraint, _Tag tag) {
      buffer.writeln(constraint);
    });
722 723 724

    return buffer.toString();
  }
725
}