Contact
CoCalc Logo Icon
StoreFeaturesDocsShareSupport News AboutSign UpSign In
| Download

Jupyter notebook SPC_3.ipynb

Project: SPC
Views: 59
Kernel: Python 2 (Ubuntu, plain)
import numpy as np import matplotlib.pyplot as plt %matplotlib inline

System Prototype

The basic Runto Run Controller

The W2W control scheme:

  • The yky_k is either form AVM or Metrology tool

  • AkA_k (or A) is typically chosen to be least square estimates of β1β_1 based on historical data.The control valuable is set to nullify the deviation from target.

CMP Example

  1. yky_k is the actual removal amount measured from the metrology tool and PostYkPostY_k is the actual post CMP thickness of run k. The specification of PostYkPostY_k is 2800±150 Angstrom (Å) with 2800 being the target value denoted by TgtPostYTgtPostY Deifned PreYkPreY_k be the one before CMP process, ARRkARR_k is the polish rate, μk\mu_k is the polish time that we can control!

The material removal model for CMP can be divided into two parts, mechanical model and chemical model. The chemical action of slurry is responsible for continuosly softening the silicon oxide. The fresh silicon oxide or metal surface is then rapidly removed by mechanical part.

  1. The AkA_k is the nominal removal rate, which is empirically simulated by a polynomial curve fitting of parts usage count between PMs (denoted by PU varying from 1 to 600)

  1. Process gain:

  1. Simulation parameters:

  1. Assumption:

  1. AkA_k is set to 300, and A is set to mean of ARRkARR_k (also 300)

  2. β0\beta_0 is set to 700 and b0b_0 is 500

  3. RI is simulated with 0.9 plus some variation, and assume controller can detect any value that is beyond UCL or LCL

  4. yk from AVM is virtual times yz(actual metrodlogy data), where virtual is a normal distribution variable

  5. ηk\eta_k is simulated as ARIMA(0 1 1)

  6. Control β\beta not η\eta

yz is atual measurement value, and A is used internal by R2R

from statsmodels.tsa.arima_model import ARIMA import sklearn.linear_model as skl_lm
# Simulated environment Error = np.random.normal(0, 25, 600) PM1 = np.random.normal(0, 50, 600) PM2 = np.random.normal(0, 3, 600) Stress1 = np.random.normal(1000, 55,600) Stress2 = np.random.normal(0, 4, 600) Rotspd1 = np.random.normal(100, 5.5, 600) Rotspd2 = np.random.normal(0, 0.4, 600) Sfuspd1 = np.random.normal(100, 5.5, 600) Sfuspd2 = np.random.normal(0, 0.4, 600) PreY = np.random.normal(3800, 70, 600) TgtPostY = np.array([2800]*600) # Assumption Virtual = np.random.normal(1, np.sqrt(0.01), 600) Ak0 = 1200 # Not sure, does not know PU Ak = np.array([1200]*600) RI = 0.9 etak0 = 0 etak = 500 TgtY = PreY - TgtPostY def ARR (stress1, stress2, rotspd1, rotspd2, sfuspd1, sfuspd2, pm1, pm2, error): return Ak0*((stress1+stress2)/1000)*((rotspd1+rotspd2)/100)*((sfuspd1+sfuspd2)/100)+pm1+pm2+error ARRk = ARR(Stress1, Stress2, Rotspd1, Rotspd2, Sfuspd1, Sfuspd2, PM1, PM2, Error) # beta1 actual process gain RI=RI*np.sqrt(np.sqrt(Virtual)) # ARIMA (0 1 1) assume error model = ARIMA(Error, order=(0, 1, 1)) results_AR = model.fit(disp=-1) alpha1 = 0.35 etak = np.zeros(len(PreY)) mu = np.zeros(len(PreY)) PostY = np.zeros(len(PreY)) yz = np.zeros(len(PreY)) mu[0] = (TgtY[0] - etak0)/Ak[0] yz[0] = ARRk[0]*mu[0]+ results_AR.fittedvalues[0] # Actual value PostY[0] = PreY[0] - yz[0] etak[1] = alpha1*(yz[0]-Ak[0]*mu[0])+(1-alpha1)*etak0 mu[1] = (TgtY[1] - etak[1])/Ak[1] yz[1] = ARRk[1]*mu[1]+ results_AR.fittedvalues[1] # Actual value PostY[1] = PreY[1] - yz[1] etak[2] = alpha1*(yz[1]-Ak[1]*mu[1])+(1-alpha1)*etak[1] mu[2] = (TgtY[2] - etak[2])/Ak[2] for k in range(2,599): yz[k] = ARRk[k]*mu[k]+ 700 + results_AR.fittedvalues[k] ## ARRk[k], Toub[k] approximate beta_1, beta_0(not tou!!) PostY[k] = PreY[k] - yz[k] etak[k+1] = alpha1*(yz[k]-Ak[k]*mu[k])+(1-alpha1)*etak[k] mu[k+1] = (TgtY[k+1] - etak[k+1])/Ak[k+1] # Regression coefficients (Ordinary Least Squares) mu = mu.reshape((600,1)) yz = yz.reshape((600,1)) regr = skl_lm.LinearRegression() regr.fit(mu,yz) print(regr.intercept_) print(regr.coef_) print(regr.score(mu,yz))
[ 755.61616784] [[ 962.33970185]] 0.531339007585
# Simulated environment Error = np.random.normal(0, 25, 600) PM1 = np.random.normal(0, 50, 600) PM2 = np.random.normal(0, 3, 600) Stress1 = np.random.normal(1000, 55,600) Stress2 = np.random.normal(0, 4, 600) Rotspd1 = np.random.normal(100, 5.5, 600) Rotspd2 = np.random.normal(0, 0.4, 600) Sfuspd1 = np.random.normal(100, 5.5, 600) Sfuspd2 = np.random.normal(0, 0.4, 600) PreY = np.random.normal(3800, 70, 600) TgtPostY = np.array([2800]*600) # Assumption Virtual = np.random.normal(1, np.sqrt(0.01), 600) Ak0 = 1200 # Not sure, does not know PU Ak = np.array([1200]*600) RI = 0.9 etak0 = 0 etak = 500 TgtY = PreY - TgtPostY def ARR (stress1, stress2, rotspd1, rotspd2, sfuspd1, sfuspd2, pm1, pm2, error): return Ak0*((stress1+stress2)/1000)*((rotspd1+rotspd2)/100)*((sfuspd1+sfuspd2)/100)+pm1+pm2+error ARRk = ARR(Stress1, Stress2, Rotspd1, Rotspd2, Sfuspd1, Sfuspd2, PM1, PM2, Error) # beta1 actual process gain # ARIMA (0 1 1) assume error model = ARIMA(Error, order=(0, 1, 1)) results_AR = model.fit(disp=-1)
RI=RI*np.sqrt(np.sqrt(Virtual))

Round 1

First two rounds, we use actual metrodlogy value

alpha1 = 0.35 etak = np.zeros(len(PreY)) mu = np.zeros(len(PreY)) PostY = np.zeros(len(PreY)) yz = np.zeros(len(PreY))
mu[0] = (TgtY[0] - etak0)/Ak[0] yz[0] = ARRk[0]*mu[0]+ results_AR.fittedvalues[0] # Actual value PostY[0] = PreY[0] - yz[0] etak[1] = alpha1*(yz[0]-Ak[0]*mu[0])+(1-alpha1)*etak0 mu[1] = (TgtY[1] - etak[1])/Ak[1] print etak[1] print yz[0] print mu[0] print mu[1] print PostY[0]
9.4087832029 1064.34565718 0.864552849547 0.872436677123 2773.11776228

Round2

yz[1] = ARRk[1]*mu[1]+ results_AR.fittedvalues[1] # Actual value PostY[1] = PreY[1] - yz[1] etak[2] = alpha1*(yz[1]-Ak[1]*mu[1])+(1-alpha1)*etak[1] mu[2] = (TgtY[2] - etak[2])/Ak[2] print etak[2] print yz[1] print mu[2] print PostY[1]
45.0780102323 1158.24487298 0.790866989964 2698.08792277

CASE1 R2R with in-situ metrology

𝛼1𝛼_1 set to 0.35, and all actual metrology data are available

for k in range(2,599): yz[k] = ARRk[k]*mu[k]+ results_AR.fittedvalues[k] ## ARRk[k], Toub[k] approximate beta_1, beta_0(not tou!!) PostY[k] = PreY[k] - yz[k] etak[k+1] = alpha1*(yz[k]-Ak[k]*mu[k])+(1-alpha1)*etak[k] mu[k+1] = (TgtY[k+1] - etak[k+1])/Ak[k+1]
for k in range(0, 7): print etak[k+1]
9.4087832029 45.0780102323 24.864613292 5.37253882486 8.19373790403 38.2825699152 65.735506947
for k in range(0, 7): print mu[k+1]
0.872436677123 0.790866989964 0.881409276899 0.822494223586 0.760932003591 0.769812147481 0.858734964209
for k in range(0, 7): print yz[k]
1064.34565718 1158.24487298 936.365835502 1026.86410424 1000.42617593 1007.28023367 1040.49410984
for k in range(0, 7): print PostY[k]
2773.11776228 2698.08792277 2857.75256269 2855.69164133 2791.9394312 2714.03190854 2721.56303705
x = range(2, 590) plt.plot(x,etak[2:590])
[<matplotlib.lines.Line2D at 0x7f2d3dbdc690>]
Image in a Jupyter notebook
# Plotting x = range(10, 590) U = np.empty(580) L = np.empty(580) UCL = 2950 LCL = 2650 U.fill(2950) L.fill(2650) plt.plot(x,U) plt.plot(x,L) plt.plot(x,PostY[10:590]) plt.ylim([2600, 3000]) np.mean(PostY[10:590])
2800.5189283424916
Image in a Jupyter notebook

CpK(Process Capability)

Cpk=np.zeros(len(PreY)) MAPEp =np.zeros(len(PreY))
for k in range(2,599): Cpk[k] = min((UCL-np.mean(PostY[1:k]))/(3*np.std(PostY[1:k], dtype=np.float64)), (np.mean(PostY[1:k])-LCL)/(3*np.std(PostY[1:k], dtype=np.float64))) plt.plot(x,Cpk[10:590]) plt.ylim([0.3, 2.0])
(0.3, 2.0)
Image in a Jupyter notebook

Mean-absolute-percentage error (MAPEp)

for k in range(2,599): MAPEp[k] = sum(np.absolute((PostY[1:k]-2800)/2800))/k*100 plt.plot(x,MAPEp[10:590]) plt.ylim([0.7, 2.5])
(0.7, 2.5)
Image in a Jupyter notebook

CASE2 R2R+VM without RI

𝛼2=𝛼1=0.35𝛼_2=𝛼_1=0.35, apply to following wafers

# Reinitialized alpha1 = 0.35 betak = np.zeros(len(PreY)) mu = np.zeros(len(PreY)) PostY = np.zeros(len(PreY)) PostYv = np.zeros(len(PreY)) y = np.zeros(len(PreY)) yz = np.zeros(len(PreY)) # Round1 mu[0] = (TgtY[0] - betak0)/Ak[0] yz[0] = ARRk[0]*mu[0]+beta0+results_AR.fittedvalues[0] # Actual value PostY[0] = PreY[0] - yz[0] betak[1] = alpha1*(yz[0]-Ak[0]*mu[0])+(1-alpha1)*betak0 mu[1] = (TgtY[1] - betak[1])/Ak[1] # Round2 yz[1] = ARRk[1]*mu[1]+beta0+results_AR.fittedvalues[1] # Actual value PostY[1] = PreY[1] - yz[1] betak[2] = alpha1*(yz[1]-Ak[1]*mu[1])+(1-alpha1)*betak[1] mu[2] = (TgtY[2] - betak[2])/Ak[2]
for k in range(2,599): yz[k] = ARRk[k]*mu[k]+results_AR.fittedvalues[k]+beta0 # Metrology data if(k%25)==1: y[k] = yz[k] else: y[k] = yz[k]*Virtual[k] # Assume VM value ok ARR=ARR_bar PostY[k] = PreY[k] - yz[k] betak[k+1] = alpha1*(y[k]-Ak[k]*mu[k])+(1-alpha1)*betak[k] mu[k+1] = (TgtY[k+1] - betak[k+1])/Ak[k+1]
x = range(2, 590) plt.plot(x,betak[2:590])
[<matplotlib.lines.Line2D at 0x7f68a7814c50>]
Image in a Jupyter notebook
x = range(10, 590) U = np.empty(580) L = np.empty(580) U.fill(2950) L.fill(2650) plt.plot(x,U) plt.plot(x,L) plt.plot(x,PostY[10:590]) np.mean(PostY[10:590]) plt.ylim([2600, 3000])
(2600, 3000)
Image in a Jupyter notebook
Cpk=np.zeros(len(PreY)) MAPEp =np.zeros(len(PreY))
for k in range(2,599): Cpk[k] = min((UCL-np.mean(PostY[1:k]))/(3*np.std(PostY[1:k], dtype=np.float64)), (np.mean(PostY[1:k])-LCL)/(3*np.std(PostY[1:k], dtype=np.float64))) plt.plot(x,Cpk[10:590]) plt.ylim([0.3, 2.0])
(0.3, 2.0)
Image in a Jupyter notebook
for k in range(2,599): MAPEp[k] = sum(np.absolute((PostY[1:k]-2800)/2800))/k*100 plt.plot(x,MAPEp[10:590]) plt.ylim([0.7, 2.5])
(0.7, 2.5)
Image in a Jupyter notebook

CASE3: R2R+VM with RI

# Reinitialized alpha1 = 0.35 betak = np.zeros(len(PreY)) mu = np.zeros(len(PreY)) PostY = np.zeros(len(PreY)) PostYv = np.zeros(len(PreY)) y = np.zeros(len(PreY)) yz = np.zeros(len(PreY)) # Round1 mu[0] = (TgtY[0] - betak0)/Ak[0] yz[0] = ARRk[0]*mu[0]+beta0+results_AR.fittedvalues[0] # Actual value PostY[0] = PreY[0] - yz[0] betak[1] = alpha1*(yz[0]-Ak[0]*mu[0])+(1-alpha1)*betak0 mu[1] = (TgtY[1] - betak[1])/Ak[1] # Round2 yz[1] = ARRk[1]*mu[1]+beta0+results_AR.fittedvalues[1] # Actual value PostY[1] = PreY[1] - yz[1] betak[2] = alpha1*(yz[1]-Ak[1]*mu[1])+(1-alpha1)*betak[1] mu[2] = (TgtY[2] - betak[2])/Ak[2]
for k in range(2,599): yz[k] = ARRk[k]*mu[k]+beta0+results_AR.fittedvalues[k] # process output how to filter out if(k%25)==1: y[k] = yz[k] alpha = alpha1 else: y[k] = yz[k]*Virtual[k] alpha = RI[k]*alpha1 PostYv[k] = PreY[k] - yz[k]*Virtual[k] if PostYv[k]>2950 or PostYv[k]<2650: alpha = 0 betak[k+1] = alpha*(y[k]-Ak[k]*mu[k])+(1-alpha)*betak[k] #Toub canceled out mu[k+1] = (TgtY[k+1] - betak[k+1])/Ak[k+1] PostY[k] = PreY[k] - yz[k] # Process output
x = range(2, 590) plt.plot(x,betak[2:590])
[<matplotlib.lines.Line2D at 0x7f68a74fb850>]
Image in a Jupyter notebook
x = range(10, 590) U = np.empty(580) L = np.empty(580) U.fill(2950) L.fill(2650) # plot(x, PostY[3:600], type="o") plt.plot(x,U) plt.plot(x,L) plt.plot(x,PostY[10:590]) np.mean(PostY[10:590]) plt.ylim([2600, 3000])
(2600, 3000)
Image in a Jupyter notebook
Cpk=np.zeros(len(PreY)) MAPEp =np.zeros(len(PreY))
for k in range(2,599): Cpk[k] = min((UCL-np.mean(PostY[1:k]))/(3*np.std(PostY[1:k], dtype=np.float64)), (np.mean(PostY[1:k])-LCL)/(3*np.std(PostY[1:k], dtype=np.float64))) plt.plot(x,Cpk[10:590]) plt.ylim([0.3, 2.0])
(0.3, 2.0)
Image in a Jupyter notebook
for k in range(2,599): MAPEp[k] = sum(np.absolute((PostY[1:k]-2800)/2800))/k*100 plt.plot(x,MAPEp[10:590]) plt.ylim([0.7, 2.5])
(0.7, 2.5)
Image in a Jupyter notebook

CASE4: R2R+VM with (1-RI)

# Reinitialized alpha1 = 0.35 betak = np.zeros(len(PreY)) mu = np.zeros(len(PreY)) PostY = np.zeros(len(PreY)) PostYv = np.zeros(len(PreY)) y = np.zeros(len(PreY)) yz = np.zeros(len(PreY)) # Round1 mu[0] = (TgtY[0] - betak0)/Ak[0] yz[0] = ARRk[0]*mu[0]+beta0+results_AR.fittedvalues[0] # Actual value PostY[0] = PreY[0] - yz[0] betak[1] = alpha1*(yz[0]-Ak[0]*mu[0])+(1-alpha1)*betak0 mu[1] = (TgtY[1] - betak[1])/Ak[1] # Round2 yz[1] = ARRk[1]*mu[1]+beta0+results_AR.fittedvalues[1] # Actual value PostY[1] = PreY[1] - yz[1] betak[2] = alpha1*(yz[1]-Ak[1]*mu[1])+(1-alpha1)*betak[1] mu[2] = (TgtY[2] - betak[2])/Ak[2]
for k in range(2,599): yz[k] = ARRk[k]*mu[k]+beta0+results_AR.fittedvalues[k] # process output how to filter out if(k%25)==1: y[k] = yz[k] alpha = alpha1 else: y[k] = yz[k]*Virtual[k] alpha = (1-RI[k])*alpha1 PostYv[k] = PreY[k] - yz[k]*Virtual[k] if PostYv[k]>2950 or PostYv[k]<2650: alpha = 0 betak[k+1] = alpha*(y[k]-Ak[k]*mu[k])+(1-alpha)*betak[k] #Toub canceled out mu[k+1] = (TgtY[k+1] - betak[k+1])/Ak[k+1] PostY[k] = PreY[k] - yz[k] # Process output
x = range(2, 590) plt.plot(x,betak[2:590])
[<matplotlib.lines.Line2D at 0x7f68a77adbd0>]
Image in a Jupyter notebook
x = range(10, 590) U = np.empty(580) L = np.empty(580) U.fill(2950) L.fill(2650) # plot(x, PostY[3:600], type="o") plt.plot(x,U) plt.plot(x,L) plt.plot(x,PostY[10:590]) np.mean(PostY[10:590]) plt.ylim([2600, 3000])
(2600, 3000)
Image in a Jupyter notebook
Cpk=np.zeros(len(PreY)) MAPEp =np.zeros(len(PreY))
for k in range(2,599): Cpk[k] = min((UCL-np.mean(PostY[1:k]))/(3*np.std(PostY[1:k], dtype=np.float64)), (np.mean(PostY[1:k])-LCL)/(3*np.std(PostY[1:k], dtype=np.float64))) plt.plot(x,Cpk[10:590]) plt.ylim([0.3, 2.0])
(0.3, 2.0)
Image in a Jupyter notebook
for k in range(2,599): MAPEp[k] = sum(np.absolute((PostY[1:k]-2800)/2800))/k*100 plt.plot(x,MAPEp[10:590]) plt.ylim([0.7, 2.5])
(0.7, 2.5)
Image in a Jupyter notebook

CASE5: R2R+VM with RI.S.(1-RI)

# Reinitialized alpha1 = 0.35 betak = np.zeros(len(PreY)) mu = np.zeros(len(PreY)) PostY = np.zeros(len(PreY)) PostYv = np.zeros(len(PreY)) y = np.zeros(len(PreY)) yz = np.zeros(len(PreY)) # Round1 mu[0] = (TgtY[0] - betak0)/Ak[0] yz[0] = ARRk[0]*mu[0]+beta0+results_AR.fittedvalues[0] # Actual value PostY[0] = PreY[0] - yz[0] betak[1] = alpha1*(yz[0]-Ak[0]*mu[0])+(1-alpha1)*betak0 mu[1] = (TgtY[1] - betak[1])/Ak[1] # Round2 yz[1] = ARRk[1]*mu[1]+beta0+results_AR.fittedvalues[1] # Actual value PostY[1] = PreY[1] - yz[1] betak[2] = alpha1*(yz[1]-Ak[1]*mu[1])+(1-alpha1)*betak[1] mu[2] = (TgtY[2] - betak[2])/Ak[2]
for k in range(2,599): yz[k] = ARRk[k]*mu[k]+beta0+results_AR.fittedvalues[k] PostYv[k] = PreY[k]- yz[k]*Virtual[k] if k>= 25: if PostYv[k]>2950 or PostYv[k]<2650: # should be replace with RI, GSI alpha = 0 else: alpha = (1-RI[k])*alpha1 else: if PostYv[k]>2950 or PostYv[k]<2650: alpha = 0 else: alpha = RI[k]*alpha1 if(k%25)==1: y[k] = yz[k] alpha = alpha1 else: y[k] = yz[k]*Virtual[k] betak[k+1] = alpha*(y[k]-Ak[k]*mu[k])+(1-alpha)*betak[k] #Toub canceled out mu[k+1] = (TgtY[k+1] - betak[k+1])/Ak[k+1] PostY[k] = PreY[k] - yz[k] # Process output
x = range(2, 590) plt.plot(x,betak[2:590])
[<matplotlib.lines.Line2D at 0x7f68a76ad550>]
Image in a Jupyter notebook
x = range(10, 590) U = np.empty(580) L = np.empty(580) U.fill(2950) L.fill(2650) # plot(x, PostY[3:600], type="o") plt.plot(x,U) plt.plot(x,L) plt.plot(x,PostY[10:590]) np.mean(PostY[10:590]) plt.ylim([2600, 3000])
(2600, 3000)
Image in a Jupyter notebook
Cpk=np.zeros(len(PreY)) MAPEp =np.zeros(len(PreY))
for k in range(2,599): Cpk[k] = min((UCL-np.mean(PostY[1:k]))/(3*np.std(PostY[1:k], dtype=np.float64)), (np.mean(PostY[1:k])-LCL)/(3*np.std(PostY[1:k], dtype=np.float64))) plt.plot(x,Cpk[10:590]) plt.ylim([0.3, 2.0])
(0.3, 2.0)
Image in a Jupyter notebook
for k in range(2,599): MAPEp[k] = sum(np.absolute((PostY[1:k]-2800)/2800))/k*100 plt.plot(x,MAPEp[10:590]) plt.ylim([0.7, 2.5])
(0.7, 2.5)
Image in a Jupyter notebook