1111-- 'optimizeFeasibleSystem' performs phase two of the two-phase simplex method.
1212-- 'twoPhaseSimplex' performs both phases of the two-phase simplex method.
1313-- 'twoPhaseSimplex'' performs both phases with variable domain support.
14- module Linear.Simplex.Solver.TwoPhase
14+ module Linear.Simplex.Solver.TwoPhase
1515 ( findFeasibleSolution
1616 , optimizeFeasibleSystem
1717 , twoPhaseSimplex
@@ -410,30 +410,22 @@ twoPhaseSimplex :: (MonadIO m, MonadLogger m) => VarDomainMap -> ObjectiveFuncti
410410twoPhaseSimplex domainMap objFunction constraints = do
411411 logMsg LevelInfo $
412412 " twoPhaseSimplex: Solving system with domain map " <> showT domainMap
413- let (transformedObj, transformedConstraints, transforms) = preprocess objFunction domainMap constraints
413+ let (transformedObj, transformedConstraints, transforms) = preprocess objFunction domainMap constraints
414414 logMsg LevelInfo $
415415 " twoPhaseSimplex: Applied transforms " <> showT transforms
416416 <> " ; Transformed objective: " <> showT transformedObj
417417 <> " ; Transformed constraints: " <> showT transformedConstraints
418- phase1Result <- findFeasibleSolution transformedConstraints
419- case phase1Result of
420- Nothing -> do
421- logMsg LevelInfo " twoPhaseSimplex: No feasible solution found in phase 1"
422- pure Nothing
423- Just feasibleSystem -> do
424- logMsg LevelInfo $
425- " twoPhaseSimplex: Feasible system found for transformed system; Feasible system: "
426- <> showT feasibleSystem
427- mOptimizedSystem <- optimizeFeasibleSystem transformedObj feasibleSystem
428- case mOptimizedSystem of
429- Nothing -> do
430- logMsg LevelInfo " twoPhaseSimplex: No optimized solution found in phase 2"
431- pure Nothing
432- Just result -> do
433- let finalResult = postprocess objFunction transforms result
434- logMsg LevelInfo $
435- " twoPhaseSimplex: Postprocessed result: " <> showT finalResult
436- pure (Just finalResult)
418+ mFeasibleSystem <- findFeasibleSolution transformedConstraints
419+ let phase1FailureLog = logMsg LevelInfo " twoPhaseSimplex: No feasible solution found in phase 1"
420+ let runPhase2 feasibleSystem = do
421+ logMsg LevelInfo $
422+ " twoPhaseSimplex: Feasible system found for transformed system; Feasible system: "
423+ <> showT feasibleSystem
424+ mOptimizedSystem <- optimizeFeasibleSystem transformedObj feasibleSystem
425+ let mFinalResult = postprocess objFunction transforms <$> mOptimizedSystem
426+ logMsg LevelInfo $ maybe " twoPhaseSimplex: No optimized solution found in phase 2" ((" twoPhaseSimplex: Postprocessed result: " <> ) . showT) mFinalResult
427+ pure mFinalResult
428+ maybe (phase1FailureLog >> pure Nothing ) runPhase2 mFeasibleSystem
437429
438430-- | Postprocess the result by unapplying variable transformations and computing
439431-- the objective value in the original space.
@@ -457,8 +449,8 @@ computeObjective objFunction varVals =
457449
458450-- | Preprocess the system by applying variable transformations based on domain information.
459451-- Returns the transformed objective, constraints, and the list of transforms applied.
460- preprocess :: ObjectiveFunction
461- -> VarDomainMap
452+ preprocess :: ObjectiveFunction
453+ -> VarDomainMap
462454 -> [PolyConstraint ]
463455 -> (ObjectiveFunction , [PolyConstraint ], [VarTransform ])
464456preprocess objFunction (VarDomainMap domainMap) constraints =
@@ -507,24 +499,24 @@ getTransform nextFreshVar var (Bounded mLower mUpper) =
507499 | l == 0 -> ([] , 0 ) -- NonNegative: no transform needed
508500 | l > 0 -> ([AddLowerBound var l], 0 ) -- Positive lower bound: add constraint
509501 | otherwise -> ([Shift var nextFreshVar l], 1 ) -- Negative lower bound: shift
510-
502+
511503 -- Handle upper bound (if present)
512504 upperTransforms = case mUpper of
513505 Nothing -> []
514506 Just u -> [AddUpperBound var u]
515-
507+
516508 -- If no lower bound (Nothing), we need Split transformation
517509 -- Split replaces the variable, so upper bound would apply to the original var
518510 -- which gets expressed as posVar - negVar
519511 (finalTransforms, finalOffset) = case mLower of
520- Nothing ->
512+ Nothing ->
521513 -- Unbounded: split the variable
522514 -- Note: upperTransforms will still be added and will apply to the original variable
523515 -- expression (posVar - negVar) via the constraint system
524516 (Split var nextFreshVar (nextFreshVar + 1 ) : upperTransforms, 2 )
525517 Just _ ->
526518 (lowerTransforms ++ upperTransforms, varOffset)
527-
519+
528520 in (finalTransforms, finalOffset)
529521
530522-- | Apply all transforms to the objective function and constraints.
@@ -539,19 +531,19 @@ applyTransform transform (objFunction, constraints) =
539531 -- AddLowerBound: Add a GEQ constraint for the variable
540532 AddLowerBound v bound ->
541533 (objFunction, GEQ (M. singleton v 1 ) bound : constraints)
542-
534+
543535 -- AddUpperBound: Add a LEQ constraint for the variable
544536 AddUpperBound v bound ->
545537 (objFunction, LEQ (M. singleton v 1 ) bound : constraints)
546-
538+
547539 -- Shift: originalVar = shiftedVar + shiftBy (where shiftBy < 0)
548540 -- Substitute: wherever we see originalVar, replace with shiftedVar
549541 -- and adjust the RHS by -coeff * shiftBy
550542 Shift origVar shiftedVar shiftBy ->
551543 ( applyShiftToObjective origVar shiftedVar shiftBy objFunction
552544 , map (applyShiftToConstraint origVar shiftedVar shiftBy) constraints
553545 )
554-
546+
555547 -- Split: originalVar = posVar - negVar
556548 -- Substitute: wherever we see originalVar with coeff c,
557549 -- replace with posVar with coeff c and negVar with coeff -c
@@ -585,13 +577,13 @@ applyShiftToObjective origVar shiftedVar _shiftBy objFunction =
585577applyShiftToConstraint :: Var -> Var -> SimplexNum -> PolyConstraint -> PolyConstraint
586578applyShiftToConstraint origVar shiftedVar shiftBy constraint =
587579 case constraint of
588- LEQ m rhs ->
580+ LEQ m rhs ->
589581 let (newMap, rhsAdjust) = substituteVarInMap origVar shiftedVar shiftBy m
590582 in LEQ newMap (rhs - rhsAdjust)
591- GEQ m rhs ->
583+ GEQ m rhs ->
592584 let (newMap, rhsAdjust) = substituteVarInMap origVar shiftedVar shiftBy m
593585 in GEQ newMap (rhs - rhsAdjust)
594- EQ m rhs ->
586+ EQ m rhs ->
595587 let (newMap, rhsAdjust) = substituteVarInMap origVar shiftedVar shiftBy m
596588 in EQ newMap (rhs - rhsAdjust)
597589 where
@@ -644,18 +636,18 @@ unapplyTransform transform result@(Result {varValMap = valMap, ..}) =
644636 case transform of
645637 -- AddLowerBound: No variable substitution was done, nothing to unapply
646638 AddLowerBound {} -> result
647-
639+
648640 -- AddUpperBound: No variable substitution was done, nothing to unapply
649641 AddUpperBound {} -> result
650-
642+
651643 -- Shift: originalVar = shiftedVar + shiftBy
652644 -- So originalVar's value = shiftedVar's value + shiftBy
653645 Shift origVar shiftedVar shiftBy ->
654646 let shiftedVal = M. findWithDefault 0 shiftedVar valMap
655647 origVal = shiftedVal + shiftBy
656648 newMap = M. insert origVar origVal (M. delete shiftedVar valMap)
657649 in result { varValMap = newMap }
658-
650+
659651 -- Split: originalVar = posVar - negVar
660652 -- So originalVar's value = posVar's value - negVar's value
661653 Split origVar posVar negVar ->
0 commit comments