// (c) 2023-2024 Fair Isaac Corporation
using System;
using System.Linq;
using System.Collections.Generic;
using Optimizer.Maps;
using Optimizer.Objects;
using ColumnType = Optimizer.ColumnType;
using static Optimizer.Objects.Utils;
using static Optimizer.XPRSprob;
using static Optimizer.Objects.ConstantExpression;
namespace XpressExamples
{
/// <summary>A production planning example.</summary>
/// <remarks>
/// <p>
/// There is a set F of factories at which products can be produced from
/// raw materials.
/// </p>
/// <p>
/// Products are sold (at a certain price) out of the factories.
/// Products can be produced (at a certain cost) at the factories from raw
/// materials. Each product has raw material requirements.
/// Raw materials can be bought at the factories.
/// Products and raw materials can be stocked.
/// There is an initial stock of products and raw materials.
/// A factory can produce only if it is open. It can buy raw material or
/// sell products even it it is closed.
/// </p>
///
/// <p>
/// <b>Objective</b>
/// Maximize the total profit which comprises of
/// <list type='bullet'>
/// <item><description>
/// the revenue from selling products,
/// </description></item>
/// <item><description>
/// the production cost,
/// </description></item>
/// <item><description>
/// the cost for keeping products on stock,
/// </description></item>
/// <item><description>
/// the cost for keeping raw materials on stock,
/// </description></item>
/// <item><description>
/// the cost for buying raw material,
/// </description></item>
/// <item><description>
/// the cost for keeping a factory open
/// </description></item>
/// </list>
/// </p>
/// </remarks>
public class ProductionPlanning_Index
{
const int PMAX = 2; // number of products
const int FMAX = 2; // number of factories
const int RMAX = 2; // number of raw material
const int TMAX = 4; // number of time periods
const double CPSTOCK = 2.0; // unit cost to store a product
const double CRSTOCK = 2.0; // unit cost to store a raw material
const double MAXRSTOCK = 300; // maximum number of raw material that can be stocked in a factory
static readonly double[,] REV = // REV[p, t] equals unit price for selling product p in period t
new double[PMAX, TMAX] {
{400, 380, 405, 350},
{410, 397, 412, 397}};
static readonly double[,] CMAK = // CMAK[p, f] unit cost for producing product p at factory f
new double[PMAX, FMAX] {
{150, 153},
{ 75, 68}};
static readonly double[,] CBUY = // CBUY[r, t] unit cost to buy raw material r in period t
new double[RMAX, TMAX] {
{100, 98, 97, 100},
{200, 195, 198, 200}};
static readonly double[] COPEN = // COPEN[f] fixed cost for factory f being open for one period
new double[FMAX] { 50000, 63000 };
static readonly double[,] REQ = // REQ[p, r] raw material requirement (in units) of r to make one unit of p
new double[PMAX, RMAX] {
{1.0, 0.5},
{1.3, 0.4}};
static readonly double[,] MAXSELL = // MAXSELL[p, t] maximum number of units that can be sold of product p in period t
new double[PMAX, TMAX] {
{650, 600, 500, 400},
{600, 500, 300, 250}};
static readonly double[] MAXMAKE = // MAXMAKE[f] maximum number of units (over all products) a factory can produce per period
new double[FMAX] { 400, 500 };
static readonly double[,] PSTOCK0 = // PSTOCK0[p, f] initial stock of product p at factory f
new double[PMAX, FMAX] {
{50, 100},
{50, 50}};
static readonly double[,] RSTOCK0 = // RSTOCK0[r, f] initial stock of raw material r at factor f
new double[RMAX, FMAX] {
{100, 150},
{ 50, 100}};
public static void Main(string[] args)
{
Console.WriteLine("Formulating the production planning problem");
using (XpressProblem prob = new XpressProblem())
{
// make[p, f, t]: Amount of product p to make at factory f in period t
Variable[,,] make = prob.AddVariables(PMAX, FMAX, TMAX)
.WithName("make_p{0}_f{1}_t{2}")
.ToArray();
// sell[p, f, t]: Amount of product p to sell from factory f in period t
Variable[,,] sell = prob.AddVariables(PMAX, FMAX, TMAX)
.WithName("sell_p{0}_f{1}_t{2}")
.ToArray();
// pstock[p, f, t]: Stock level of product p at factor f at start of period t
Variable[,,] pstock = prob.AddVariables(PMAX, FMAX, TMAX)
.WithName("pstock_p{0}_f{1}_t{2}")
.ToArray();
// buy[r, f, t]: Amount of raw material r bought for factory f in period t
Variable[,,] buy = prob.AddVariables(RMAX, FMAX, TMAX)
.WithName("buy_r{0}_f{1}_t{2}")
.ToArray();
// rstock[r, f, t]: Stock level of raw material r at factory f at start of period t
Variable[,,] rstock = prob.AddVariables(RMAX, FMAX, TMAX)
.WithName("rstock_r{0}_f{1}_t{2}")
.ToArray();
// openm[f, t]: If factory f is open in period t
Variable[,] openm = prob.AddVariables(FMAX, TMAX)
.WithType(ColumnType.Binary)
.WithUB(1)
.WithName("openm_f{0}_t{1}")
.ToArray();
// ## Objective:
// Maximize total profit
// revenue from selling products
// + REV[p, t] * sell[p, f, t]
LinExpression revenue = LinExpression.Create();
foreach (Index idx in LoopThrough3D(PMAX, FMAX, TMAX))
{
int p = idx.i1, f = idx.i2, t = idx.i3;
revenue.AddTerm(REV[p, t] * sell[p, f, t]);
}
// cost for making products (must be subtracted from profit)
// - CMAK[p, f] * make[p, f, t]
LinExpression prodCost = LinExpression.Create();
foreach (Index idx in LoopThrough3D(PMAX, FMAX, TMAX))
{
int p = idx.i1, f = idx.i2, t = idx.i3;
prodCost.AddTerm(-CMAK[p, f] * make[p, f, t]);
}
// cost for storing products (must be subtracted from profit)
// - CPSTOCK * pstock[p, f, t]
LinExpression prodStorageCost = LinExpression.Create();
foreach (Index idx in LoopThrough3D(PMAX, FMAX, TMAX))
{
int p = idx.i1, f = idx.i2, t = idx.i3;
prodStorageCost.AddTerm(-CPSTOCK * pstock[p, f, t]);
}
// cost for opening a factory in a time period
// - openm[f, t] * COPEN[f]
LinExpression factoryCost = LinExpression.Create();
foreach (Index idx in LoopThrough2D(FMAX, TMAX))
{
int f = idx.i1, t = idx.i2;
factoryCost.AddTerm(-COPEN[f] * openm[f, t]);
}
// cost for buying raw material in time period t
// - buy[r, f, t] * CBUY[r, t]
LinExpression rawMaterialBuyCost = LinExpression.Create();
foreach (Index idx in LoopThrough3D(RMAX, FMAX, TMAX))
{
int r = idx.i1, f = idx.i2, t = idx.i3;
rawMaterialBuyCost.AddTerm(-CBUY[r, t] * buy[r, f, t]);
}
// cost for storing raw material (must be subtracted from profit)
// - rstock[r, f, t] * CRSTOCK
// an alternative way of setting an objective Expression is the below nested sum expression
// this construction does not use an index loop
Expression rawMaterialStorageCost = Sum(FMAX,
f => Sum(RMAX,
r => Sum(TMAX, t => rstock[r, f, t] * -CRSTOCK)
)
);
// sum up the 6 individual contributions to the overall profit
Expression profit = Sum(revenue, prodCost, prodStorageCost, factoryCost, rawMaterialStorageCost, rawMaterialBuyCost);
// set maximization of profit as objective function
prob.SetObjective(profit, Optimizer.ObjSense.Maximize);
// constraints
// Product stock balance
foreach (Index idx in LoopThrough3D(PMAX, FMAX, TMAX))
{
int p = idx.i1, f = idx.i2, t = idx.i3;
// for each time period except the last time period, surplus is available as stock for the next time period
if (t < TMAX - 1)
{
prob.AddConstraint(
pstock[p, f, t] + make[p, f, t] == sell[p, f, t] + pstock[p, f, t + 1])
.SetName($"prod_stock_balance_p{p}_f{f}_t{t}");
}
else
{
prob.AddConstraint(
pstock[p, f, t] + make[p, f, t] >= sell[p, f, t])
.SetName($"prod_stock_balance_p{p}_f{f}_t{t}");
}
}
// Raw material stock balance
foreach (Index idx in LoopThrough3D(RMAX, FMAX, TMAX))
{
int r = idx.i1, f = idx.i2, t = idx.i3;
if (t < TMAX - 1)
{
prob.AddConstraint(
rstock[r, f, t] + buy[r, f, t] == rstock[r, f, t + 1] + Sum(PMAX, p => make[p, f, t] * REQ[p, r]))
.SetName($"raw_material_stock_balance_r{r}_f{f}_t{t}");
}
else
{
prob.AddConstraint(
rstock[r, f, t] + buy[r, f, t] >= Sum(PMAX, p => make[p, f, t] * REQ[p, r]))
.SetName($"raw_material_stock_balance_r{r}_f{f}_t{t}");
}
}
// Limit on the amount of product p to be sold
// exemplifies how to loop through multiple ranges
prob.AddConstraints(PMAX, TMAX,
(p, t) =>
(Sum(FMAX, f => sell[p, f, t]) <= MAXSELL[p, t])
.SetName($"maxsell_p{p}_t{t}")
);
// Capacity limit at factory f
// exemplifies how to loop through multiple ranges
prob.AddConstraints(FMAX, TMAX,
(f, t) =>
(Sum(PMAX, p => make[p, f, t]) <= openm[f, t].Mul(MAXMAKE[f]))
.SetName($"capacity_f{f}_t{t}")
);
// Raw material stock limit
// exemplifies how to loop through multiple ranges
prob.AddConstraints(FMAX, TMAX,
(f, t) =>
(Sum(RMAX, r => rstock[r, f, t]) <= MAXRSTOCK)
.SetName($"raw_material_stock_limit_f{f}_t{t}")
);
// Initial product storage
// loop through a custom IEnumerable, in our case our custom index object.
prob.AddConstraints(LoopThrough2D(PMAX, FMAX),
idx =>
// pstock is indexed (p, f, t), PSTOCK0 is indexed (p, f)
(pstock[idx.i1, idx.i2, 0] == PSTOCK0[idx.i1, idx.i2])
.SetName($"initial_product_stock_p{idx.i1}_f{idx.i2}")
);
// Initial raw material storage
// classic for loop
foreach (Index idx in LoopThrough2D(RMAX, FMAX))
{
int r = idx.i1, f = idx.i2;
prob.AddConstraint((rstock[r, f, 0] == RSTOCK0[r, f])
.SetName($"initial_raw_material_stock_r{r}_f{f}"));
}
// write the problem in LP format for manual inspection
Console.WriteLine("Writing the problem to 'ProductionPlanning.lp'");
prob.WriteProb("ProductionPlanning.lp", "l");
// Solve the problem
Console.WriteLine("Solving the problem");
prob.Optimize();
Console.WriteLine("Problem finished with SolStatus {0}", prob.SolStatus);
if (prob.SolStatus != Optimizer.SolStatus.Optimal)
{
throw new Exception("Problem not solved to optimality");
}
Console.WriteLine("Solution has objective value (profit) of {0}", prob.MIPObjVal);
Console.WriteLine("");
Console.WriteLine("*** Solution ***");
double[] sol = prob.GetSolution();
// Is factory f open at time period t?
for (int f = 0; f < FMAX; f++)
{
for (int t = 0; t < TMAX; t++)
{
Console.WriteLine($"Factory {f}\tTime Period {t} : Open = {openm[f, t].GetValue(sol)}");
}
}
Console.WriteLine("");
// Production plan for producing
WriteSol3D(sol, make, PMAX, FMAX, TMAX, new string[] { "Product", "Factory", "Time Period" }, "Make");
Console.WriteLine("");
// Production plan for selling products
WriteSol3D(sol, sell, PMAX, FMAX, TMAX, new string[] { "Product", "Factory", "Time Period" }, "Sell");
Console.WriteLine("");
// Production plan for keeping products in stock
WriteSol3D(sol, pstock, PMAX, FMAX, TMAX, new string[] { "Product", "Factory", "Time Period" }, "Pstock");
Console.WriteLine("");
// Production plan for keeping raw material in stock
WriteSol3D(sol, rstock, RMAX, FMAX, TMAX, new string[] { "Material", "Factory", "Time Period" }, "Rstock");
Console.WriteLine("");
// Buying plan for raw material
WriteSol3D(sol, buy, RMAX, FMAX, TMAX, new string[] { "Material", "Factory", "Time Period" }, "Buy");
Console.WriteLine("");
}
}
/// <summary>Helper class for convenient index iteration</summary>
public sealed class Index
{
public readonly int i1;
public readonly int i2;
public readonly int i3;
public Index(int i1, int i2, int i3)
{
this.i1 = i1;
this.i2 = i2;
this.i3 = i3;
}
public override string ToString() { return $"Index {this.i1} {this.i2} {this.i3}"; }
}
/// <summary>
/// Loop through a 3-dimensional range, always starting at 0.
/// </summary>
/// <param name="max1">First index varies between 0 and max1</param>
/// <param name="max2">Second index varies between 0 and max2</param>
/// <param name="max3">Third index varies between 0 and max3</param>
/// <returns>an Enumerable that can be used to loop over the cross product of the three ranges</returns>
public static IEnumerable<Index> LoopThrough3D(int max1, int max2, int max3)
{
for (int i1 = 0; i1 < max1; i1++)
{
for (int i2 = 0; i2 < max2; i2++)
{
if (max3 > 0)
{
for (int i3 = 0; i3 < max3; i3++)
{
Index idx = new Index(i1, i2, i3);
yield return idx;
}
}
else
{
yield return new Index(i1, i2, -1);
}
}
}
}
/// <summary>
/// Loop through a 2-dimensional range, starting at 0
/// </summary>
/// <param name="max1">First index varies between 0 and max1</param>
/// <param name="max2">Second index varies between 0 and max2</param>
/// <returns>an Enumerable that can be used to loop over the cross product of the two ranges. The third
/// entry of the returned Index object will always be -1</returns>
public static IEnumerable<Index> LoopThrough2D(int max1, int max2)
{
return LoopThrough3D(max1, max2, 0);
}
/// <summary>
/// Convenience function for printing solution values stored in a 3-dimensional Variable array
/// </summary>
/// <param name="sol">solution object, obtained via prob.getSolution()</param>
/// <param name="array">3-dimensional array of Xpress Variables</param>
/// <param name="max1">First index varies between 0 and max1</param>
/// <param name="max2">Second index varies between 0 and max2</param>
/// <param name="max3">Third index varies between 0 and max3</param>
/// <param name="dimNames">An array with a name for every dimension</param>
/// <param name="name">The name of the array</param>
public static void WriteSol3D(double[] sol, Variable[,,] array, int max1, int max2, int max3, string[] dimNames, string name)
{
foreach (Index idx in LoopThrough3D(max1, max2, max3))
{
int i1 = idx.i1, i2 = idx.i2, i3 = idx.i3;
Console.WriteLine($"{dimNames[0]} {i1}\t{dimNames[1]} {i2}\t{dimNames[2]} {i3} : {name} = {array[i1, i2, i3].GetValue(sol)}");
}
}
}
}
|
