Добавил:
Опубликованный материал нарушает ваши авторские права? Сообщите нам.
Вуз: Предмет: Файл:

VT Mod 3

.py
Скачиваний:
1
Добавлен:
05.03.2020
Размер:
4.22 Кб
Скачать
import random as r
import math as m
from scipy import stats
from scipy import special
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt


exp_lamda = 0.5

vg_a=1
vg_c=5

gam_lam=1
gam_a=2

#функции
def expF(x):
    return 1-m.exp(-exp_lamda*x)

def vgF(x):
    return 1-m.exp(-(x/vg_a)**vg_c)

def gamF(x):
    return stats.gamma.cdf(x,gam_a,0,gam_lam)


#плотности   
def exp(x):
    return exp_lamda*m.exp(-exp_lamda*x)

def vg(x):
    return (vg_c/vg_a)*((x/vg_a)**(vg_c-1))*m.exp(-(x/vg_a)**vg_c)


def gam(x):
    return stats.gamma.pdf(x,gam_a,0,gam_lam)
    #return ((gam_lam**gam_a)*(x**(gam_a-1))*exp(-gam_lam*x))/stats.gamma(gam_a)

def lognormalCheck(x, muC, sigmaC):
    sq2pi = m.sqrt(2*m.pi)
    return 1/(sigmaC*sq2pi)*m.exp(-m.pow((x-muC),2)/2/sigmaC/sigmaC)


#генерация чисел
def expG(x):
    return -1/exp_lamda*m.log(x)

def vgG(x):
    return vg_a*((-m.log(x))**(1/vg_c))

def gamG(x):
    return stats.gamma.rvs(gam_a,0,gam_lam)

def plotE(bins, vals):
    left,right = bins[:-1],bins[1:]
    X = np.array([left,right]).T.flatten()
    Y = np.array([vals,vals]).T.flatten()

    plt.plot(X,Y)
    plt.show()
    return

bins = []
for i in range(600):
    bins.append(i/100)
bins = np.array(bins)

width = 0.8 * (bins[1] - bins[0])
center = (bins[:-1] + bins[1:]) / 2 - width/2

val = []
for value in center:
    val.append(expF(value))
val = np.array(val)

plotE(bins, val)

val = []
for value in center:
    val.append(exp(value))
val = np.array(val)

plotE(bins, val)

val = []
for value in center:
    val.append(vgF(value))
val = np.array(val)

plotE(bins, val)

val = []
for value in center:
    val.append(vg(value))
val = np.array(val)

plotE(bins, val)

val = []
for value in center:
    val.append(gamF(value))
val = np.array(val)

plotE(bins, val)

val = []
for value in center:
    val.append(gam(value))
val = np.array(val)

plotE(bins, val)

tau = []
L = int(1 + 3.222 * m.log(5000, 10))

for i in range(5000):
    y1 = expG(r.random())
    y2 = gamG(r.random())
    y3 = vgG(r.random())
    y4 = gamG(r.random())
    y5 = gamG(r.random())
    y6 = expG(r.random())
    y7 = vgG(r.random())
    y8 = gamG(r.random())
    y9 = expG(r.random())
    y10 = vgG(r.random())

    tau.append(max(min(y1,y6,y9),min(y3,y7,y10),min(y1,y6,y8,y10),min(y3,y7,y8,y9),min(y2,y4,y6,y9),min(y2,y5,y7,y10),min(y2,y4,y6,y8,y10),min(y2,y5,y7,y8,y9),min(y3,y4,y5,y6,y9),min(y1,y4,y5,y7,y10)))


p, x = np.histogram(tau, bins=L)
pf = []
for i in range(L):
    pf.append(float(p[i]) / 5000.0)
plotE(x,pf)

cumulative = []
cumulative.append(pf[0]);
for i in range(1,len(pf)):
    cumulative.append(pf[i]+cumulative[i-1])
    
plotE(x, cumulative)
    
width = 0.9 * (x[1] - x[0])
center = (x[:-1] + x[1:]) / 2 - width/2
plt.bar(center, pf, align='center', width=width)
plt.show()

mTau = 0
dTau = 0
for i in range(L):
    mTau += x[i]*p[i]
mTau /= 5000
    
for i in range(L):
    dTau += (x[i] - mTau) * (x[i] - mTau)
dTau/=L
dTau=m.sqrt(dTau)

print(L)
print("X:", mTau)
print("sigma:",dTau)

tal99=2.5582
tal975=2.1791
tal95=1.6821
tal90=1.2715
tal80=0.7649
tal60=0.2298
tal55=0.1134

print("P=0.80:")
print("left:", mTau-tal90*dTau/m.sqrt(L-1))
print("right:", mTau+tal90*dTau/m.sqrt(L-1))

expected = []
for i in range(L):
    expected.append(lognormalCheck(x[i], mTau, dTau))
    
print("Chi-square:",stats.chisquare(pf, expected))
iofunc = []
x = []
y = []
y.append(expG(r.random()))
y.append(gamG(r.random()))
y.append(vgG(r.random()))
y.append(vgG(r.random()))
y.append(gamG(r.random()))
y.append(gamG(r.random()))
y.append(expG(r.random()))
y.append(expG(r.random()))
y.append(vgG(r.random()))
y.append(gamG(r.random()))

for i in range(1,10001):
    broken=0
    for j in range(9):
        if (i/2000 > y[j]):
            broken+=1
    if broken != 0:
        x.append(i/2000)
        if broken != 9:
            iofunc.append(broken/i*2000/(9-broken))
        else:
            iofunc.append(4)
            
plt.plot(x,iofunc)
plt.show()
Соседние файлы в предмете Моделирование систем