CS计算机代考程序代写 In [1]:

In [1]:
import numpy as np
import matplotlib.pyplot as plt
import math
import scipy.special as spsp

import scipy.stats as spst

import simpy
%matplotlib inline
In [2]:
def arrival(env,server,interarrival,servicetime,departtime,stop):
for i in range(20):
yield env.timeout(interarrival[i])
env.process(service(env, server, servicetime, departtime,i,stop))

def service(env, server, servicetime, departtime,i,stop):
rqt=server.request()
yield rqt
yield env.timeout(servicetime[i])
server.release(rqt)
departtime.append(env.now)

def system():
env=simpy.Environment()
server=simpy.Resource(env)
stop=simpy.Event(env)
u1=np.random.rand(20)
u2=np.random.rand(20)
lmbda1,lmbda2=10,11
interarrival=-1/lmbda1*np.log(u1)
servicetime=-1/lmbda2*np.log(u2)
departtime=[]
env.process(arrival(env,server,interarrival,servicetime,departtime,stop))
arrivaltimes=np.cumsum(interarrival)
env.run()
average_time=np.mean(np.array(departtime)-arrivaltimes)
return average_time

N=1000
ATs=[system() for i in range(N)]
ATs=np.array(ATs)
mean=np.mean(ATs)
std=np.std(ATs,ddof=1)
lower=mean-spst.norm.ppf(0.975)*std/N**0.5
upper=mean+spst.norm.ppf(0.975)*std/N**0.5
print(“95% confidence interval”, lower,upper)
print(“mean”, mean)
print(“width of confidence interval”,upper-lower)

95% confidence interval 0.26603852706140946 0.28611443438848966
mean 0.27607648072494956
width of confidence interval 0.0200759073270802

Now, let’s apply the control variate method
• Step 1: Generate a set of $x$ samples and also a set of $z$ samples. 
• Step 2: estimate $c=-\frac{Cov(x,z)}{ Var(z)}$ based on the samples we have 
• Step 3: Generate $y$ samples usign $x+c(z-\mu_z)$ 
• Step 4: Construct the confidence interval using the $y$ samples 
In [0]:
def arrival(env,server,interarrival,servicetime,departtime,stop):
for i in range(20):
yield env.timeout(interarrival[i])
env.process(service(env, server, servicetime, departtime,i,stop))

def service(env, server, servicetime, departtime,i,stop):
rqt=server.request()
yield rqt
yield env.timeout(servicetime[i])
server.release(rqt)
departtime.append(env.now)

N=1000
result=[system() for i in range(N)]
result=np.array(result)
x=result[:,0]
z=result[:,1]
c=-np.cov(x,z,ddof=1)[0,1]/np.var(z,ddof=1)
y=x+c*(z-1/11)

#ATs=np.array(ATs)
mean=np.mean(y)
std=np.std(y,ddof=1)
lower=mean-spst.norm.ppf(0.975)*std/N**0.5
upper=mean+spst.norm.ppf(0.975)*std/N**0.5
print(“95% confidence interval”, lower,upper)
print(“mean”, mean)
print(“width of confidence interval”,upper-lower)
In [0]:
def arrival(env,server,interarrival,servicetime,departtime,stop):
for i in range(20):
yield env.timeout(interarrival[i])
env.process(service(env, server, servicetime, departtime,i,stop))

def service(env, server, servicetime, departtime,i,stop):
rqt=server.request()
yield rqt
yield env.timeout(servicetime[i])
server.release(rqt)
departtime.append(env.now)

N=1000
result=[system() for i in range(N)]
result=np.array(result)
x=result[:,0]
z1=result[:,1]
z2=result[:,2]
c1=-np.cov(x,z1,ddof=1)[0,1]/np.var(z1,ddof=1)
c2=-np.cov(x,z2,ddof=1)[0,1]/np.var(z2,ddof=1)

y=x+c1*(z1-1/11)+c2*(z2-1/10)

Related Posts