Originally posted by: Sulivan

I am formulating a resource constrainted scheduling problem: you have a series of tasks to be completed, which consume resource. The task operator(s) have limited capacity of resource, which can be smaller or larger than the requirement of resource of a single task, starting with full resource. Whenever the resource is depleted or not enough to complete any task, the task operator has to replenish the resource before it can operate again. For some tasks whose resource consumption are low, the task operator can complete serveral of them between two replenishments; for some tasks whose resource consumption are high (larger than the capacity), the task operator has to replenish at least once to complete them.

Here is a simplified version of my model with only 1 task operator and 5 tasks:

using CP;
 
{int} Tasks={20,20,40,130,80};
int Capacity=50;
int RepTime=5;                  //Replenish time

dvar interval opr[t in Tasks];  //Operation of tasks
dvar int+ rep[t in Tasks];      //Number of replenishments before/during operating the task
dvar sequence seq in all(t in Tasks) opr[t];
cumulfunction resource = step (0,Capacity) + sum(t in Tasks) stepAtStart(opr[t],0,rep[t]*RepTime) - sum(t in Tasks) stepAtStart(opr[t],t);

minimize max(t in Tasks) endOf(opr[t]);
subject to{
ct1: resource >= 0;
ct2: forall(t in Tasks) sizeOf(opr[t]) == t+rep[t]*RepTime;
ct3: forall(t in Tasks) noOverlap(seq);
}

I get errors in the cumulative function. Seems that rep[t] cannot be used as the parameter, since it is a decision variable. How can I deal with it?

Besides, I'd like to ask the details about difference between stepAtStart(a, h) and stepAtStart(a, hmin, hmax). For example, consider the resource value is 10 before interval a. If I use stepAtStart(a, 50), it well burst to 60 (increase by 50) at the start of a; if I use stepAtStart(a, 0, 50), it will only burst to 50 (increase up to 50). Is that correct?

Originally posted by: PhilippeLaborie

Hello,
I think you almost got it right. The only thing is that the basic cumul functions like startAtStart do not take a decision variable (dvar int) as parameter but only a range.
So you can model resource production by a stepAtStart(prod,0,MaxProduction): it will increase the level of the resource at the start time of interval variable 'prod' by a quantity in range [0,MaxProduction] fixed by the engine. You can get this quantity as an integer expression by using heightAtStart(prod,resource) that returns the contribution of interval variable 'prod' to cumul function 'resource' at its start time.
Note that it is possible to say as in your model that it is the same activity 'opr' that consumes the resource and potentially produces it, but in this case, heightAtStart(opr,resource) will return the sum of the two contribution (consumption and production). Here is a version of your model with additional production on the operator, one for each 'opr' activity. In my model I denote Qty[t] the quantity of resource consumed by operation t (not sure why it was equal to t in your model); beware that in your model, Tasks is a set, so you will only have one tasks of consumption/duration 20.

using CP;
range Tasks = 1..5;
int Qty[t in Tasks] = [20,20,40,130,80];
int Capacity = 50;
int RepTime  = 5;   // Replenish time
dvar interval opr[t in Tasks];
dvar sequence seq in opr;
cumulfunction resource = 
    step (0,Capacity)                              // Initial level
  + sum(t in Tasks) stepAtStart(opr[t], 0, 1000)   // Production
  - sum(t in Tasks) stepAtStart(opr[t], Qty[t]);   // Consumption
dexpr int rep[t in Tasks] = heightAtStart(opr[t],resource) + Qty[t];      
minimize max(t in Tasks) endOf(opr[t]);
subject to {
  ct1: resource <= Capacity;
  forall(t in Tasks) {
   ct2: sizeOf(opr[t]) == Qty[t] + rep[t]*RepTime;
   ct3: noOverlap(seq);
  }
}

Originally posted by: Sulivan

Dear Philippe, thank you very much for your answer! In fact I can define rep[t] as intervals whose length equal to the total length or replenish time. Then rep[t] would be the production of resource. However the problem is, the quatity of resource produced each time in not determined but depend on the current quantity, as well as the quantity required by the next task.

Originally posted by: PhilippeLaborie

If you do nothing special on this produced resource quantity, the engine will choose any value that satisfies the constraints and optimizes the objective function. If in the model you state the constraint that the resource capacity cannot exceed Capacity (and is non-negative), then the engine will ensure that activities rep[t] produce enough for the consumption activities (resource level >= 0) but not too much so that the resource never overflows (resource level <= capacity). And the objective function which is to finish as soon as possible will ensure that you do not produce fr nothing (as producing increasing the duration of the producing activities).

If you need to post additional constraints on the quantity produced that depends on the actual level at the time of the activities, you could use an alternative model that do not use a cumul function but tracks the current resource level before/after all operations in the sequence, as suggested in this thread: https://www.ibm.com/developerworks/community/forums/html/topic?id=5c799735-2717-4456-ba41-80a08730a4c6&ps=25

Philippe

Originally posted by: Sulivan

Thank you for your suggestion. I have been attempting using the method your suggested in the thread. This time I've made a 10-task 2-operator example. Same task assigned to multiple opertors is allowed, as long as the sum of operation time reaches the required.

using CP;
range Tasks=1..10;
range Operators=1..2;
int Qty[t in Tasks]=...;
int Capacity[o in Operators]=...;
int RepTime[o in Operators]=...;

dvar interval opr[o in Operators][t in Tasks] optional;
dvar sequence seq[o in Operators] in all(t in Tasks) opr[o][t] types all(t in Tasks) t;
dvar int+ rep[o in Operators][t in Tasks];
dvar int+ resSt[o in Operators][t in Tasks];
dvar int+ resEd[o in Operators][t in Tasks];

minimize max(o in Operators, t in Tasks) endOf(opr[o][t]);
subject to{
  forall(o in Operators, t in Tasks){
    ct1: presenceOf(opr[o][t]) == 0 => rep[o][t] == 0;
    ct2: presenceOf(opr[o][t]) == 0 => resSt[o][t] == 0;
    ct3: presenceOf(opr[o][t]) == 0 => resEd[o][t] == 0;
    ct4: presenceOf(opr[o][t]) == 1 => resEd[o][t] == resSt[o][t] + rep[o][t]*Capacity[o] - (sizeOf(opr[o][t]) - rep[o][t]*RepTime[o]);   //Ed_resource = St_resource + total replenishment - total consumption.
    ct5: presenceOf(opr[o][t]) == 1 => resSt[o][t] == resEd[o][typeOfPrev(seq[o],opr[o][t],Capacity[o],0)];                               //St_resource = Ed_resource of previous task.
    ct6: noOverlap(seq[o]);
  }
  ct7: forall(t in Tasks) sum(o in Operators) (sizeOf(opr[o][t]) - rep[o][t]*RepTime[o]) == Qty[t];                                       //All tasks are operated for enough duration
}

The data is as follows:

Qty=[20,20,40,130,80,70,30,60,40,100];
Capacity=[50,70];
RepTime=[5,7];

However ct5 itself yields a conflict. Am I improperly using typeOfPrev?

Originally posted by: PhilippeLaborie

Yes indeed, the problem coms from this constraint ct5:

ct5: presenceOf(opr[o][t]) == 1 => resSt[o][t] == resEd[o][typeOfPrev(seq[o],opr[o][t],Capacity[o],0)];                               //St_resource = Ed_resource of previous task.

This constraint says that the level at the start time of a task is equal to the level at the send time of the previous task (which is what you want), but for the first interval in the sequence, the value of the expression will be Capacity[o], but Capacity[o] is a capacity value and it cannot be used to index an array indexed on task identifier. 

You could have a special index (typically 0) for a (non-existing) task that is before the first actual task. For that you just need to index the variables resSt and resEd starting from 0 instead of 1. And then specify a return value of 0 for typeOfPrev(seq[o],opr[o][t]) when the interval is first.

Additionally, I think it would be more efficient to avoid as much as possible the meta-constraints presenceOf(opr[o][t]) =>  XXX and instead leverage the notion of optionality. When you know that the right side of the implication is satisfied when the interval is absent, you do not need to write the left side condition  presenceOf(opr[o][t]) =>.

You can also limit the domain of the variables if you know a reasonable bound. And for instance here I suppose you want that the level of the resource never exceeds Capacity. Same for the horizon of the schedule if you know a reasonable one.

Here is how I would reformulate your model:

using CP;
range Tasks=1..10;
range Operators=1..2;
int Qty[Tasks]=[20,20,40,130,80,70,30,60,40,100];
int Capacity[Operators]=[50,70];
int RepTime[Operators]=[5,7];
int Horizon = 1000;

dvar interval opr[o in Operators][t in Tasks] optional in 0..Horizon size 1..Horizon;
dvar sequence seq[o in Operators] in all(t in Tasks) opr[o][t] types all(t in Tasks) t;
dvar int+ rep[o in Operators][t in Tasks];
dvar int+ resSt[o in Operators][t in 0..10] in 0..Capacity[o];
dvar int+ resEd[o in Operators][t in 0..10] in 0..Capacity[o];

minimize max(o in Operators, t in Tasks) endOf(opr[o][t]);
subject to{
  forall(o in Operators) {
    resSt[o][0] == 0;
    resEd[o][0] == 0;
  }
  forall(o in Operators, t in Tasks){
    ct1: !presenceOf(opr[o][t]) => ( (rep[o][t] == 0) && (resSt[o][t] == 0) && (resEd[o][t] == 0));
    ct4: resEd[o][t] == resSt[o][t] + rep[o][t]*Capacity[o] - (sizeOf(opr[o][t]) - rep[o][t]*RepTime[o]);   
    ct5: resSt[o][t] == resEd[o][typeOfPrev(seq[o],opr[o][t],0,0)];
    ct6: noOverlap(seq[o]);
  }
  forall(t in Tasks) {
    ct7: Qty[t] == sum(o in Operators) (sizeOf(opr[o][t]) - rep[o][t]*RepTime[o]);  //All tasks are operated for enough duration  
  }                                     
}