It’s time to get down to business and write a SchedulingModel class that encapsulates the resource constrained project scheduling model I wrote about in my previous two posts. Let’s break it down into three phases: defining the state, writing the code that constructs the model, and writing a function that solves the model and returns the results. To define the state, let’s review the important data associated with the model:

- Each event has a start date (t). We’ll call this Decision
**start**. - Each task is either active or inactive for each event (z). Let’s call this Decision
**isActive**. - The project duration is represented by C_max. We’ll call this
**makespan**. - The
**duration**of each task (p) is an input Parameter. - The above Decisions and Parameters are indexed by the sets of tasks and events.

We will have private member fields for this information since they will be used by the model throughout. Here is the start of our SchedulingModel class:

using System; using System.Collections.Generic; using System.Diagnostics; using System.Linq; using System.Text; using Microsoft.SolverFoundation.Common; using Microsoft.SolverFoundation.Services; using SolverFoundation.Plugin.Gurobi; namespace ProjectScheduling { /// <summary>A Solver Foundation model for solving resource constrained project scheduling problems. /// </summary> /// <remarks> /// This model produces a schedule with minimum duration with no overallocations. /// /// The source for this model is: /// "Event-based MILP models for resource-constrained project scheduling problems" /// by Oumar Kone, Christian Artigues, Pierre Lopez, Marcel Mongeau in /// Computers and Operations Research, December 2009. This model is the "On/Off Event-Based Formulation" /// presented in section 3.2. /// </remarks> public class SchedulingModel { private bool verbose; private SolverContext context; private Model model; private Set tasks, events, events1ToN; // in the paper: sets A, R, E, and E - 0. private Decision makespan; // C_max private Decision isActive; // z private Decision start; // t private Parameter duration; // p /// <summary>Create a new instance. /// </summary> /// <param name="verbose">Indicates whether to dump intermediate output.</param> public SchedulingModel(bool verbose) { this.verbose = verbose; } /// <summary>Create the model. /// </summary> /// <param name="project">The project to be scheduled.</param> public void Initialize(Project project) { context = SolverContext.GetContext(); context.ClearModel(); model = context.CreateModel(); int eventCount = project.Tasks.Count; // we will fill these in in the remainder of the post. InitializeSets(project.Tasks.Count, eventCount, project.Resources.Count); InitializeParameters(project.Tasks, project.Resources); InitializeDecisions(); InitializeGoal(); InitializeConstraints(project, eventCount); } }

Now let’s do the important part: write the model. This amounts to transcribing equations (34) – (45) on Page 9 of the paper. You might find the remainder of the post easier to follow if you have the paper with you as you read through the code and descriptions. The first thing I like to do is think about the Sets, Parameters (inputs) and Decisions (outputs) I will need, and then write down the constraints and goals. Let’s start with Sets. Looking over the equations we can observe the following:

- Most of the equations involve a “foreach” over the set E of events or the set A of tasks (or both). In these cases we will want to use the Model.Foreach construct.
- Equations (39) – (43) use summations. The natural SFS construct to use is Model.Sum. But note some subtle variations:
- In equation (39) we sum over the events from 0 to e-1, where e itself can be any event from 1 to Events.Count. The summations in (40) and (42) similar structure.
- The summation in (41) is simple: sum over all events.
- (43) appears simple: we just sum from i = 0 to n – 1.

If you review the documentation for Model.Sum you will see that the first argument is a Set, and the second is a Func that maps a Term (which represents an element of the set) to a Term (which represents the operand to Sum). This means that it won’t work too well for the sums in equations (39), (40), (42) because the number of terms in the sum varies depending on the value of e. We can still do what we want by just computing the sum ourselves using a standard C# for-loop. Therefore we need three sets: the set of tasks, the set of events, and the set of events from 1 to n. There are two ways to establish the elements of a Set. The first works for the case where the elements are numeric: specify the start, end, and step. I will use this constructor for the event-based sets. The second is more commonly used: define a Parameter or Decision that uses the Set, and bind the Parameter or Decision to a data source using SetBinding. The elements of the Set will be determined automatically. I like this way of defining sets because it is less brittle: if the size of my problem changes I generally don’t need to touch my code. I will use this for the task set, since I know I need to bind the task duration to a Parameter.

private void InitializeSets(int taskCount, int eventCount, int resourceCount) { tasks = new Set(Domain.Real, "tasks"); events = new Set(0, eventCount, 1); events1ToN = new Set(1, eventCount, 1); // starting from event 1. }

Next come Parameters. You can identify these by looking at the “input variables” of the equations. Here they are p_I (the task duration), b_ik (amount of work on task I performed by resource k) and B_k (resource availability). B_k and b_ik appear in equation (43). Stare at it for a second, and think about the data structures we defined last time. We defined an Assignment class that represents work performed by a resource on a task. Each task has a list of assignments. Therefore we do not really have a matrix b_ik, at least not in a single collection. We’ll need to iterate over the assignments for each task to collect the resource usage and capacity – so it doesn’t really buy us anything to define Parameters for these quantities. So we’ll just define a parameter for task duration – and remember that this will also establish the elements of the tasks set we defined earlier.

private void InitializeParameters(ICollection<Task> t, ICollection<Resource> r) { duration = new Parameter(Domain.RealNonnegative, "duration", tasks); duration.SetBinding(t, "Duration", "ID"); model.AddParameters(duration); }

Now for the Decisions. Those come straight out of the paper. When we define the decisions we need to specify the Sets over which they are indexed, and the range of possible values (the Domain) for each.

private void InitializeDecisions() { makespan = new Decision(Domain.RealNonnegative, "makespan"); isActive = new Decision(Domain.IntegerRange(0, 1), "isActive", tasks, events); start = new Decision(Domain.RealNonnegative, "start", events); model.AddDecisions(makespan, isActive, start); }

The Goal is incredibly simple: minimize the makespan.

private void InitializeGoal() { Goal goal = model.AddGoal("goal", GoalKind.Minimize, makespan); // (34) }

Now the fun part: the constraints. Of course, we have already examined most of the constraints so that we could define the other parts of our model. This is not unusual! Let’s work from top to bottom, starting with (35). This constraint establishes the makespan for the project: the latest finish time for any task. The equation says the constraint holds for all events and tasks, but that isn’t quite true: notice that we compare the start of event e with the start of event (e-1). So we actually want to start from event 1.

private void InitializeConstraints(Project project, int eventCount) { // Establish the makespan: the maximum finish time over all activities. model.AddConstraint("c35", Model.ForEach(events1ToN, e => Model.ForEach(tasks, i => makespan >= start[e] + (isActive[i, e] - isActive[i, e - 1]) * duration[i] )));

(36) is trivial:

// The first event starts at time 0. model.AddConstraint("c_36", start[0] == 0);

(37) is a foreach over each task and unique pair of events. It relates the isActive and start decisions. To model the fact that we want to consider unique pairs, we place a ForeachWhere inside a ForEach. The condition on the ForeachWhere ensures that the inner event index f exceeds outer event index e.

// Link the isActive decision to the starts of each event and task durations. model.AddConstraint("c_37", Model.ForEach(events1ToN, e => Model.ForEachWhere(events, f => Model.ForEach(tasks, i => start[f] >= start[e] + ((isActive[i, e] - isActive[i, e - 1]) - (isActive[i, f] - isActive[i, f - 1]) - 1) * duration[i] ), f => f > e)));

As I alluded to in my first post, the original paper does not consider the first event, so let us add one more constraint (37a).

// Link the isActive decision with the first event. This constraint is missing in the original // paper. model.AddConstraint("c_37_e0", Model.ForEach(events1ToN, f => Model.ForEach(tasks, i => start[f] >= start[0] + (isActive[i, 0] - (isActive[i, f] - isActive[i, f - 1]) - 1) * duration[i] )));

(38) is easy.

// Order the events. model.AddConstraint("c_38", Model.ForEach(events1ToN, e => start[e] >= start[e - 1]));

(39) and (40) are very similar, so let’s focus on (39). We looked at these when we were figuring out which Sets to define. Since the sum ranges from event 0 to event e, where e varies, we decided to use nested for-loops instead of the Model.Foreach construct. The outer loop is over the event e. If we look at the summation, we sum over events 0 .. e-1. We will use the SumTermBuilder class to accumulate these terms (I blogged about this class awhile ago – it is part of Solver Foundation 3.0).

// Ensure adjacency of events for an in-progress task. SumTermBuilder sum = new SumTermBuilder(eventCount); for (int i = 0; i < project.Tasks.Count; i++) { for (int e = 1; e < eventCount; e++) { sum.Add(isActive[i, e - 1]); model.AddConstraint("c_39_" + i + "_" + e, sum.ToTerm() <= e * (1 - (isActive[i, e] - isActive[i, e - 1]))); } sum.Clear(); } sum = new SumTermBuilder(eventCount); for (int i = 0; i < project.Tasks.Count; i++) { for (int e = 1; e < eventCount; e++) { for (int e1 = e; e1 < eventCount; e1++) { sum.Add(isActive[i, e1]); // (it's inefficient to reconstruct this for each value of e.) } model.AddConstraint("c_40_" + i + "_" + e, sum.ToTerm() <= e * (1 + (isActive[i, e] - isActive[i, e - 1]))); sum.Clear(); } }

(41) is easy.

// All activities must be active during the project. model.AddConstraint("c_41", Model.ForEach(tasks, i => Model.Sum(Model.ForEach(events, e => isActive[i, e])) >= 1));

(42) is not too hard at this point either. We need to iterate over the set of task dependencies, and add a constraint involving a summation for each. We once again use the SumTermBuilder class for performance.

// A link (i, j) means that the start of task j must be after the finish of task i. int c42 = 0; foreach (TaskDependency link in project.Dependencies) { int i = link.Source.ID; int j = link.Destination.ID; sum = new SumTermBuilder(eventCount); for (int e = 0; e < eventCount; e++) { sum.Add(isActive[j, e]); // sum now has isActive[j, 0] .. isActive[j, e]. model.AddConstraint("c_42_" + c42++, isActive[i, e] + sum.ToTerm() <= 1 + e * (1 - isActive[i, e])); } }

(43) is the last constraint – we already took care of (44) and (45) in our Decision definitions. We considered this constraint in some detail when we were figuring out our Parameter definitions. We want to express the fact that the sum of resource work on all tasks is no greater than its availability. This condition needs to hold for each event. The key is to iterate over all the assignments for all the tasks, and group the b_ik * z_ie terms by resource, using a Dictionary. Since most resources only work on a small percentage of the tasks, there will be very few nonzero b_ik terms. If we had stored the b_ik in a sparse matrix, this constraint would have been easier to write, on the other hand this would be more inconvenient in other parts of the system.

// Resource usage during each event must not exceed resource capacity. Dictionary<Resource, int> resToId = new Dictionary<Resource, int>(project.Resources.Count); for (int k = 0; k < project.Resources.Count; k++) { resToId[project.Resources[k]] = k; } SumTermBuilder[] totalResWork = new SumTermBuilder[project.Resources.Count]; int c43 = 0; for (int e = 0; e < eventCount; e++) { for (int taskID = 0; taskID < project.Tasks.Count; taskID++) { foreach (var assn in project.Tasks[taskID].Assignments) { int resID = resToId[assn.Resource]; if (totalResWork[resID] == null) { totalResWork[resID] = new SumTermBuilder(4); // most tasks will have <= 4 assignments. } totalResWork[resID].Add(assn.Units * isActive[taskID, e]); } } for (int resID = 0; resID < totalResWork.Length; resID++) { if (totalResWork[resID] != null) { model.AddConstraint("c43_" + c43++, totalResWork[resID].ToTerm() <= project.Resources[resID].MaxUnits); totalResWork[resID] = null; // note: memory churn... } } }

Whew. The rest is easy. We will expose a wrapper for SolverContext.Solve(). The return type for this function is a dictionary that maps task IDs to scheduled start date. We use a LINQ query to produce this dictionary from the Solver Foundation results.

/// <summary>Solve the model. /// </summary> public IDictionary<int, double> Solve() { GurobiDirective directive = new GurobiDirective(); directive.OutputFlag = verbose; //SimplexDirective directive = new SimplexDirective(); // use me if you do not have Gurobi. var solution = context.Solve(directive); if (verbose) { Console.WriteLine(solution.GetReport()); } return GetStartTimesByTask(); } private IDictionary<int, double> GetStartTimesByTask() { // The "start" decision has the start times for each *event*. Index 0 has the start and index 1 has the event. // Recall that the event set was created using Rationals, so we must convert to int here. var eventToStart = start.GetValues().ToDictionary(o => (int)((Rational)o[1]).ToDouble(), o => (double)o[0]); // Group the isActive decision by task. var isActiveByTask = isActive.GetValues().GroupBy(o => o[1]); // Find the first event where each task is active, and transform the output into (task, event) pairs. var firstActiveByTask = isActiveByTask .Select(g => g.First(a => (double)a[0] > 0.5)) .Select(o => new { TaskID = (int)((double)o[1]), EventID = (int)((Rational)o[2]).ToDouble() }); // Now we can use the (task, event) pairs to get a dictionary that maps task IDs to start dates. return firstActiveByTask.ToDictionary(f => f.TaskID, f => eventToStart[f.EventID]); } }

That’s a lot of code and it’s possible (perhaps probable) that I have made a mistake. So in the next post I will write a simple test program to see how it works.