Wait long, buy low

Let $\Phi(x)$ denote the area under the standard normal distribution in the range $[-\infty,x]$ . Assume that we set a threshold of $x$ for buying the gas. Hence our expected cost would be $C(x)=x\bigg(1-\big(1-\Phi\left(\frac{x-\bar{g}}{\sigma_G}\right)\big)^{29}\bigg)+ \bar{g}\bigg(1-\Phi\left(\frac{x-\bar{g}}{\sigma_G}\right)\bigg)^{29}$ Here the first term gives the expected cost when the gas price falls below $x$ in the first $29$ days and the second term gives the expected cost when it does not. Now the answer may be obtained by optimizing $C(x)$ over $x$ graphically.

Here is one way to get to the answer without programming.

Please excuse me for posting an image since I am not good enough with LaTeX yet.

p.s. One way to turn this into a true programming problem would be to buy gas on such a day that it is unlikely to witness a lower price on the "remaining" days, which means we will accept higher prices as we get nearer to the end of the month. Also, this strategy seems to be better, judging from the outputs of my program centered around 2.561.

If someone would like to run some massive simulation, I am attaching my code.

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171

#include <iostream>
#include <random>
#include <chrono>
#include <thread>
#include <numeric>

template <typename T>
class FuelPriceSimulator
{
    private:
        static unsigned seedNum;
        mutable std::default_random_engine random_generator;
        mutable std::normal_distribution<T> price_distribution;

        unsigned maxDayNum;
        unsigned maxIter;
    public:
        FuelPriceSimulator(const T& meanPrice, const T& priceSd, unsigned maxDayNum, unsigned maxIter) :
        random_generator(std::chrono::system_clock::now().time_since_epoch().count()+FuelPriceSimulator::seedNum),
        price_distribution(meanPrice, priceSd),
        maxDayNum(maxDayNum),
        maxIter(maxIter)
        {
            ++FuelPriceSimulator::seedNum;
        }

        T GetPrice() const
        {
            return this->price_distribution(this->random_generator);
        }

        T SingleRun() const
        {
            unsigned curDay = 1;
            T curPrice = this->GetPrice();
            while(!shouldBuy(curDay, curPrice))
            {
                curPrice = this->GetPrice();
                ++curDay;
            }
            return curPrice;
        }

    private:
        bool shouldBuy(const unsigned dayNumber, const T& inPrice) const
        {
            bool res;
            if(dayNumber >= this->maxDayNum)
            {
                res = true;
            }
            else
            {
                unsigned withLower = 0;
                for(unsigned iSim = 0; iSim < this->maxIter; ++iSim)
                {
                    if(this->simulateNextDays(inPrice))
                    {
                        ++withLower;
                    }
                }
                if(T(withLower)/T(this->maxIter) > 0.5)
                {
                    res = false;
                }
                else
                {
                    res = true;
                }
            }
            return res;
        }

        bool simulateNextDays(const T& inPrice) const
        {
            bool res = false;
            for(unsigned curDayNum = 0; curDayNum < this->maxDayNum; ++curDayNum)
            {
                if(inPrice > this->GetPrice())
                {
                    res = true;
                    break;
                }
            }
            return res;
        }
};

template<typename T>
unsigned FuelPriceSimulator<T>::seedNum = 0;

template <typename T>
T mcFuelPriceSingle(const T& meanPrice, const T& priceSd, const unsigned maxDayNum, const unsigned maxIter, const unsigned totalIter)
{
    FuelPriceSimulator<double> fpg(meanPrice, priceSd, maxDayNum, maxIter);
    T sum = 0;
    for(unsigned i = 0; i < totalIter; ++i)
    {
        sum += fpg.SingleRun();
    }
    return sum/T(totalIter);
}

template <typename T>
class Worker
{
    private:
        T meanPrice, priceSd;
        unsigned maxDayNum, maxIter, localIter;
        T resultValue;
    public:
        Worker(const T& inMeanPrice, const T& inPriceSd, const unsigned inMaxDayNum, const unsigned inMaxIter, const unsigned inLocalIter) :
        meanPrice(inMeanPrice),
        priceSd(inPriceSd),
        maxDayNum(inMaxDayNum),
        maxIter(inMaxIter),
        localIter(inLocalIter),
        resultValue(-1.0)
        {}
        void operator()()
        {
            this->resultValue = mcFuelPriceSingle<T>(this->meanPrice, this->priceSd, this->maxDayNum,this->maxIter, this->localIter);
        }
        T GetResultValue() const
        {
            return this->resultValue;
        }
};

template <typename T>
T mcFuelPricePar(const T& meanPrice, const T& priceSd, const unsigned maxDayNum, const unsigned maxIter, const unsigned totalIter, const unsigned numberOfThreads)
{
    const unsigned localIter = totalIter/numberOfThreads;
    std::vector<std::thread> threadVector;
    threadVector.reserve(numberOfThreads);
    std::vector<Worker<T> > workerVector;
    workerVector.reserve(numberOfThreads);
    for(unsigned i = 0; i < numberOfThreads; ++i)
    {
        workerVector.push_back(Worker<T>(meanPrice, priceSd, maxDayNum, maxIter, localIter));
        threadVector.push_back(std::thread(std::ref(workerVector[i])));
    }
    for(auto & th : threadVector)
    {
        th.join();
    }

    auto sumFun = [](const T& tmpSum, const Worker<T>& inWorker)
    {
        return tmpSum+inWorker.GetResultValue();
    };
    T sum = std::accumulate(workerVector.begin(), workerVector.end(), T(0.0), sumFun);
    return sum/T(numberOfThreads);
}

int main()
{
    const double meanPrice = 3.0;
    const double priceSd = 0.25;
    const unsigned maxDayNum = 30;
    const unsigned maxIter = 400; //used for deciding, if one should buy fuel on a given day
    const unsigned totalIter = 400000; //total number of Monte-Carlo iterations
    const unsigned threadNum = std::thread::hardware_concurrency()/2 > 0 ? std::thread::hardware_concurrency()/2 : 1;
    auto startTime = std::chrono::steady_clock::now();
    std::cout<<"Result: "<<mcFuelPricePar(meanPrice, priceSd, maxDayNum, maxIter, totalIter, threadNum)<<std::endl;
    auto endTime = std::chrono::steady_clock::now();
    std::chrono::duration<double> runTime = endTime-startTime;
    std::cout<<"Runtime: "<<runTime.count()<<" [s]"<<std::endl;

    return 0;
}