蒙特卡罗容差模拟

Question

我是新的工作蒙特卡罗模拟，我试图找到最佳的公差适合两个部分。该模型应产生符合规格的部件的概率。这应该很简单，但我不知道如何使用truncnorm.rvs：

n = 1000000 # number of trials

Component specifications
sigma = 3       # assume 99.7% of components in tolerance spec
dp = 3.900       # (Part Number 721634) diameter
tol_dpmax = 3.9015   # (Part Number 721634) diameter tolerance + (0.0015)
tol_dpmin = 3.9000   # 721634 diameter tolerance - (-0.0)
do = 3.900          # d ring diameter 
tol_domax = 3.90   # d ring diameter tolerance + (+0.0)
tol_domin = 3.8995   # d ring diameter tolerance - (-0.0005)

dp_std= (tol_dpmax-tol_dpmin/2/3)
do_std= (tol_domax-tol_domin/2/3)
def calculateInterference(dp, do):
    """
    in: diameters
    out: o ring interference 
    """
    
    interference = dp - do 
    return interference

def checkInterference(interference, lower, upper):
    """
    in: o ring interference
    out: Bool, T for pass, F for fail
    """
    
    if interference < lower or interference > upper:
        return False
    else:
        return True

# set interference limits
lower = 0.0 # smallest permissible interference 
upper = 0.0005 # largest permissible interference
      
interference_list = []  # initialize a list to store interference values per trial 
results = []            # initialize a list to store our trial results

# Main simulation loop
for trial in range(n):
    # Sample PDF for each component
   
    piston_sample = scipy.stats.truncnorm.rvs((tol_dpmin-dp)/sigma,(tol_dpmax-dp)/sigma)
    #piston_sample = scipy.stats.truncnorm.rvs((tol_dpmin-dp)/sigma,(tol_dpmax-dp)/sigma,dp,dp_std) #good one ish

    
    #oring_sample = np.random.normal(do, tol_domin/sigma)
    oring_sample = scipy.stats.truncnorm.rvs((tol_domin-do)/sigma,(tol_domax-do)/sigma)
    #oring_sample = scipy.stats.truncnorm.rvs((tol_domin-do)/sigma,(tol_domax-do)/sigma,do,do_std) #good one ish
    #cylinder_sample = np.random.normal(dc, tol_dc/sigma)
    
    # compute and store interference 
    interference = calculateInterference(piston_sample, oring_sample)
    interference_list.append(interference)
    
    # Trial check: log results of interference check
    results.append(checkInterference(interference, lower, upper))

# results check
goodAssemblyCount = sum(results)
failurePercentage =  100*((n - goodAssemblyCount) / n)

对此我将不胜感激。

谢谢!

Answer 1

要正确使用truncnorm，请尝试生成大量的示例并检查它们的范围：

sigma = 3
dp = 3.900
tol_dpmax = 3.9015
tol_dpmin = 3.9000

ntest = 10000
samp = scipy.stats.truncnorm.rvs((tol_dpmin-dp)/sigma,(tol_dpmax-dp)/sigma, size=ntest)
# scale & translate to original (non-normalized) location
samp_non_norm = (samp*sigma)+dp
print(f'minimum value: {samp_non_norm.min()}')
print(f'maximum value: {samp_non_norm.max()}')

import matplotlib.pyplot as plt
plt.hist(samp_non_norm, 100)

​

​

从直方图中我们可以看到截断和预期的一样工作，但是分布看起来是一致的。这会引起人们对sigma参数的注意，这是过程中的标准差。评论

# assume 99.7% of components in tolerance spec

提示sigma = 3与0.0015的耐受范围不一致。考虑重新调整模型的中心位置

dp_new = (tol_dpmin + tol_dpmax)/2

调整标准差，使单面范围为三西格玛：

dp_std = (tol_dpmax - dp_new)/sigma  # top boundary is "three sigmas" from center

并重新建模：

ntest = 10000
samp = scipy.stats.truncnorm.rvs((tol_dpmin-dp_new)/dp_std,(tol_dpmax-dp_new)/dp_std, size=ntest)
# scale & translate to original (non-normalized) location
samp_non_norm = (samp*dp_std)+dp_new
print(f'minimum value: {samp_non_norm.min()}')
print(f'maximum value: {samp_non_norm.max()}')

plt.hist(samp_non_norm, 100)

​

​

这看起来好多了！这是一个正态分布，在上下边界有轻微的截断。

..。然后重复这个过程，得到d环模型的参数。