RMOL Logo  1.00.0
C++ library of Revenue Management and Optimisation classes and functions
 All Classes Namespaces Files Functions Variables Typedefs Friends Pages
MCOptimiser.cpp
Go to the documentation of this file.
1 // //////////////////////////////////////////////////////////////////////
2 // Import section
3 // //////////////////////////////////////////////////////////////////////
4 // STL
5 #include <cassert>
6 #include <string>
7 #include <sstream>
8 #include <algorithm>
9 #include <cmath>
10 // StdAir
11 #include <stdair/stdair_basic_types.hpp>
12 #include <stdair/bom/BomManager.hpp>
13 #include <stdair/bom/LegCabin.hpp>
14 #include <stdair/bom/SegmentCabin.hpp>
15 #include <stdair/bom/BookingClass.hpp>
16 #include <stdair/bom/VirtualClassStruct.hpp>
17 #include <stdair/service/Logger.hpp>
18 #include <stdair/basic/RandomGeneration.hpp>
19 #include <stdair/basic/BasConst_General.hpp>
20 // RMOL
22 #include <rmol/bom/MCOptimiser.hpp>
23 
24 namespace RMOL {
25 
26  // // //////////////////////////////////////////////////////////////////////
27  void MCOptimiser::
28  optimalOptimisationByMCIntegration (stdair::LegCabin& ioLegCabin) {
29  // Retrieve the segment-cabin
30  const stdair::SegmentCabinList_T lSegmentCabinList =
31  stdair::BomManager::getList<stdair::SegmentCabin> (ioLegCabin);
32  stdair::SegmentCabinList_T::const_iterator itSC = lSegmentCabinList.begin();
33  assert (itSC != lSegmentCabinList.end());
34  const stdair::SegmentCabin* lSegmentCabin_ptr = *itSC;
35  assert (lSegmentCabin_ptr != NULL);
36 
37  // Retrieve the class list.
38  const stdair::BookingClassList_T lBookingClassList =
39  stdair::BomManager::getList<stdair::BookingClass> (*lSegmentCabin_ptr);
40 
41  // Retrieve the remaining cabin capacity.
42  const stdair::Availability_T& lCap = ioLegCabin.getAvailabilityPool();
43  const int lCapacity = static_cast<const int> (lCap);
44  const stdair::UnsignedIndex_T lCapacityIndex =
45  static_cast<const stdair::UnsignedIndex_T> ((lCapacity+abs(lCapacity))/2);
46 
47  // Retrieve the virtual class list.
48  stdair::VirtualClassList_T& lVCList = ioLegCabin.getVirtualClassList();
49 
50  // Parse the virtual class list and compute the protection levels.
51  stdair::VirtualClassList_T::iterator itCurrentVC = lVCList.begin();
52  assert (itCurrentVC != lVCList.end());
53  stdair::VirtualClassList_T::iterator itNextVC = itCurrentVC; ++itNextVC;
54 
55  // Initialise the partial sum holder with the demand sample of the first
56  // virtual class.
57  stdair::VirtualClassStruct& lFirstVC = *itCurrentVC;
58  stdair::GeneratedDemandVector_T lPartialSumHolder =
59  lFirstVC.getGeneratedDemandVector();
60 
61  // Initialise the booking limit for the first class, which is equal to
62  // the remaining capacity.
63  lFirstVC.setCumulatedBookingLimit (lCap);
64 
65  // Initialise bid price vector with the first element (index 0) equal to
66  // the highest yield.
67  ioLegCabin.emptyBidPriceVector();
68  stdair::BidPriceVector_T& lBPV = ioLegCabin.getBidPriceVector();
69  //const stdair::Yield_T& y1 = lFirstVC.getYield ();
70  //lBPV.push_back (y1);
71  stdair::UnsignedIndex_T idx = 1;
72 
73  for (; itNextVC != lVCList.end(); ++itCurrentVC, ++itNextVC) {
74  // Get the yields of the two classes.
75  stdair::VirtualClassStruct& lCurrentVC = *itCurrentVC;
76  stdair::VirtualClassStruct& lNextVC = *itNextVC;
77  const stdair::Yield_T& yj = lCurrentVC.getYield ();
78  const stdair::Yield_T& yj1 = lNextVC.getYield ();
79 
80  // Consistency check: the yield/price of a higher class/bucket
81  // (with the j index lower) must be higher.
82  assert (yj > yj1);
83 
84  // Sort the partial sum holder.
85  std::sort (lPartialSumHolder.begin(), lPartialSumHolder.end());
86  const stdair::UnsignedIndex_T K = lPartialSumHolder.size ();
87 
88  // Compute the optimal index lj = floor {[y(j)-y(j+1)]/y(j) . K}
89  const double ljdouble = std::floor (K * (yj - yj1) / yj);
90  stdair::UnsignedIndex_T lj =
91  static_cast<stdair::UnsignedIndex_T> (ljdouble);
92 
93  // Consistency check.
94  assert (lj >= 1 && lj < K);
95 
96  // The optimal protection: p(j) = 1/2 [S(j,lj) + S(j, lj+1)]
97  const double sjl = lPartialSumHolder.at (lj - 1);
98  const double sjlp1 = lPartialSumHolder.at (lj + 1 - 1);
99  const double pj = (sjl + sjlp1) / 2;
100 
101  // Set the cumulated protection level for the current class.
102  lCurrentVC.setCumulatedProtection (pj);
103  // Set the cumulated booking limit for the next class.
104  lNextVC.setCumulatedBookingLimit (lCap - pj);
105 
110  const stdair::UnsignedIndex_T pjint = static_cast<const int> (pj);
111  stdair::GeneratedDemandVector_T::iterator itLowerBound =
112  lPartialSumHolder.begin();
113  for (; idx <= pjint && idx <= lCapacityIndex; ++idx) {
114  itLowerBound =
115  std::lower_bound (itLowerBound, lPartialSumHolder.end(), idx);
116  const stdair::UnsignedIndex_T pos =
117  itLowerBound - lPartialSumHolder.begin();
118 
119  const stdair::BidPrice_T lBP = yj * (K - pos) / K;
120  lBPV.push_back (lBP);
121  }
122 
123  // Update the partial sum holder.
124  const stdair::GeneratedDemandVector_T& lNextPSH =
125  lNextVC.getGeneratedDemandVector();
126  assert (K <= lNextPSH.size());
127  for (stdair::UnsignedIndex_T i = 0; i < K - lj; ++i) {
128  lPartialSumHolder.at(i) = lPartialSumHolder.at(i + lj) + lNextPSH.at(i);
129  }
130  lPartialSumHolder.resize (K - lj);
131  }
132 
137  stdair::VirtualClassStruct& lLastVC = *itCurrentVC;
138  const stdair::Yield_T& yn = lLastVC.getYield();
139  stdair::GeneratedDemandVector_T::iterator itLowerBound =
140  lPartialSumHolder.begin();
141  for (; idx <= lCapacityIndex; ++idx) {
142  itLowerBound =
143  std::lower_bound (itLowerBound, lPartialSumHolder.end(), idx);
144  const stdair::UnsignedIndex_T pos =
145  itLowerBound - lPartialSumHolder.begin();
146  const stdair::UnsignedIndex_T K = lPartialSumHolder.size();
147  const stdair::BidPrice_T lBP = yn * (K - pos) / K;
148  lBPV.push_back (lBP);
149  }
150  }
151 
152  // ///////////////////////////////////////////////////////////////////
153  stdair::GeneratedDemandVector_T MCOptimiser::
154  generateDemandVector (const stdair::MeanValue_T& iMean,
155  const stdair::StdDevValue_T& iStdDev,
156  const stdair::NbOfSamples_T& K) {
157  stdair::GeneratedDemandVector_T oDemandVector;
158  if (iStdDev > 0) {
159  stdair::RandomGeneration lGenerator (stdair::DEFAULT_RANDOM_SEED);
160  for (unsigned int i = 0; i < K; ++i) {
161  stdair::RealNumber_T lDemandSample =
162  lGenerator.generateNormal (iMean, iStdDev);
163  oDemandVector.push_back (lDemandSample);
164  }
165  } else {
166  for (unsigned int i = 0; i < K; ++i) {
167  oDemandVector.push_back (iMean);
168  }
169  }
170  return oDemandVector;
171  }
172 
173  // /////////////////////////////////////////////////////////////////////////
174  void MCOptimiser::
175  optimisationByMCIntegration (stdair::LegCabin& ioLegCabin) {
176  // Number of MC samples
177  stdair::NbOfSamples_T K = DEFAULT_NUMBER_OF_DRAWS_FOR_MC_SIMULATION;
178 
179  const stdair::YieldLevelDemandMap_T& lYieldDemandMap =
180  ioLegCabin.getYieldLevelDemandMap();
181  assert (!lYieldDemandMap.empty());
182 
183  std::ostringstream oStr;
184  oStr << "Yield list ";
185  for (stdair::YieldLevelDemandMap_T::const_iterator itYD =
186  lYieldDemandMap.begin();
187  itYD != lYieldDemandMap.end(); ++itYD) {
188  const stdair::Yield_T& y = itYD->first;
189  oStr << y << " ";
190  }
191 
192  STDAIR_LOG_DEBUG (oStr.str());
193  ioLegCabin.emptyBidPriceVector();
194  stdair::BidPriceVector_T& lBidPriceVector =
195  ioLegCabin.getBidPriceVector();
196  const stdair::Availability_T& lAvailabilityPool =
197  ioLegCabin.getAvailabilityPool();
198  // Initialise the minimal bid price to 1.0 (just to avoid problems
199  // of division by zero).
200  const stdair::BidPrice_T& lMinBP = 1.0;
201 
202  stdair::YieldLevelDemandMap_T::const_reverse_iterator itCurrentYD =
203  lYieldDemandMap.rbegin();
204  stdair::YieldLevelDemandMap_T::const_reverse_iterator itNextYD = itCurrentYD;
205  ++itNextYD;
206 
207  // Initialise the partial sum holder
208  stdair::MeanStdDevPair_T lMeanStdDevPair = itCurrentYD->second;
209  stdair::GeneratedDemandVector_T lPartialSumHolder =
210  generateDemandVector(lMeanStdDevPair.first, lMeanStdDevPair.second, K);
211 
212  stdair::UnsignedIndex_T idx = 1;
213  for (; itNextYD!=lYieldDemandMap.rend(); ++itCurrentYD, ++itNextYD) {
214  const stdair::Yield_T& yj = itCurrentYD->first;
215  const stdair::Yield_T& yj1 = itNextYD->first;
216  // Consistency check: the yield/price of a higher class/bucket
217  // (with the j index lower) must be higher.
218  assert (yj > yj1);
219  // Sort the partial sum holder.
220  std::sort (lPartialSumHolder.begin(), lPartialSumHolder.end());
221  // STDAIR_LOG_DEBUG ("Partial sums : max = " << lPartialSumHolder.back()
222  // << " min = " << lPartialSumHolder.front());
223  K = lPartialSumHolder.size ();
224  // Compute the optimal index lj = floor {[y(j)-y(j+1)]/y(j) . K}
225  const double ljdouble = std::floor (K * (yj - yj1) / yj);
226  stdair::UnsignedIndex_T lj =
227  static_cast<stdair::UnsignedIndex_T> (ljdouble);
228  // Consistency check.
229  assert (lj >= 1 && lj < K);
230  // The optimal protection: p(j) = 1/2 [S(j,lj) + S(j, lj+1)]
231  const double sjl = lPartialSumHolder.at (lj - 1);
232  const double sjlp1 = lPartialSumHolder.at (lj + 1 - 1);
233  const double pj = (sjl + sjlp1) / 2;
238  const stdair::UnsignedIndex_T pjint = static_cast<const int> (pj);
239  stdair::GeneratedDemandVector_T::iterator itLowerBound =
240  lPartialSumHolder.begin();
241  for (; idx <= pjint && idx <= lAvailabilityPool; ++idx) {
242  itLowerBound =
243  std::lower_bound (itLowerBound, lPartialSumHolder.end(), idx);
244  const stdair::UnsignedIndex_T pos =
245  itLowerBound - lPartialSumHolder.begin();
246 
247  const stdair::BidPrice_T lBP = yj * (K - pos) / K;
248  lBidPriceVector.push_back (lBP);
249  }
250  // Update the partial sum holder.
251  lMeanStdDevPair = itNextYD->second;
252  const stdair::GeneratedDemandVector_T& lNextDV =
253  generateDemandVector (lMeanStdDevPair.first,
254  lMeanStdDevPair.second, K - lj);
255  for (stdair::UnsignedIndex_T i = 0; i < K - lj; ++i) {
256  lPartialSumHolder.at(i) = lPartialSumHolder.at(i + lj) + lNextDV.at(i);
257  }
258  lPartialSumHolder.resize (K - lj);
259  }
268  // STDAIR_LOG_DEBUG ("Partial sums : max = " << lPartialSumHolder.back()
269  // << " min = " << lPartialSumHolder.front());
270 
271  std::sort (lPartialSumHolder.begin(), lPartialSumHolder.end());
272  const stdair::Yield_T& yn = itCurrentYD->first;
273  stdair::GeneratedDemandVector_T::iterator itLowerBound =
274  lPartialSumHolder.begin();
275  K = lPartialSumHolder.size();
276 
277  bool lMinBPReached = false;
278  for (; idx <= lAvailabilityPool; ++idx) {
279  itLowerBound =
280  std::lower_bound (itLowerBound, lPartialSumHolder.end(), idx);
281 
282  if (!lMinBPReached) {
283  const stdair::UnsignedIndex_T pos =
284  itLowerBound - lPartialSumHolder.begin();
285  stdair::BidPrice_T lBP = yn * (K - pos) / K;
286 
287  if (lBP < lMinBP) {
288  lBP = lMinBP; lMinBPReached = true;
289  }
290 
291  lBidPriceVector.push_back (lBP);
292 
293  } else {
294  lBidPriceVector.push_back (lMinBP);
295  }
296 
297  }
298 
299  // Updating the bid price values
300  ioLegCabin.updatePreviousBidPrice();
301  ioLegCabin.setCurrentBidPrice (lBidPriceVector.back());
302 
303  // Compute and display the bid price variation after optimisation
304  const stdair::BidPrice_T lPreviousBP = ioLegCabin.getPreviousBidPrice();
305  stdair::BidPrice_T lNewBP = ioLegCabin.getCurrentBidPrice();
306  // Check
307  assert (lPreviousBP != 0);
308  stdair::BidPrice_T lBidPriceDelta = lNewBP - lPreviousBP;
309 
310  double lBidPriceVariation = 100*lBidPriceDelta/lPreviousBP;
311 
312  STDAIR_LOG_DEBUG ("Bid price: previous value " << lPreviousBP
313  << ", new value " << lNewBP
314  << ", variation " << lBidPriceVariation << " %"
315  <<", BPV size " << lBidPriceVector.size());
316  }
317 
318 }