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
|
# Copyright (C) 2018 The Regents of the University of California, Michael Ludkovski and Aditya Maheshwari
# All Rights Reserved
# This code is published under the GNU Lesser General Public License (GNU LGPL)
from __future__ import division
import numpy as np
import math
import matplotlib.pyplot as plt
import microgridDEA.parameters as bv
import microgridDEA.parameters as bp
def residualDemand(param, nstep, nbsimulOpt, T, X0, typeRMC, sampleType="insample"):
lambd = param.lambd
sigma = param.sigma
A = param.ampl
mu = param.mu
dt = T/nstep
if sampleType=="insample":
np.random.seed(seed=5)
elif sampleType=="outsample":
np.random.seed(seed=96)
if typeRMC == 'regress now 2D':
loadMatrix = np.zeros((nbsimulOpt,nstep+1))
loadMatrix[:,0] = X0
for i in range(1,nstep+1):
dW = np.random.normal(0,1,nbsimulOpt)*np.sqrt(dt)
loadMatrix[:,i] = loadMatrix[:,i-1] + lambd*(mu - loadMatrix[:,i-1])*dt + sigma*dW #sigma*loadMatrix[:,i-1]*dW
loadMatrix = np.where(loadMatrix>param.DemandUB,param.DemandUB,loadMatrix)
# inventory = np.random.uniform(param.I_minMax[0],param.I_minMax[1],nbsimulOpt)
return loadMatrix#, inventory
if typeRMC == 'gd':
loadMatrix = np.zeros((nbsimulOpt,nstep+1))
loadMatrix[:,0] = X0
for i in range(1,nstep+1):
dW = np.random.normal(0,1,nbsimulOpt)*np.sqrt(dt)
loadMatrix[:,i] = loadMatrix[:,i-1] + lambd*(mu - loadMatrix[:,i-1])*dt + sigma*dW #sigma*loadMatrix[:,i-1]*dW
# inventory = np.linspace(param.I_minMax[0],param.I_minMax[1],param.meshI+1)
loadMatrix = np.where(loadMatrix>param.DemandUB,param.DemandUB,loadMatrix)
return loadMatrix#, inventory
if __name__ == '__main__':
param = bp.basicVariables()
nbsimulOpt =5
X0 = np.linspace(-param.ampl,param.ampl,nbsimulOpt)
X0[:]=0
# print len(X0)
demand, inventory = residualDemand(param, param.nstep, nbsimulOpt, param.maturity, X0, 'gd')
print(len(inventory))
# plt.figure(2)
for simcnt in range(nbsimulOpt):
plt.plot(demand[simcnt,:])
plt.grid(True)
plt.show()
# param = bv.basicVariables()
# print param.B_minMax, param.I_minMax[0],param.I_minMax[1], param.d_minMax, param.K
# X0 = np.arange(-10,10,20/param.nbsimulOpt)
# # matrix of prices with rows as (number of simulations) and columns as (number of time steps+1)
# demand, inventory, control = residualDemand(param, param.nstep, param.nbsimulOpt, param.maturity, X0, 'regress now 2D')
# for i in range(0,100):
# plt.plot(demand[i,:])
# plt.show()
|
