CN113189967A - Control performance diagnosis method for semiconductor process batch control system - Google Patents

Control performance diagnosis method for semiconductor process batch control system Download PDF

Info

Publication number
CN113189967A
CN113189967A CN202110491048.7A CN202110491048A CN113189967A CN 113189967 A CN113189967 A CN 113189967A CN 202110491048 A CN202110491048 A CN 202110491048A CN 113189967 A CN113189967 A CN 113189967A
Authority
CN
China
Prior art keywords
model
window
disturbance
representing
data
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202110491048.7A
Other languages
Chinese (zh)
Other versions
CN113189967B (en
Inventor
凌丹
栗朝松
王妍
孙军伟
王英聪
郭群力
王延峰
张勋才
刘鹏
姜素霞
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhengzhou University of Light Industry
Original Assignee
Zhengzhou University of Light Industry
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Zhengzhou University of Light Industry filed Critical Zhengzhou University of Light Industry
Priority to CN202110491048.7A priority Critical patent/CN113189967B/en
Publication of CN113189967A publication Critical patent/CN113189967A/en
Application granted granted Critical
Publication of CN113189967B publication Critical patent/CN113189967B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G05—CONTROLLING; REGULATING
    • G05B—CONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
    • G05B23/00—Testing or monitoring of control systems or parts thereof
    • G05B23/02—Electric testing or monitoring
    • G05B23/0205—Electric testing or monitoring by means of a monitoring system capable of detecting and responding to faults
    • G05B23/0218—Electric testing or monitoring by means of a monitoring system capable of detecting and responding to faults characterised by the fault detection method dealing with either existing or incipient faults
    • G05B23/0243—Electric testing or monitoring by means of a monitoring system capable of detecting and responding to faults characterised by the fault detection method dealing with either existing or incipient faults model based detection method, e.g. first-principles knowledge model
    • G—PHYSICS
    • G05—CONTROLLING; REGULATING
    • G05B—CONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
    • G05B2219/00—Program-control systems
    • G05B2219/20—Pc systems
    • G05B2219/24—Pc safety
    • G05B2219/24065—Real time diagnostics

Abstract

The invention provides a control performance diagnosis method of a semiconductor process batch control system, which comprises the following steps of firstly, estimating reference working condition white noise, a reference working condition disturbance model, reference working condition model residual error, actual working condition white noise, an actual working condition disturbance model and actual working condition model residual error according to closed loop data of historical normal working conditions and closed loop data of actual working conditions; secondly, acquiring four diagnostic monitoring quantities according to the reference working condition white noise and the model residual thereof, the actual working condition white noise and the model residual thereof, and acquiring the control performance monitoring quantity and the control upper limit/lower limit of the four diagnostic monitoring quantities according to the closed loop data of the historical normal working condition; and finally, evaluating the control performance of the batch control system according to the control performance monitoring quantity and the control upper limit/lower limit thereof. The method not only can correctly evaluate the control performance of the closed-loop system under the conventional operating condition of the semiconductor manufacturing process, but also can effectively diagnose the reason causing the performance reduction of the system, thereby reducing the cost of system maintenance and improving the safety of the system.

Description

Control performance diagnosis method for semiconductor process batch control system
Technical Field
The invention belongs to the field of control performance evaluation and diagnosis of a process industrial control system, and particularly relates to a control performance diagnosis method of a semiconductor manufacturing process batch control system, which is used for evaluating the control performance of the semiconductor manufacturing process batch control system and diagnosing the reason of the degradation of the control performance.
Background
In recent years, with the adoption of a great market demand, abundant manpower resources and other advantages, the semiconductor industry in China has rapidly developed and has become the focus of attention of the global integrated circuit industry. The industrial level of semiconductors in China is still extremely disproportionate to the status of the big country in China. For example, high-tech enterprises in the middle and high-tech village take the majority of computer market share in China, but the profit margin is only about 3%, while the pure profit margin of the central processor manufacturer intel for years is more than 30%.
The semiconductor industry is primarily concerned with the industry of manufacturing semiconductor solid state devices and circuits, a process known as wafer fabrication. The production of wafers is a rather long and complex process involving hundreds of process steps, which cannot be performed perfectly at each time, and the occurrence of defects in any one process affects the final quality of the product. An Exponentially Weighted Moving Average (EWMA) batch-to-batch controller is an advanced controller commonly used in semiconductor manufacturing processes, and the controller adjusts manipulated variables of the next batch according to input and output data of historical batches, so that influences of process disturbances in forms of offset, drift and the like on product quality are reduced, and yield is improved. Wafer processing is often affected by internal and external uncertainties of the system, which may cause degradation of system control performance, reduction of energy, material utilization, and economic benefits, and even safety hazards. Therefore, it is important to research the performance monitoring and diagnosis technology of the EWMA batch control system in the semiconductor manufacturing process to improve the quality of the wafer product.
After the system is put into operation, the control quality of the EWMA batch controller is influenced by factors such as model mismatch, external disturbance change, controller adjustment and the like, and the system performance is reduced due to the factors. The product quality, operational safety, material energy consumption rate and economic performance of semiconductor manufacturing processes are closely related to the control performance of batch controllers. The loss of control performance from any one process can be quite dramatic, with losses of up to $ 5000 for 300 mm wafers.
Disclosure of Invention
The invention provides a control performance diagnosis method of a semiconductor manufacturing process batch control system, which aims to monitor the control performance of an EWMA batch controller and diagnose the reason of performance degradation, and solves the technical problem that the existing control performance monitoring method is difficult to determine the specific factors of control performance degradation.
The technical scheme of the invention is realized as follows:
a control performance diagnostic method for a semiconductor process batch control system comprises the following steps:
s1, obtaining N groups of input data and output data of historical normal working conditions of the EWMA batch control system in the semiconductor manufacturing process as sample data, dividing the sample data into S +1 windows, and calculating a white noise estimation value epsilon of the reference working condition by using the input data and the output data of the first windowrAnd reference working condition disturbance model estimation value
Figure BDA0003052110160000011
Obtaining process-disturbance combined model residual e (t), process data set
Figure BDA0003052110160000012
And output data set Y1 0And based on the process data set
Figure BDA0003052110160000013
Obtaining a matrix X of reference conditions1From the output data set Y1 0Obtaining a matrix Y of reference conditions1Construction of Hankel matrix XpAnd XfAnd projection matrix Jq；
S2, constructing a process data set using the input data and the output data of the k' th window
Figure BDA0003052110160000014
And outputting the data set
Figure BDA0003052110160000015
From a process data set
Figure BDA0003052110160000016
Obtaining matrix X of actual working conditionk'From the output data set
Figure BDA0003052110160000017
Obtaining matrix Y of actual working conditionsk'Obtaining the control performance monitoring quantity of the EWMA batch-to-batch control system corresponding to the kth window by adopting a typical variable analysis method
Figure BDA0003052110160000018
Wherein k ═ 2,3, …, S + 1;
s3, according to the reference working condition disturbance model estimated value in the step S1
Figure BDA0003052110160000019
Obtaining reference working condition model residual error v from input data and output data of kth windowr,k'(t) and model mismatch/perturbation dynamics joint monitoring volume
Figure BDA00030521101600000110
S4, input data using k' th windowObtaining the white noise estimation value epsilon of the actual working condition corresponding to the kth window by the output datah,k'(t) and actual condition disturbance model estimation value
Figure BDA00030521101600000111
And according to the actual working condition disturbance model estimated value
Figure BDA00030521101600000112
Obtaining the residual error v of the actual working condition modelh,k'(t)；
S5, according to the reference working condition white noise estimated value epsilon in the step S1rAnd the white noise estimate of the actual condition ε in step S4h,k'(t) calculating an interference characteristic monitoring quantity
Figure BDA0003052110160000021
According to the reference working condition model residual v in the step S3r,k'(t) and the actual condition model residual v in step S4h,k'(t) calculating the amount of monitoring of the disturbance model
Figure BDA0003052110160000022
According to the actual working condition white noise estimated value epsilon in the step S4h,k'(t) and the actual condition model residual vh,k'(t) obtaining a model mismatch monitoring quantity
Figure BDA0003052110160000023
S6, respectively obtaining the control upper limit UP of the control performance monitoring quantity by adopting a kernel density estimation methodCPIUpper control limit UP of model mismatch/disturbance dynamic joint monitoring quantityPDMIUpper limit UP of control of the monitored quantity of interference characteristicsDVILower limit of control LP of disturbance characteristic monitoring amountDVILower control limit LP of disturbance model monitoring quantityDMIAnd the upper control limit UP of the model mismatch monitoring quantityPMI；
S7, collecting N1 groups of input data and output data of the EWMA batch control system of the semiconductor process on line, and obtaining the control performance monitoring quantity eta of the current working conditionCPIAccording to the upper limit of control in step S6UPCPIJudging whether the control performance of the closed-loop system is reduced or not in real time;
s8, if the control performance monitoring quantity etaCPIControl upper limit UP greater than control performance monitoring amountCPIIf the performance of the EWMA batch control system is reduced, executing the step S9, otherwise, the EWMA batch control system has good performance;
s9, using the reference working condition disturbance model estimated value in the step S1
Figure BDA0003052110160000024
Obtaining the model mismatch/disturbance dynamic combined monitoring quantity eta of the current working condition from the N1 groups of input data and output data in the step S7PDMI；
S10, if the model mismatching/disturbance dynamic joint monitoring quantity etaPDMIGreater than the upper control limit UP of the combined monitoring quantity of model mismatch/disturbance dynamics in step S6PDMIIf the performance of the EWMA batch control system is reduced due to model mismatch or disturbance dynamics, executing the step S11, otherwise, the performance of the EWMA batch control system is reduced due to the change of the controller parameters;
s11, calculating the disturbance characteristic monitoring quantity eta of the current working conditionDVIMonitoring quantity eta of disturbance modelDMIAnd model mismatch monitor ηPMI(ii) a Disturbance characteristic monitoring quantity eta when current working conditionDVIGreater than the upper control limit UP of the disturbance characteristic monitoring amount in step S6DVIOr less than the lower control limit LP of the disturbance characteristic monitored quantityDVIWhen the interference characteristic changes, the performance of the EWMA batch control system is reduced; when the disturbance model monitors the quantity etaDMILower control limit LP less than disturbance model monitoring quantityDMIMeanwhile, the performance of the EWMA batch control system is reduced due to the change of the disturbance model; when model mismatch monitor ηPMIControl upper limit UP greater than model mismatch monitoring amountPMIModel mismatch results in performance degradation of the EWMA batch-to-batch control system.
Matrix X of the reference condition1And Y1Respectively expressed as:
Figure BDA0003052110160000025
wherein u isiRepresents the ith process input vector, i 1, …, mu，muDenotes the number of input vectors, t denotes the sampling time, ui"represents the i-th process input vector after normalization, i.e.
Figure BDA0003052110160000026
Figure BDA0003052110160000027
Means representing the ith process input vector, (σ)u,r)iRepresents the standard deviation of the ith process input vector; y isjRepresents the jth process output vector, j ═ 1, …, my，myIndicating the number of output vectors, yj"represents the normalized j-th process output vector, i.e.
Figure BDA0003052110160000028
Figure BDA0003052110160000029
Means for representing the jth process output vector, (σ)y,r)jRepresents the standard deviation of the jth process output vector; e.g. of the typej(t) denotes the jth process-perturbation joint model residual vector, where ej＝(ej(t))， ej(t)＝yj(t)-(Gm(q-1))ju(t)，(Gm(q-1))jRepresenting the process model G used in the first windowm(q-1) J-th line of (1), yj(t) denotes the jth process output at the tth sampling instant, u (t) denotes the column vector formed by all process inputs at the tth sampling instant in the first window, e ″jRepresents the normalized j process-perturbation joint model residual vector,
Figure BDA00030521101600000210
Figure BDA00030521101600000211
means for representing the jth process-perturbation combined model residual, (σ)e,r)jRepresenting the standard deviation of the jth process-perturbation combined model residual error;
the projection matrix JqComprises the following steps:
Figure BDA00030521101600000212
wherein, VqRepresenting a matrix constructed from the first q columns of the unitary matrix V,
Figure BDA00030521101600000213
represents the Hankel matrix XpThe covariance matrix of (2).
The process data set
Figure BDA00030521101600000214
And outputting the data set
Figure BDA00030521101600000215
Respectively expressed as:
Figure BDA00030521101600000216
wherein the content of the first and second substances,
Figure BDA00030521101600000217
the process input vector representing the kth window,
Figure BDA00030521101600000218
the process-perturbation joint model residual vector representing the kth window,
Figure BDA0003052110160000031
a process output vector representing the kth window;
matrix X of the actual working conditionk'And Yk'Respectively expressed as:
Figure BDA0003052110160000032
wherein u isi,k'The process input vector, u ", representing the k' th windowi,k'The process input vector representing the normalized k' th window, i.e.
Figure BDA0003052110160000033
yj,k'Represents the process output vector, y ″, of the k' th windowj,k'The process output vector representing the normalized k' th window, i.e.
Figure BDA0003052110160000034
ej,k'Process-perturbation joint model residual vector representing the kth window, ej,k'＝(ej,k'(t))， ej,k'(t)＝yj,k'(t)-(Gm,k′(q-1))juk'(t)，yj,k'(t) denotes the jth process output, u, at the tth sampling instant in the kth windowk'(t) represents the column vector formed by all process inputs at the tth sampling instant in the kth window, (G)m,k′(q-1))jRepresenting the Process model G used by the EWMA batch-to-batch controller in the k' th windowm，k′(q-1) Line j, e ″)j,k'Process-perturbation joint model residuals representing the normalized k' th window, i.e.
Figure BDA0003052110160000035
Obtaining control performance monitoring quantity of EWMA batch control system by adopting typical variable analysis method
Figure BDA0003052110160000036
Comprises the following steps:
Figure BDA00030521101600000323
wherein the content of the first and second substances,
Figure BDA0003052110160000037
indicating correspondence of the kth windowAmount of control Performance monitoring of EWMA inter-batch control System, JqRepresenting a projection matrix, pk,k'Representing matrix Xp,k'The k-th column of (1).
The model mismatch/disturbance dynamic joint monitoring quantity
Figure BDA0003052110160000038
Comprises the following steps:
Figure BDA0003052110160000039
wherein the content of the first and second substances,
Figure BDA00030521101600000310
represents the model mismatch/disturbance dynamic joint monitoring quantity corresponding to the k' th window,
Figure BDA00030521101600000311
has a value range of (0, 1)]，εr(t) the white noise estimation value of the reference condition at the t-th sampling moment, vr,k'(t) represents the reference condition model residual at the t-th sampling instant in the kth' window, | cov (εr(t)) | denotes εrDeterminant of covariance matrix of (t), | cov (v)r,k'(t)) | denotes vr,k'(t) determinant of covariance matrix;
the reference working condition model residual vr,k'(t) is:
Figure BDA00030521101600000312
wherein v isr,k'(t) represents the reference condition model residual error corresponding to the k' th window,
Figure BDA00030521101600000313
disturbance model estimation value representing reference working condition
Figure BDA00030521101600000314
Inverse of (y)k'(t) represents the process output vector at the t-th sampling instant of the k' -th window, uk'(t) represents the process input vector at the t-th sampling instant in the k' -th window, Gm,k′(q-1) Denotes the k' th windowA process model for use with a intraoral EWMA batch controller.
The actual working condition model residual error vh,k'(t) is:
Figure BDA00030521101600000315
wherein v ish,k'(t) represents the real model residual at the t-th sampling instant of the k' -th window,
Figure BDA00030521101600000316
actual condition disturbance model estimated value representing k' th window
Figure BDA00030521101600000317
Inverse of (y)k'(t) represents the process output vector at the t-th sampling instant of the k' -th window, uk'(t) represents the process input vector at the t-th sampling instant of the k' -th window, Gm,k′(q-1) Representing the process model used by the EWMA batch controller in the k' th window.
Step S5 includes the following steps:
s5.1, calculating interference characteristic monitoring quantity
Figure BDA00030521101600000318
Wherein epsilonr(t) represents the white noise estimation value of the reference condition at the t-th sampling moment, epsilonh,k'(t) represents the white noise estimate for the actual condition at the t-th sampling instant in the kth' window, | cov (εh,k'(t)) | denotes εh,k'Determinant of covariance matrix of (t), | cov (εr(t)) | denotes εr(t) determinant of covariance matrix;
s5.2, calculating the monitoring quantity of the disturbance model
Figure BDA00030521101600000319
Wherein the content of the first and second substances,
Figure BDA00030521101600000320
represents the perturbation model monitoring quantity of the k' th window,
Figure BDA00030521101600000321
has a value range of (0, 1)]If, if
Figure BDA00030521101600000322
Near or equal to 1, the perturbation model is unchanged, vr,k'(t) represents the reference operating condition model residual at the t-th sampling moment in the k' -th window, vh,k'(t) represents the actual condition model residual at the t-th sampling instant in the k' -th window, J (f (v)r,k'(t)),g(vh,k'(t))) represents vh,k'(t) and vr,k'(t) KL divergence, expressed as:
Figure BDA0003052110160000041
wherein, f (v)r,k'(t)) represents v constructed using the nuclear density estimation methodr,k'(t) probability density function, f (v)h,k'(t)) represents v obtained by the nuclear density estimation methodh,k'(t) probability density function, i ═ 1,2, …, Nf，NfRepresenting the total number of independent identically distributed sample points resulting from the kernel density estimation;
s5.3, calculating model mismatch monitoring quantity
Figure BDA0003052110160000042
Wherein epsilonh,k'(t) represents the white noise estimation value of the actual working condition at the t sampling moment in the k' th window, vh,k'(t) represents the real-world model residual at the t-th sampling instant in the kth' window, | cov (εh,k'(t)) | denotes εh,k'Determinant of covariance matrix of (t), | cov (v)h,k'(t)) | denotes vh,k'(t) determinant of covariance matrix.
Step S6 includes the following steps:
s6.1, calculating the control upper limit UP of the control performance monitoring quantityCPI：
Figure BDA0003052110160000043
Wherein the content of the first and second substances,
Figure BDA0003052110160000044
the control performance monitoring quantity corresponding to the kth window is K ', K' is 2, 3.., S +1, h represents the kernel bandwidth, and K (·) represents the kernel density estimation with the kernel as a radial basis function;
s6.2, calculating the control upper limit UP of the model mismatch/disturbance dynamic combined monitoring quantityPDMI：
Figure BDA0003052110160000045
Wherein the content of the first and second substances,
Figure BDA0003052110160000046
model mismatch/disturbance dynamic combined monitoring quantity corresponding to the kth' window;
s6.3, calculating the control upper limit UP of the interference characteristic monitoring quantityDVIAnd a lower control limit LPDVI：
Figure BDA0003052110160000047
Wherein the content of the first and second substances,
Figure BDA0003052110160000048
is the interference characteristic monitoring quantity corresponding to the kth window;
s6.4, calculating a control lower limit LP of the monitoring quantity of the disturbance modelDMI：
Figure BDA0003052110160000049
Wherein the content of the first and second substances,
Figure BDA00030521101600000410
is the disturbance model monitoring quantity corresponding to the kth window;
s6.5, calculating the control upper limit UP of the model mismatch monitoring quantityPMI：
Figure BDA00030521101600000411
Wherein the content of the first and second substances,
Figure BDA00030521101600000412
is the interference corresponding to the k' th windowFeature monitoring amount.
Step S7 includes the following steps:
s7.1, collecting N1 groups of input data and output data of the EWMA batch control system of the semiconductor process on line as actual working condition sample data, and dividing the actual working condition sample data into S1An online window, each window containing n groups of samples, wherein S1＝N1-n+1；
S7.2, calculating the residual error of the online process-disturbance combined model: e.g. of the typel(t)＝yl(t)-Gm,l(q-1)ul(t) wherein el(t) represents the on-line process-perturbation joint model residual error at the tth sampling instant of the ith window, yl(t) represents the process output vector at the tth sampling instant of the ith window, ul(t) denotes the input vector at the tth sampling instant of the l-th window, t denotes the sampling instant, Gm,l(q-1) Representing a Process model, q, used by an EWMA batch-to-batch controller in an on-line regime-1Denotes a difference operator, where l is 1,2, …, S1；
S7.3, acquiring a process input data set of the current working condition by using the process output vector, the process input vector and the process-disturbance combined model residual error which are acquired in real time
Figure BDA00030521101600000413
And output data set Yl 0：
Figure BDA00030521101600000414
Wherein the content of the first and second substances,
Figure BDA00030521101600000415
representing the on-line process input vector and,
Figure BDA00030521101600000416
representing an online process-perturbation joint model residual vector,
Figure BDA00030521101600000417
representing an online process output vector;
s7.4, constructing a matrix XlAnd Yl：
Figure BDA0003052110160000051
Wherein u isi,lRepresents the process input vector, u ″, in the ith window of the online datai,lThe process input vector representing the normalized l-th window, i.e.
Figure BDA0003052110160000052
yj,lRepresents the process output vector, y ″, of the ith window of online dataj,lThe process output vector representing the normalized l-th window, i.e.
Figure BDA0003052110160000053
ej,lProcess-perturbation joint model residual vector representing the ith window of online data, ej,l＝(ej,l(t))， ej,l(t)＝yj,l(t)-(Gm,l(q-1))jul(t)，(Gm,l(q-1))jProcess model G representing EWMA inter-batch controller usage in the ith window of online datam，l(q-1) J-th line of (1), yj,l(t) denotes the jth process output at the tth sampling instant of the ith window of online data, ul(t) represents the vector formed by all process inputs at the tth sampling time of the ith window of online data, e ″j,lRepresenting the normalized process-perturbation joint model residual,
Figure BDA0003052110160000054
s7.5, constructing a Hankel matrix Xp,lAnd Xf,l：
Figure BDA0003052110160000055
Wherein p isp+1,lRepresenting a history data vector at the p +1 th sampling instant in the ith window of line data, i.e. pp+1,l＝[xl(p)T … xl(1)T yl(p)T … yl(1)T]T，xl(p) a representation matrix
Figure BDA0003052110160000056
P column of (2), xl(1) Representation matrix
Figure BDA0003052110160000057
The p-th column of (2),
Figure BDA0003052110160000058
representation matrix XlTranspose of (y)l(p) represents a matrix Yl TP column, yl(1) Representation matrix Yl TColumn 1, Yl TRepresentation matrix YlTransposing; f. ofp+1,lThe vector of predicted data representing the p +1 th sampling instant in the l window of line data, i.e. fp+1,l＝[yl(p+1)T yl(p+2)T … yl(p+f-1)T]T，yl(p +1) represents Yl TP +1 column, yl(p +2) represents Yl TP +2 column, yl(p + f-1) represents Yl TP + f-1 column (b);
s7.6, acquiring control performance monitoring quantity eta of the current working condition by using online dataCPI：
Figure BDA0003052110160000059
Wherein the content of the first and second substances,
Figure BDA00030521101600000510
indicating the control performance monitoring quantity, p, of the ith window of the on-line datak,lRepresentation matrix Xp,lThe k thA column vector;
s7.7, control Upper Limit UP according to control Performance monitoring amountCPIJudgment of
Figure BDA00030521101600000511
Whether or not less than UPCPIIf, if
Figure BDA00030521101600000512
Less than UPCPIThen the system performance is better, if
Figure BDA00030521101600000513
Greater than UPCPIThe system performance is degraded.
Step S9 includes the following steps:
s9.1, obtaining a reference working condition model residual error by using online data:
Figure BDA00030521101600000514
wherein v isr,l(t) represents the reference operating mode model residual error at the tth sampling moment of the ith window of the online data,
Figure BDA00030521101600000515
representing the reference condition disturbance model estimate, yl(t) an online data process output vector, u, at the tth sampling instant of the ith windowl(t) an on-line data process input vector, G, representing the t-th sampling instant of the l-th windowm,l(q-1) Representing a process model used by the EWMA batch controller in the ith window of the online data;
s9.2, calculating model mismatching/disturbance dynamic combined monitoring quantity eta by utilizing online collected input data and output dataPDMI：
Figure BDA00030521101600000516
Wherein the content of the first and second substances,
Figure BDA00030521101600000517
has a value range of (0, 1)]，εr(t) reference for the t-th sampling instantEstimate of operating condition white noise, | cov (εr(t)) | denotes εr(t) determinant of covariance matrix, vr,l(t) represents the on-line data reference case model residual at the tth sampling instant of the ith window, | cov (v)r,l(t)) | denotes vr,l(t) determinant of covariance matrix;
s9.3, combining control upper limit UP of model mismatch/disturbance dynamic joint monitoring quantityPDMIThe factors for controlling the performance degradation are preliminarily diagnosed if the current working condition is met
Figure BDA00030521101600000518
Greater than UPPDMIIf the model mismatch or the external disturbance dynamics is changed, the system performance is reduced due to the model mismatch or the external disturbance dynamics; if the current operating condition is
Figure BDA00030521101600000519
Less than UPPDMIThen changes in the controller tuning parameters are a factor that causes degradation in system performance.
Step S11 includes the following steps:
s11.1, calculating interference characteristic monitoring quantity eta of online dataDVI：
Figure BDA0003052110160000061
Wherein the content of the first and second substances,
Figure BDA0003052110160000062
representing the interference characteristic monitoring quantity of the online data corresponding to the ith window,
Figure BDA0003052110160000063
for characterizing the change of the white noise variance, ε, in the actual operating conditionsr(t) represents the white noise estimate of the reference condition at the tth sampling time, εh,l(t) represents the white noise estimation value of the actual working condition at the t sampling moment of the ith window in the online monitoring stage, | cov (epsilon)h,l(t)) | denotes εh,lDeterminant of covariance matrix of (t), | cov (εr(t)) | denotes εr(t) determinant of covariance matrix;
s11.2, calculating the actual working condition model residual error of the online data:
Figure BDA0003052110160000064
wherein v ish,l(t) represents the model residual error of the online data actual condition corresponding to the ith window,
Figure BDA0003052110160000065
on-line data actual condition disturbance model estimation value y representing ith windowl(t) an online data process output vector, u, at the tth sampling instant of the ith windowl(t) an online data process input vector, G, for the tth sampling instant of the l-th windowm,l(q-1) Representing a process model used by the EWMA batch controller in the ith window of the online data;
s11.3, calculating disturbance model monitoring quantity eta of online dataDMI：
Figure BDA0003052110160000066
Wherein the content of the first and second substances,
Figure BDA0003052110160000067
representing the monitored quantity, v, of the disturbance model corresponding to the l-th windowr,l(t) represents the on-line data reference condition model residual at the tth sampling instant of the ith window, vh,l(t) represents the on-line data real-time model residual error at the tth sampling time of the ith window, J (f (v)r,l(t)),g(vh,l(t))) represents vh,l(t) and vr,l(t) KL divergence;
s11.4, calculating model mismatch monitoring quantity eta of online dataPMI：
Figure BDA0003052110160000068
Wherein the content of the first and second substances,
Figure BDA0003052110160000069
model mismatch monitoring quantity, epsilon, of the on-line data corresponding to the ith window of the on-line datah,l(t) denotes the ith of the on-line dataWhite noise estimation value of practical working condition at the t-th sampling moment of window vh,l(t) represents the real model residual at the tth sampling instant of the l windows of online data, cov (εh,l(t)) represents εh,l(t) determinant of covariance matrix, cov (v)h,l(t)) represents vh,l(t) determinant of covariance matrix;
s11.5, combining the lower control limit LP of the interference characteristic monitoring quantityDVIUpper limit UP of control of the monitored quantity of interference characteristicsDVILower control limit LP of disturbance model monitoring quantityDMIAnd upper control limit UP of model mismatch monitoring quantityPMIFurther diagnosing and identifying the factors causing the performance degradation of the EWMA batch-to-batch control system, and a plurality of disturbance characteristic monitoring quantities
Figure BDA00030521101600000610
Greater than UPDVIOr less than LPDVIIf so, the white noise variance changes to cause the performance of the EWMA batch-to-batch control system to be reduced; if the model monitoring quantity is disturbed
Figure BDA00030521101600000611
Less than LPDMIIf so, the performance of the EWMA batch control system is reduced due to the change of the external disturbance model; monitoring quantity if model mismatching
Figure BDA00030521101600000612
Greater than UPPMIThen the model mismatch causes the EWMA inter-batch control system performance to degrade.
Compared with the prior art, the invention has the following beneficial effects:
1) when monitoring the control performance of a semiconductor process batch control system, the control performance monitoring quantity constructed by adopting a typical variable analysis method is closely related to controller regulation, model mismatch and disturbance dynamics according to the relation between process output and process input and a process-disturbance combined model residual error; compared with other control performance evaluation indexes, the control performance monitoring quantity provided by the invention directly reflects the influence of controller regulation, model mismatch and disturbance dynamics on the control performance.
2) When the control performance of a semiconductor process batch control system is diagnosed, four diagnosis statistics are established, wherein the four diagnosis statistics are as follows: model mismatch/disturbance dynamic monitoring quantity, interference characteristic monitoring quantity, disturbance model monitoring quantity and model mismatch monitoring quantity; the model mismatch/disturbance dynamic monitoring quantity is not influenced by the adjusting parameters of the controller, the disturbance characteristic monitoring quantity can separate the white noise variance change from other factors influencing the control performance, the disturbance model monitoring quantity can separate the disturbance model change from other factors influencing the control performance, and the model mismatch monitoring quantity separates the model mismatch from other factors influencing the control performance.
3) The method can effectively monitor the control performance change condition of the system and can also effectively distinguish the influence of three factors of model mismatch, disturbance dynamics and controller adjusting parameters on the control performance, when the control performance is reduced, the closed-loop system needs to be shut down, and a system identification experiment is carried out on the closed-loop system to diagnose the reason of the reduction of the control performance of the system.
Drawings
In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings used in the description of the embodiments or the prior art will be briefly described below, it is obvious that the drawings in the following description are only some embodiments of the present invention, and for those skilled in the art, other drawings can be obtained according to the drawings without creative efforts.
FIG. 1 is a flow chart of the present invention;
FIG. 2 is a schematic diagram of a CMP process of the present invention;
FIG. 3 is a schematic diagram of the control performance monitoring index and each diagnostic index when the model is mismatched according to the present invention, wherein (a) is a monitoring graph of the control performance monitoring amount, (b) is a monitoring graph of the model mismatch/disturbance dynamics monitoring amount, (c) is a monitoring graph of the disturbance characteristic monitoring amount, (d) is a monitoring graph of the disturbance model monitoring amount, and (e) is a monitoring graph of the model mismatch monitoring amount;
FIG. 4 is a schematic diagram of the control performance monitoring index and each diagnostic index when the white noise variance changes according to the present invention, wherein (a) is a monitoring graph of the monitoring amount of control performance, (b) is a monitoring graph of the monitoring amount of model mismatch/disturbance dynamics, (c) is a monitoring graph of the monitoring amount of disturbance characteristic, (d) is a monitoring graph of the monitoring amount of disturbance model, and (e) is a monitoring graph of the monitoring amount of model mismatch;
FIG. 5 is a diagram of a control performance monitoring index and diagnostic indices when a controller parameter changes according to the present invention, wherein (a) is a monitoring graph of a control performance monitoring amount, (b) is a monitoring graph of a model mismatch/disturbance dynamics monitoring amount, (c) is a monitoring graph of a disturbance characteristic monitoring amount, (d) is a monitoring graph of a disturbance model monitoring amount, and (e) is a monitoring graph of a model mismatch monitoring amount;
FIG. 6 is a schematic diagram of the control performance monitoring index and each diagnostic index when the disturbance model is changed, wherein (a) is a monitoring graph of the control performance monitoring amount, (b) is a monitoring graph of the model mismatch/disturbance dynamics monitoring amount, (c) is a monitoring graph of the disturbance characteristic monitoring amount, (d) is a monitoring graph of the disturbance model monitoring amount, and (e) is a monitoring graph of the model mismatch monitoring amount, according to the present invention;
FIG. 7 is a schematic diagram of the control performance monitoring index and each diagnostic index when the model is mismatched and the white noise variance is changed, wherein (a) is a monitoring graph of the monitoring quantity of the control performance, (b) is a monitoring graph of the monitoring quantity of the model mismatch/disturbance dynamics, (c) is a monitoring graph of the monitoring quantity of the disturbance characteristic, (d) is a monitoring graph of the monitoring quantity of the disturbance model, and (e) is a monitoring graph of the monitoring quantity of the model mismatch;
fig. 8 is a schematic diagram of the control performance monitoring index and each diagnostic index when the model is mismatched and the disturbance model is changed according to the present invention, where (a) is a monitoring map of the control performance monitoring amount, (b) is a monitoring map of the model mismatch/disturbance dynamics monitoring amount, (c) is a monitoring map of the disturbance characteristic monitoring amount, (d) is a monitoring map of the disturbance model monitoring amount, and (e) is a monitoring map of the model mismatch monitoring amount.
Detailed Description
The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments of the present invention, and it is obvious that the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. All other embodiments, which can be obtained by a person skilled in the art without inventive effort based on the embodiments of the present invention, are within the scope of the present invention.
As shown in fig. 1, an embodiment of the invention provides a method for diagnosing control performance of a semiconductor manufacturing lot-to-lot control system, which includes the following steps:
s1, obtaining N groups of input data and output data of historical normal working conditions of the EWMA batch control system in the semiconductor manufacturing process as sample data, dividing the sample data into S +1 windows, and calculating a white noise estimation value epsilon of the reference working condition by using the input data and the output data of the first windowrAnd reference working condition disturbance model estimation value
Figure BDA0003052110160000071
Obtaining process-disturbance combined model residual e (t), process data set
Figure BDA0003052110160000072
And output data set Y1 0Collecting the process data
Figure BDA0003052110160000073
Is subtracted from the ith columnProcess data set
Figure BDA0003052110160000074
Mean of all samples in column i divided by the process data set
Figure BDA0003052110160000075
Obtaining a matrix X of the reference working condition by the standard deviation of the sample of the ith column1Will output data set Y1 0Is subtracted from the output data set Y1 0Is divided by the output data set Y after averaging all samples in column i1 0Obtaining a matrix Y of the reference working condition by the standard deviation of the sample of the ith column1Construction of Hankel matrix XpAnd XfAnd projection matrix Jq。
Step S1 includes:
s1.1, acquiring input data and output data of N groups of historical normal working conditions in a semiconductor manufacturing process with an EWMA batch controller as sample data, dividing the sample data into S +1 windows, wherein each window comprises N groups of samples, and S is N-N. The calculation process of the window S is as follows: s ═ N/d ], where d is the window moving step length, [ · ] denotes rounding the data down, in this example, d ═ 1.
S1.2, regarding a first window in the S +1 windows acquired in the step S1.1 as a reference working condition, and respectively performing centralization processing on input data and output data of the first window:
Figure BDA0003052110160000076
wherein u' (t) represents the process input vector at the t-th sampling time after the centering process, u (t) represents the input vector at the t-th sampling time of the first window,
Figure BDA0003052110160000077
representing the mean vector, m, of the input vectors at the t-th sampling instant of the first windowuRepresenting the number of input vectors, t representing the sampling instant, n representing the firstThe number of samples in the window, y' (t) represents the process output vector at the t-th sampling instant after the centering process, y (t) represents the output vector at the t-th sampling instant of the first window,
Figure BDA0003052110160000078
mean vector, m, representing the output vector at the t-th sampling instant of the first windowyIndicating the number of output vectors.
S1.3, obtaining an estimated value of white noise of a reference working condition:
εr(t)＝y′P(t){I-YP′(t)T[YP′(t)YP′(t)T]-1YP′(t)}，
wherein epsilonr(t) represents a reference condition white noise estimation value y 'at the t-th sampling time'P(t)＝[y′(t) … y′(t-P+1)]Representing a data vector constructed from a centralised process output vector Y' (t), YP′(t)＝[y′P(t-1)T… y′P(t-L)T]TIs represented by a data vector y'P(t) constructing a data matrix, P representing the size of a data window for estimating white noise, L representing the data matrix YP' (T) and T denotes a transpose operation.
S1.4, constructing an autoregressive moving average model of the process output data after centralization processing and with out-of-band source input:
Figure BDA0003052110160000081
wherein the content of the first and second substances,
Figure BDA0003052110160000082
denotes a delay factor, phi1，…，φICoefficient, phi, representing the autoregressive partI+1，…，φI+KCoefficient of exogenous input, phiI+K+1， …，φI+K+JCoefficient representing moving average, I represents maximum order of autoregressive part, and K represents exogenous inputMaximum order of entry; j denotes the maximum order of the moving average.
S1.5, acquiring coefficient vectors and orders of an autoregressive moving average model with out-of-band source input according to an adaptive minimum absolute value contraction and operator selection method:
Figure BDA0003052110160000083
wherein the content of the first and second substances,
Figure BDA0003052110160000084
YW(t)＝[y′(t) … y′(t-W+1)]T(ii) a W denotes the window size of the sampled data, YW(t) denotes a W × 1-dimensional data matrix, HW(t) represents a W × (I + K + J) dimensional data matrix;
Figure BDA0003052110160000085
expressing the adaptive minimum absolute value shrinkage and selection operator method to estimate the coefficient vector of the autoregressive moving average model, i.e.
Figure BDA0003052110160000086
Denotes the coefficients of the estimated autoregressive moving average model, phi denotes the coefficient vector of the actual but unknown autoregressive moving average type, phi ═ phi1,…,φI+K+J]λ denotes an adjustable parameter controlling the penalty factor in the adaptive minimum absolute value shrinkage and selection operator method,
Figure BDA0003052110160000087
j-th coefficient phi in coefficients representing autoregressive moving average modeljThe weighting factor of (2).
S1.6, calculating a process-disturbance combined model residual error according to an input vector and an output vector of a first window:
e(t)＝y(t)-Gm(q-1)u(t)，
wherein e (t) represents the process-perturbation joint model residual, Gm(q-1) Indicating EWMA inter-batch controller usage in the first windowProcess model of (a), q-1A difference operator is represented.
S1.7, acquiring a process data set of a reference working condition by utilizing a process input vector, a process output vector and a process-disturbance combined model residual error of a first window
Figure BDA0003052110160000088
And output data set Y1 0：
Figure BDA0003052110160000089
Wherein the content of the first and second substances,
Figure BDA00030521101600000810
which represents the input vector of the process,
Figure BDA00030521101600000811
representing the process-perturbation joint model residual vector,
Figure BDA00030521101600000812
representing the process output vector.
S1.8, calculating a mean vector and a standard deviation vector of process input vectors of a first window:
Figure BDA00030521101600000813
Figure BDA00030521101600000814
wherein the content of the first and second substances,
Figure BDA00030521101600000815
is the mean vector, σ, of the process input vectors at the t-th sampling instant of the first windowu,rThe process of the first window inputs the vector of standard deviation of the vector;
calculating a mean vector and a standard deviation vector of the output vector of the first window process:
Figure BDA00030521101600000816
Figure BDA00030521101600000817
wherein the content of the first and second substances,
Figure BDA0003052110160000091
is the mean vector, σ, of the process output vectors at the t-th sampling instant of the first windowy,rThe process of the first window outputs a vector of standard deviations of the vectors;
calculating a mean vector and a standard deviation vector of a first window process-disturbance combined model residual:
Figure BDA0003052110160000092
Figure BDA0003052110160000093
wherein the content of the first and second substances,
Figure BDA0003052110160000094
mean vector, σ, of the process-perturbation combined model residuals for the first windowe,rThe process-perturbation for the first window is the standard deviation vector of the joint model residual.
S1.9, constructing a matrix X1And Y1：
Figure BDA0003052110160000095
Wherein u isiRepresents the ith process input vector, i 1, …, mu，muDenotes the number of input vectors, t denotes the sampling time, ui"represents the i-th process input vector after normalization, i.e.
Figure BDA0003052110160000096
Figure BDA0003052110160000097
Means representing the ith process input vector, (σ)u,r)iRepresents the standard deviation of the ith process input vector; y isjRepresents the jth process output vector, j ═ 1, …, my，myIndicating the number of output vectors, yj"represents the normalized j-th process output vector, i.e.
Figure BDA0003052110160000098
Figure BDA0003052110160000099
Means for representing the jth process output vector, (σ)y,r)jRepresents the standard deviation of the jth process output vector; e.g. of the typej(t) denotes the jth process-perturbation joint model residual vector, where ej＝(ej(t))， ej(t)＝yj(t)-(Gm(q-1))ju(t)，(Gm(q-1))jRepresenting the process model G used in the first windowm(q-1) J-th line of (1), yj(t) denotes the jth process output at the tth sampling instant, u (t) denotes the column vector formed by all process inputs at the tth sampling instant in the first window, e ″jRepresents the normalized j process-perturbation joint model residual vector,
Figure BDA00030521101600000910
Figure BDA00030521101600000911
means for representing the jth process-perturbation combined model residual, (σ)e,r)jRepresents the standard deviation of the jth process-perturbation combined model residual.
S1.10, constructionHankel matrix XpAnd Xf：
Figure BDA00030521101600000912
Wherein m is mu+2my，pp+1Representing the historical data vector at the p +1 th sampling instant, i.e. pp+1＝[x(p)T … x(1)T y(p)T … y(1)T]TX (p) represents a matrix
Figure BDA00030521101600000913
P column of (2), x (1) represents a matrix
Figure BDA00030521101600000914
The number 1. sup. st column of (a),
Figure BDA00030521101600000915
representing matrix X1Y (p) denotes a matrix Y1 TY (1) denotes the matrix Y1 TColumn 1, Y1 TRepresentation matrix Y1Transposing; f. ofp+1The vector of prediction data representing the p +1 th sampling instant, i.e. fp+1＝[y(p+1)T y(p+2)T … y(p+f-1)T]TY (p +1) represents Y1 TColumn p +1, Y (p +2) represents Y1 TIn the p +2 th column, Y (p + f-1) represents Y1 TM denotes the number of vectors of the history data vector and the prediction data vector constructed from n sets of observation data, and M is n-p-f + 1.
S1.11, calculating a unitary matrix V: h ═ U ∑ VTWherein V represents a unitary matrix constructed when singular value decomposition is performed on the matrix H,
Figure BDA00030521101600000916
is represented by XpAnd XfThe constructed matrix Hankel matrix is used as a matrix,
Figure BDA00030521101600000917
represents the sample covariance of the historical data vector,
Figure BDA00030521101600000918
represents the sample covariance of the prediction data vector,
Figure BDA00030521101600000919
representing the cross covariance of the historical data vector and the predicted data vector.
S1.12, calculating a projection matrix Jq：
Figure BDA00030521101600000920
Wherein, VqRepresenting a matrix constructed from the first q columns of the unitary matrix V,
Figure BDA00030521101600000921
represents the Hankel matrix XpThe covariance matrix of (2).
S2, constructing a process data set using the input data and the output data of the k' th window
Figure BDA00030521101600000922
And outputting the data set
Figure BDA00030521101600000923
Collecting process data
Figure BDA00030521101600000924
Column i minus the process data set
Figure BDA0003052110160000101
Mean of all samples in column i divided by the process data set
Figure BDA0003052110160000102
Obtaining a matrix X of actual working conditions by the standard deviation of the sample of the ith columnk'Will output the data set
Figure BDA0003052110160000103
I th of (1)Column-minus output dataset
Figure BDA0003052110160000104
Is divided by the output data set after averaging all samples in column i
Figure BDA0003052110160000105
Obtaining a matrix Y of actual working conditions by the standard deviation of the sample of the ith columnk'Acquiring control performance monitoring quantity of an EWMA batch-to-batch control system by adopting a typical variable Analysis (CVA)
Figure BDA0003052110160000106
Where k' is 2,3, …, S + 1.
Step S2 includes:
s2.1, obtaining a process-disturbance combined model residual error of a k 'th window by using input data and output data of the k' th window: e.g. of the typek'(t)＝yk'(t)-Gm,k'(q-1)uk' (t) wherein ek'(t) denotes the process-perturbation joint model residual at the t-th sampling instant of the k' -th window, yk'(t) denotes the kth' window process output vector at the tth sampling instant uk'(t) represents the process input vector for the kth window at the tth sampling instant, Gm,k'(q-1) Representing the process model used by the EWMA batch controller in the current operating regime.
S2.2, acquiring a process input data set under the current working condition according to the process input vector and the output vector of the kth window and the process-disturbance combined model residual error of the kth window
Figure BDA0003052110160000107
And outputting the data set
Figure BDA0003052110160000108
Figure BDA0003052110160000109
Wherein the content of the first and second substances,
Figure BDA00030521101600001010
the process input vector representing the kth window,
Figure BDA00030521101600001011
the process-perturbation joint model residual vector representing the kth window,
Figure BDA00030521101600001012
the process output vector representing the k' th window.
S2.3, constructing a matrix Xk'And Yk'：
Figure BDA00030521101600001013
Wherein u isi,k'The process input vector, u ", representing the k' th windowi,k'The process input vector representing the normalized k' th window, i.e.
Figure BDA00030521101600001014
yj,k'Represents the process output vector, y ″, of the k' th windowj,k'The process output vector representing the normalized k' th window, i.e.
Figure BDA00030521101600001015
ej,k'Process-perturbation joint model residual vector representing the kth window, ej,k'＝(ej,k'(t))， ej,k'(t)＝yj,k'(t)-(Gm,k′(q-1))juk'(t)，yj,k'(t) denotes the jth process output, u, at the tth sampling instant in the kth windowk'(t) represents the column vector formed by all process inputs at the tth sampling instant in the kth window, (G)m,k′(q-1))jRepresents Gm，k′(q-1) Line j, e ″)j,k'Representing the process-perturbation combined model residual of the normalized k' th window, i.e.
Figure BDA00030521101600001016
S2.4, constructing a Hankel matrix Xp,k'And Xf,k'：
Figure BDA00030521101600001017
Wherein p isp+1kRepresenting the historical data vector at the p +1 th sampling instant in the k' th window, i.e. pp+1,k'＝[xk'(p)T … xk'(1)T yk'(p)T … yk'(1)T]T，xk'(p) a representation matrix
Figure BDA00030521101600001018
P column of (2), xk'(1) Representation matrix
Figure BDA00030521101600001019
The number 1. sup. st column of (a),
Figure BDA00030521101600001020
representation matrix Xk'Transpose of (y)k'(p) a representation matrix
Figure BDA00030521101600001021
P column, yk'(1) Representation matrix
Figure BDA00030521101600001022
The number 1. sup. st column of (a),
Figure BDA00030521101600001023
representation matrix Xk'Transposing; f. ofp+1kDenotes the vector of prediction data at the p +1 th sampling instant in the k' th window, i.e. fp+1,k'＝[y(k')(p+1)T y(k')(p+2)T … y(k')(p+f-1)T]T，y(k')(p +1) represents
Figure BDA00030521101600001024
P +1 column, y(k')(p +2) represents
Figure BDA00030521101600001025
P +2 column, y(k')(p + f-1) represents
Figure BDA00030521101600001026
Column p + f-1.
S2.5, obtaining control performance monitoring quantity of the EWMA batch control system by adopting a typical variable analysis method
Figure BDA00030521101600001027
Comprises the following steps:
Figure BDA00030521101600001028
wherein the content of the first and second substances,
Figure BDA00030521101600001029
represents the amount of control performance monitoring, J, of the EWMA batch-to-batch control system corresponding to the kth windowqRepresenting a projection matrix, pk,k'Representation matrix Xp,k'The k-th column of (1).
S3, according to the reference working condition disturbance model estimated value in the step S1
Figure BDA00030521101600001030
Obtaining reference working condition model residual error v from input data and output data of k' th windowr,k'(t) and model mismatch/perturbation dynamics joint monitoring volume
Figure BDA00030521101600001031
Step S3 includes:
s3.1, disturbing the model estimated value according to the reference working condition
Figure BDA0003052110160000111
Obtaining reference working condition model residual error v from input data and output data of k' th windowr,k'(t)：
Figure BDA0003052110160000112
Wherein v isr,k'(t) reference behavior model residual, y, at the t-th sampling time of the k' -th windowk'(t) represents the process output vector at the t-th sampling instant of the k' -th window, uk'(t) represents the process input vector at the t-th sampling instant of the k' -th window, Gm,k'(q-1) Representing the process model used by the EWMA batch controller in the k' th window.
S3.2, according to the model residual error v of the reference working conditionr,k'(t) obtaining model mismatch/disturbance dynamic combined monitoring quantity by white noise estimated value of reference working condition
Figure BDA0003052110160000113
Wherein the content of the first and second substances,
Figure BDA0003052110160000114
represents the model mismatch/disturbance dynamic joint monitoring quantity corresponding to the k' th window,
Figure BDA0003052110160000115
has a value range of (0, 1)]，εr(t) the white noise estimation value of the reference condition at the t-th sampling moment, vr,k'(t) represents the reference condition model residual at the t-th sampling instant in the kth' window, | cov (εr(t)) | denotes εrDeterminant of covariance matrix of (t), | cov (v)r,k'(t)) | denotes vr,k'(t) determinant of covariance matrix.
S4, obtaining the white noise estimation value epsilon of the practical working condition corresponding to the k 'th window by using the input data and the output data of the k' th windowh,k'(t) and actual condition disturbance model estimation value
Figure BDA0003052110160000116
And according to the actual working condition disturbance model estimated value
Figure BDA0003052110160000117
Obtaining an actual condition modelResidual vh,k'(t)。
Step S4 includes:
s4.1, carrying out centralized processing on input data and output data of a k' th window:
Figure BDA0003052110160000118
wherein u'k'(t) represents the process input vector at the t-th sampling instant of the k' -th window after the centering process, uk'(t) represents the process input vector at the t-th sampling instant of the k' -th window that is not centered,
Figure BDA0003052110160000119
mean vector, y ', of process input vectors representing the t-th sampling instant of the k ' -th window 'k'(t) represents the process output vector at the t-th sampling instant of the k' -th window after the centering process, yk'(t) represents the process output vector at the t-th sampling instant of the k' -th window that is not centered,
Figure BDA00030521101600001110
the mean vector of the process output vectors at the t-th sampling instant of the k' -th window is represented.
S4.2, calculating an estimated value of white noise under an actual working condition: epsilonh,k'(t)＝y′P,k'(t){I-Y′P,k'(t)T[Y′P,k'(t)Y′P,k'(t)T]-1Y′P,k'(t) }, in which εh,k'(t) represents an actual working condition white noise estimated value y 'of the acquired kth sampling moment of the kth window'P,k'(t)＝[y′k'(t) … y′k'(t-P+1)]Represents the process output vector y ' using the centered k ' th window 'k'(t) constructed data vector, Y'P,k'(t)＝[y′P,k'(t-1)T … y′P,k'(t-L)T]TIs represented by a data vector y'P,k'(t) the constructed data matrix.
S4.3, constructing an autoregressive moving average model of the out-of-band source input of the process output vector of the kth window:
Figure BDA00030521101600001111
wherein the content of the first and second substances,
Figure BDA00030521101600001112
which is indicative of the delay factor, is,
Figure BDA00030521101600001113
the coefficients of the auto-regressions are represented,
Figure BDA00030521101600001114
the coefficient representing the input of the external source,
Figure BDA00030521101600001115
coefficient representing moving average, Ik' represents the maximum order of autoregressive, Kk' represents the maximum order of exogenous input; j. the design is a squarek' denotes the maximum order of the moving average.
S4.4, acquiring coefficient vector and order of the autoregressive moving average model with out-of-band source input according to the adaptive minimum absolute value contraction and operator selection method:
Figure BDA00030521101600001116
wherein the content of the first and second substances,
Figure BDA00030521101600001117
Figure BDA0003052110160000121
Wk' denotes the window size of the sampled data,
Figure BDA0003052110160000122
represents WkA data matrix of dimension' x 1,
Figure BDA0003052110160000123
represents Wk'×(Ik'+Kk'+Jk') A dimensional data matrix;
Figure BDA0003052110160000124
the method representing the adaptive minimum absolute value shrinkage and selection operator estimates the coefficient vector of the autoregressive moving average model, i.e.
Figure BDA0003052110160000125
Coefficient, phi, representing the estimated autoregressive moving average modelk'A coefficient vector representing a real but unknown type of auto-regressive moving average,
Figure BDA0003052110160000126
Figure BDA0003052110160000127
j-th coefficient phi in coefficients representing autoregressive moving average modelj,k'The weighting factor of (2).
S4.5, obtaining an estimated value of the actual working condition disturbance model according to the autoregressive moving average model with out-of-band source input:
Figure BDA0003052110160000128
wherein the content of the first and second substances,
Figure BDA0003052110160000129
and representing the estimated value of the actual working condition disturbance model of the k' th window.
S4.6, calculating the residual error of the actual working condition model:
Figure BDA00030521101600001210
wherein v ish,k'(t) actual condition model of t sampling time of k' th windowThe residual error is a residual error that is,
Figure BDA00030521101600001211
actual condition disturbance model estimation value representing k' th window
Figure BDA00030521101600001212
Inverse of (y)k'(t) represents the process output vector at the t-th sampling instant of the k' -th window, uk'(t) represents the process input vector at the t-th sampling instant of the k' -th window, Gm,k′(q-1) Representing the process model used by the EWMA batch controller in the k' th window.
S5, according to the reference working condition white noise estimated value epsilon in the step S1rAnd the white noise estimate of the actual condition ε in step S4h,k'(t) calculating an interference characteristic monitoring quantity
Figure BDA00030521101600001213
According to the reference working condition model residual v in the step S3r,k'(t) and the actual condition model residual v in step S4h,k'(t) calculating the amount of monitoring of the disturbance model
Figure BDA00030521101600001214
According to the actual working condition white noise estimated value epsilon in the step S4h,k'(t) and the actual condition model residual vh,k'(t) obtaining a model mismatch monitoring quantity
Figure BDA00030521101600001215
Step S5 includes:
s5.1, calculating interference characteristic monitoring quantity
Figure BDA00030521101600001216
Wherein epsilonr(t) represents the white noise estimation value of the reference condition at the t-th sampling moment, epsilonh,k'(t) represents the white noise estimate for the actual condition at the t-th sampling instant in the kth' window, | cov (εh,k'(t)) | denotes εh,k'Determinant of covariance matrix of (t), | cov (εr(t)) | denotes εr(t) determinant of covariance matrix.
S5.2, calculating the monitoring quantity of the disturbance model
Figure BDA00030521101600001217
Wherein the content of the first and second substances,
Figure BDA00030521101600001218
represents the perturbation model monitoring quantity of the k' th window,
Figure BDA00030521101600001219
has a value range of (0, 1)]If, if
Figure BDA00030521101600001220
Near or equal to 1, the perturbation model is unchanged, vr,k'(t) represents the reference operating condition model residual at the t-th sampling moment in the k' -th window, vh,k'(t) represents the actual condition model residual at the t-th sampling instant in the k' -th window, J (f (v)r,k'(t)),g(vh,k'(t))) represents vh,k'(t) and vr,k'(t) KL divergence, expressed as:
Figure BDA00030521101600001221
wherein, f (v)r,k'(t)) represents v constructed using the nuclear density estimation methodr,k'(t) probability density function, f (v)h,k'(t)) represents v acquired by the kernel density estimation methodh,k'(t) probability density function, i ═ 1,2, …, Nf，NfRepresenting the total number of independent identically distributed sample points produced by the kernel density estimate.
S5.3, calculating model mismatch monitoring quantity
Figure BDA00030521101600001222
Wherein epsilonh,k'(t) represents the white noise estimation value of the actual working condition at the t sampling moment in the k' th window, vh,k'(t) represents the t sample in the k' windowActual condition model residual error of moment, | cov (εh,k'(t)) | denotes εh,k'Determinant of covariance matrix of (t), | cov (v)h,k'(t)) | denotes vh,k'(t) determinant of covariance matrix.
S6, respectively obtaining the control upper limit UP of the control performance monitoring quantity by adopting a kernel density estimation methodCPIUpper control limit UP of model mismatch/disturbance dynamic joint monitoring quantityPDMIUpper limit UP of control of the monitored quantity of interference characteristicsDVILower limit of control LP of disturbance characteristic monitoring amountDVILower control limit LP of disturbance model monitoring quantityDMIAnd the upper control limit UP of the model mismatch monitoring quantityPMI。
Step S6 includes:
s6.1, calculating the control upper limit UP of the control performance monitoring quantityCPI：
Figure BDA0003052110160000131
Wherein the content of the first and second substances,
Figure BDA0003052110160000132
is the control performance monitoring amount corresponding to the kth window, K' is 2, 3.., S +1, h represents the kernel bandwidth, and K (·) represents the kernel density estimation with the kernel being a radial basis function.
S6.2, calculating the control upper limit UP of the model mismatch/disturbance dynamic combined monitoring quantityPDMI：
Figure BDA0003052110160000133
Wherein the content of the first and second substances,
Figure BDA0003052110160000134
is the model mismatch/disturbance dynamic joint monitoring quantity corresponding to the k' th window.
S6.3, calculating the control upper limit UP of the interference characteristic monitoring quantityDVIAnd a lower control limit LPDVI：
Figure BDA0003052110160000135
Wherein the content of the first and second substances,
Figure BDA0003052110160000136
is the interference characteristic monitoring quantity corresponding to the k' th window.
S6.4, calculating a control lower limit LP of the monitoring quantity of the disturbance modelDMI：
Figure BDA0003052110160000137
Wherein the content of the first and second substances,
Figure BDA0003052110160000138
is the disturbance model monitoring quantity corresponding to the k' th window.
S6.5, calculating the control upper limit UP of the model mismatch monitoring quantityPMI：
Figure BDA0003052110160000139
Wherein the content of the first and second substances,
Figure BDA00030521101600001310
is the interference characteristic monitoring quantity corresponding to the k' th window.
S7, collecting N1 groups of input data and output data of the EWMA batch control system of the semiconductor process on line, and obtaining the control performance monitoring quantity eta of the current working conditionCPIAccording to the control upper limit UP in step S6CPIAnd judging whether the control performance of the closed-loop system is reduced or not in real time.
Step S7 includes:
s7.1, collecting N1 groups of input data and output data of the EWMA batch control system of the semiconductor process on line as actual working condition sample data, and dividing the actual working condition sample data into S1An online window, each window containing n groups of samples, wherein S1＝N1-n+1。
S7.2, calculating the residual error of the online process-disturbance combined model: e.g. of the typel(t)＝yl(t)-Gm,l(q-1)ul(t) wherein el(t) represents the on-line process-perturbation joint model residual error at the tth sampling instant of the ith window, yl(t) represents the process output vector at the tth sampling instant of the ith window, ul(t) denotes an input vector at the tth sampling time of the ith window, and t denotesSampling time, Gm,l(q-1) Representing a Process model, q, used by an EWMA batch-to-batch controller in an on-line regime-1Denotes a difference operator, where l is 1,2, …, S1。
S7.3, acquiring a process input data set of the current working condition by using the process output vector, the process input vector and the process-disturbance combined model residual error which are acquired in real time
Figure BDA00030521101600001311
And output data set Yl 0：
Figure BDA00030521101600001312
Wherein the content of the first and second substances,
Figure BDA00030521101600001313
representing the on-line process input vector and,
Figure BDA00030521101600001314
representing an online process-perturbation joint model residual vector,
Figure BDA00030521101600001315
representing an online process output vector.
S7.4, constructing a matrix XlAnd Yl：
Figure BDA00030521101600001316
Wherein u isi,lRepresents the process input vector, u ″, in the ith window of the online datai,lThe process input vector representing the normalized l-th window, i.e.
Figure BDA00030521101600001317
yj,lRepresents the process output vector, y ″, of the ith window of online dataj,lProcess output representing normalized ith windowAmount of, i.e.
Figure BDA00030521101600001318
ej,lProcess-perturbation joint model residual vector representing the ith window of online data, ej,l＝(ej,l(t))， ej,l(t)＝yj,l(t)-(Gm,l(q-1))jul(t)，(Gm,l(q-1))jProcess model G representing EWMA inter-batch controller usage in the ith window of online datam，l(q-1) J-th line of (1), yj,l(t) denotes the jth process output at the tth sampling instant of the ith window of online data, ul(t) represents the vector formed by all process inputs at the tth sampling time of the ith window of online data, e ″j,lRepresenting the normalized process-perturbation joint model residual,
Figure BDA0003052110160000141
s7.5, constructing a Hankel matrix Xp,lAnd Xf,l：
Figure BDA0003052110160000142
Wherein p isp+1,lRepresenting a history data vector at the p +1 th sampling instant in the ith window of line data, i.e. pp+1,l＝[xl(p)T … xl(1)T yl(p)T … yl(1)T]T，xl(p) a representation matrix
Figure BDA0003052110160000143
P column of (2), xl(1) Representation matrix
Figure BDA0003052110160000144
The p-th column of (2),
Figure BDA0003052110160000145
representation matrix XlIs rotatedPosition yl(p) represents a matrix Yl TP column, yl(1) Representation matrix Yl TColumn 1, Yl TRepresentation matrix YlTransposing; f. ofp+1,lThe vector of predicted data representing the p +1 th sampling instant in the l window of line data, i.e. fp+1,l＝[yl(p+1)T yl(p+2)T … yl(p+f-1)T]T，yl(p +1) represents Yl TP +1 column, yl(p +2) represents Yl TP +2 column, yl(p + f-1) represents Yl TColumn p + f-1.
S7.6, acquiring control performance monitoring quantity eta of the current working condition by using online dataCPI：
Figure BDA0003052110160000146
Wherein the content of the first and second substances,
Figure BDA0003052110160000147
control performance monitor quantity, p, representing the current operating condition corresponding to the ith windowk,lRepresentation matrix Xp,lThe kth column vector of (1).
S7.7, control Upper Limit UP according to control Performance monitoring amountCPIJudgment of
Figure BDA0003052110160000148
Whether or not less than UPCPIIf, if
Figure BDA0003052110160000149
Less than UPCPIThen the system performance is better, if
Figure BDA00030521101600001410
Greater than UPCPIThe system performance is degraded.
S8, if the control performance monitoring quantity etaCPIControl upper limit UP greater than control performance monitoring amountCPIIf the performance of the EWMA inter-batch control system is degraded, the step S9 is executed, otherwise, the EWMA inter-batch control system is executedThe performance is good.
S9, using the reference working condition disturbance model estimated value in the step S1
Figure BDA00030521101600001411
Obtaining the model mismatch/disturbance dynamic combined monitoring quantity eta of the current working condition from the N1 groups of input data and output data in the step S7PDMI。
Step S9 includes:
s9.1, obtaining a reference working condition model residual error by using online data:
Figure BDA00030521101600001412
wherein v isr,l(t) represents the reference operating mode model residual error at the tth sampling moment of the ith window of the online data,
Figure BDA00030521101600001413
and representing the reference working condition disturbance model estimation value.
S9.2, calculating model mismatch/disturbance dynamic combined monitoring quantity etaPDMI：
Figure BDA00030521101600001414
Wherein the content of the first and second substances,
Figure BDA00030521101600001415
has a value range of (0, 1)]，εr(t) represents the white noise estimate for the baseline condition, | cov (εr(t)) | denotes εr(t) determinant of covariance matrix, vr,l(t) represents the reference condition model residual at the tth sampling instant, | cov (v)r,l(t)) | denotes vr,l(t) determinant of covariance matrix.
S9.3, combining control upper limit UP of model mismatch/disturbance dynamic joint monitoring quantityPDMIThe factors for controlling the performance degradation are preliminarily diagnosed if the current working condition is met
Figure BDA00030521101600001416
Greater than UPPDMIThe model mismatch or external disturbance dynamics changes, andresulting in system performance degradation; if the current operating condition is
Figure BDA00030521101600001417
Less than UPPDMIThen changes in the controller parameters are a factor that causes degradation in system performance.
S10, if the model mismatching/disturbance dynamic joint monitoring quantity etaPDMIGreater than the upper control limit UP of the combined monitoring quantity of model mismatch/disturbance dynamics in step S6PDMIIf the model mismatch or disturbance dynamics results in performance degradation of the EWMA batch control system, step S11 is executed, otherwise, the controller parameter change results in performance degradation of the EWMA batch control system.
S11, calculating the disturbance characteristic monitoring quantity eta of the current working conditionDVIMonitoring quantity eta of disturbance modelDMIAnd model mismatch monitor ηPMI(ii) a Disturbance characteristic monitoring quantity eta when current working conditionDVIGreater than the upper control limit UP of the disturbance characteristic monitoring amount in step S6DVIOr less than the lower control limit LP of the disturbance characteristic monitored quantityDVIWhen the interference characteristic changes, the performance of the EWMA batch control system is reduced; when the disturbance model monitors the quantity etaDMILower control limit LP less than disturbance model monitoring quantityDMIMeanwhile, the performance of the EWMA batch control system is reduced due to the change of the disturbance model; when model mismatch monitor ηPMIControl upper limit UP greater than model mismatch monitoring amountPMIModel mismatch results in performance degradation of the EWMA batch-to-batch control system.
Step S11 includes:
s11.1, calculating interference characteristic monitoring quantity eta of online dataDVI：
Figure BDA0003052110160000151
Wherein the content of the first and second substances,
Figure BDA0003052110160000152
representing the interference characteristic monitoring quantity of the on-line data corresponding to the ith window, representing the change of white noise variance in actual working conditionr(t) the white noise estimation value of the reference condition at the tth sampling time，εh,l(t) represents the white noise estimate for the actual condition at the tth sampling instant, | cov (εh,l(t)) | denotes εh,l(t) covariance matrix determinant, | cov (εr(t)) | denotes εr(t) covariance matrix determinant.
S11.2, calculating the actual working condition model residual error of the online data:
Figure BDA0003052110160000153
wherein v ish,l(t) represents the model residual error of the online data actual condition corresponding to the ith window,
Figure BDA0003052110160000154
on-line data actual condition disturbance model estimation value y representing ith windowl(t) an online data process output vector, u, at the tth sampling instant of the ith windowl(t) an online data process input vector, G, for the tth sampling instant of the l-th windowm,l(q-1) Representing a process model used by the EWMA batch controller in the ith window of the online data;
s11.3, calculating disturbance model monitoring quantity eta of online dataDMI：
Figure BDA0003052110160000155
Wherein the content of the first and second substances,
Figure BDA0003052110160000156
representing the monitored quantity, v, of the disturbance model corresponding to the l-th windowr,l(t) represents the on-line data reference condition model residual at the tth sampling instant of the ith window, vh,l(t) represents the on-line data real-time model residual error at the tth sampling time of the ith window, J (f (v)r,l(t)),g(vh,l(t))) represents vh,l(t) and vr,l(t) KL divergence.
S11.4, calculating model mismatch monitoring quantity eta of online dataPMI：
Figure BDA0003052110160000157
Wherein the content of the first and second substances,
Figure BDA0003052110160000158
model mismatch monitoring quantity, epsilon, of the on-line data corresponding to the ith window of the on-line datah,l(t) an actual condition white noise estimate, v, at the tth sampling time of the ith window of the on-line datah,l(t) represents the real model residual at the tth sampling instant of the l windows of online data, cov (εh,l(t)) represents εh,l(t) determinant of covariance matrix, cov (v)h,l(t)) represents vh,l(t) determinant of covariance matrix.
S11.5, combining the lower control limit LP of the interference characteristic monitoring quantityDVIUpper limit UP of control of the monitored quantity of interference characteristicsDVILower control limit LP of disturbance model monitoring quantityDMIAnd upper control limit UP of model mismatch monitoring quantityPMIFurther diagnosing and identifying the factors causing the performance degradation of the EWMA batch-to-batch control system if the characteristic monitoring quantity is disturbed
Figure BDA0003052110160000159
Greater than UPDVIOr less than LPDVIIf so, the white noise variance changes to cause the performance of the EWMA batch-to-batch control system to be reduced; if the model monitoring quantity is disturbed
Figure BDA00030521101600001510
Less than LPDMIIf so, the performance of the EWMA batch control system is reduced due to the change of the external disturbance model; monitoring quantity if model mismatching
Figure BDA00030521101600001511
Greater than UPPMIThen the model mismatch causes the EWMA inter-batch control system performance to degrade.
Specific example 1
Chemical-mechanical planarization (CMP) is the planarization method of choice in the semiconductor industry today, as shown in figure 2. In 1989, Davari et al developed a CMP process for global planarization of interlevel dielectrics. In the early 90 s of the 20 th century, IBM successfully addressed the problem of slow dry etching of copper using the CMP process. In addition, CMP has been integrated with the dual-damascen process to simplify the fabrication steps of conductive lines between successive layers of a multi-layer integrated circuit. As wafer sizes increase to 300 mm or more, equipment and process design improvements address these issues. In addition, as the demand for high-performance products is increasing in the fierce semiconductor market, the quality of semiconductor products is also more and more important as manufacturers increase. This example uses a simulated CMP process. The discrete model of the CMP process is described as follows:
Figure BDA00030521101600001512
where θ is 0.65, and ∈ (t) ∈ N (0,100). According to the above formula, the process has a manipulated variable and a controlled variable, wherein the manipulated variable u (t) is the platen rotation speed, and the controlled variable y (t) is the removal rate. The removal rate y (t) represents the average of the difference between the initial and final phosphosilicate glass thickness on the wafer. An EWMA inter-batch controller is used to control the removal rate during CMP. The EWMA batch controller consists of an EWMA filter and a minimum beat control rule, and is as follows:
Figure BDA00030521101600001513
wherein, a (t) represents the output of the EWMA filter, t represents the sampling time, lambda represents the adjusting parameter of the EWMA batch controller, b represents the coefficient of the process model in the EWMA batch controller, tau represents the set value of the controlled variable y (t), and u (t-1) represents the manipulated variable at the t-1 sampling time.
The control performance diagnosis method of the semiconductor process batch control system based on the model residual error provided by the invention is used for diagnosing the control performance of the CMP process, and comprises the following steps:
(1) the method comprises the steps of collecting 1000 groups of input data and output data of historical normal working conditions in a semiconductor manufacturing process of an EWMA batch controller as standard working condition sample data, and obtaining process input data and process output data after centralized processing.
The reference working condition sample data is divided into S +1 windows, and each window comprises 500 groups of samples, wherein S is 500. Firstly, the mean value of the input data of the process is obtained, and the calculation formula is
Figure BDA0003052110160000161
Then, 500 sets of process input data after centralization are obtained:
Figure BDA0003052110160000162
according to 500 sets of collected process output data, firstly, obtaining a mean value of process output, wherein a calculation formula is as follows:
Figure BDA0003052110160000163
then, 500 sets of process output data after centralization are obtained, and the calculation formula is as follows:
Figure BDA0003052110160000164
(2) acquiring an estimated value of white noise of a reference working condition;
firstly, outputting a data vector constructed by y' (t) according to the centralization processing process acquired in the step (1): y'P(t)＝[y′(t) … y′(t-P+1)](ii) a Then, by vector y'P(t) constructed data matrix: y isP′(t)＝[y′P(t-1)T … y′P(t-L)T]TWhere P represents the data window size of the estimated white noise. Obtaining an estimated value of the white noise of the reference working condition according to the following formula: epsilonr(t)＝y′P(t){I-YP′(t)T[YP′(t)YP′(t)T]-1YP′(t)}。
(3) Obtaining a reference working condition disturbance model estimation value;
(3.1) obtaining an autoregressive moving average model of the out-of-band source input of the centralized process output according to the following formula:
Figure BDA0003052110160000165
wherein the content of the first and second substances,
Figure BDA0003052110160000166
denotes a delay factor, phi1，…，φICoefficient of autoregressive representationI+1，…，φI+KCoefficient of exogenous input, phiI+K+1， …，φI+K+JRepresenting the coefficient of the moving average, I representing the maximum order of autoregressive, and K representing the maximum order of exogenous input; j represents the maximum order of the moving average, y '(t) represents the process output after the centering process obtained in step (1), u' (t) represents the process input after the centering process obtained in step (1), εrAnd (t) representing the reference working condition white noise estimated value obtained in the step (2).
(3.2) according to the self-adaptive minimum absolute value contraction and operator selection method, obtaining the coefficient vector and the order of the autoregressive moving average model with out-of-band source input:
Figure BDA0003052110160000167
wherein the content of the first and second substances,
Figure BDA0003052110160000168
YW(t)＝[y′(t) … y′(t-W+1)]T(ii) a W denotes the window size of the sampled data, YW(t) denotes a W × 1-dimensional data matrix, HW(t) represents a W × (I + K + J) dimensional data matrix; y '(t) represents the process output after the centering process, u' (t) represents the process input after the centering process, εr(t) represents the acquired white noise estimate of the reference condition,
Figure BDA0003052110160000169
the coefficient vector representing the adaptive minimum absolute value shrinkage and selection operator method for estimating the autoregressive moving average model, i.e.
Figure BDA00030521101600001610
Figure BDA00030521101600001611
Representing the coefficients of the estimated autoregressive moving average model, I representing the maximum order of the autoregressive portion, and K representing the maximum order of the exogenous input portion; j denotes the maximum order of the moving average part, phi denotes the coefficient vector of the actual but unknown autoregressive moving average type, phi ═ phi1,…,φI+K+J]λ represents an adjustable parameter for controlling penalty factors in the adaptive minimum absolute value shrinkage and selection operator method,
Figure BDA00030521101600001612
j-th coefficient phi in coefficients representing autoregressive moving average modeljThe weighting factor of (2).
(3.3) according to the autoregressive moving average model of the out-of-band source input estimated in the step (3.2), acquiring a reference working condition disturbance model estimated value according to the following formula:
Figure BDA00030521101600001613
wherein the content of the first and second substances,
Figure BDA00030521101600001614
coefficients representing the above estimated autoregressive portion and moving average portion, coefficients of an autoregressive moving average model,
Figure BDA00030521101600001615
and representing the reference working condition disturbance model estimation value.
(4) Obtaining process-disturbance combined model residual
Obtaining a process-perturbation combined model residual error of the first window according to the following formula: e (t) ═ y (t) — Gm(q-1) u (t); wherein e (t) represents the process-perturbation combined model residual error, y (t) represents the process output at the tth sampling moment, u (t) represents the process input at the tth sampling moment, Gm(q-1) Representing the process model used by the EWMA batch controller in the first window, q-1A difference operator is represented.
(5) Construction of Hankel matrix XpAnd Xf
(5.1) acquiring a process data set of the reference working condition according to the following formula by using the process output, the process input and the process-disturbance combined model residual error of the first window
Figure BDA0003052110160000171
And output data set Y1 0：
Figure BDA0003052110160000173
Wherein the content of the first and second substances,
Figure BDA0003052110160000174
representing the process input vector, muDenotes the number of input variables, mu＝1，
Figure BDA0003052110160000175
Representing the process-perturbation joint model residual vector,
Figure BDA0003052110160000176
representing the process output vector, myDenotes the number of output variables, my＝1。
(5.2) obtaining the process input mean vector and the standard deviation vector of the first window according to the following formula:
Figure BDA0003052110160000177
Figure BDA0003052110160000178
wherein the content of the first and second substances,
Figure BDA0003052110160000179
is the over-sampling of the tth sampling time of the first windowMean vector of range inputs, σu,rThe standard deviation vector of the process input, m, for the first windowuDenotes the number of input variables, mu＝1，myDenotes the number of input variables, myT denotes the sampling time, n denotes the number of samples of the first window, n is 500, u1(t) represents a first input variable,
Figure BDA00030521101600001710
denotes the m-thuAn input variable; obtaining a mean vector and a standard deviation vector output by the first window process according to the following formula:
Figure BDA00030521101600001711
Figure BDA00030521101600001712
wherein the content of the first and second substances,
Figure BDA00030521101600001713
is the mean vector, σ, of the process outputs at the t-th sampling instant of the first windowy,rThe standard deviation vector of the process output of the first window, myDenotes the number of output variables, myN denotes the number of samples of the first window, n is 500, y1(t) represents a first output variable,
Figure BDA00030521101600001714
denotes the m-thyAn output variable; obtaining a mean vector and a standard deviation vector of a first window process-disturbance combined model residual error according to the following formula:
Figure BDA00030521101600001715
Figure BDA00030521101600001716
wherein the content of the first and second substances,
Figure BDA00030521101600001717
mean vector, σ, of the process-perturbation combined model residuals for the first windowe,rThe standard deviation vector of the process-perturbation combined model residuals of the first window, myDenotes the number of output variables, my＝1，e1(t) represents the first process-perturbation joint model residual,
Figure BDA00030521101600001718
denotes the m-thyProcess-perturbation joint model residuals.
(5.3) constructing matrix X according to the following equation1And Y1
Figure BDA00030521101600001719
Wherein u isiRepresents the ith process input vector, i 1, …, mu，ui"represents the i-th process input vector after normalization, i.e.
Figure BDA00030521101600001720
Figure BDA00030521101600001721
Means representing the ith process input vector, (σ)u,r)iRepresents the standard deviation of the ith process input vector; y isjRepresents the jth process output vector, j ═ 1, …, my，yj"represents the normalized j-th process output vector, i.e.
Figure BDA00030521101600001722
Figure BDA00030521101600001723
Represents the mean of the jth process output vector (σ)y,r)jRepresents the standard deviation of the jth process output vector; e.g. of the typej(t) denotes the j-th process-perturbation combined model residual vector, j is 1, …, my，ej＝(ej(t))，ej(t)＝yj(t)-(Gm(q-1))juj(t)，(Gm(q-1))jRepresentation process model Gm(q-1) J-th line of (1), yj(t) denotes the jth process output at the tth sampling instant uj(t) denotes the process input at the tth sampling instant, e ″jRepresents the normalized j process-perturbation joint model residual vector,
Figure BDA00030521101600001724
Figure BDA00030521101600001725
means for representing the jth process-perturbation combined model residual, (σ)e,r)jRepresents the standard deviation of the jth process-perturbation combined model residual.
(5.4) construction of Hankel matrix X according to the following formulapAnd Xf：
Figure BDA0003052110160000181
Wherein m is mu+2my，pp+1Representing the historical data vector at the p +1 th sampling instant, i.e. pp+1＝[x(p)T … x(1)T y(p)T … y(1)T]TX (p) represents a matrix
Figure BDA0003052110160000182
P column of (2), x (1) represents a matrix
Figure BDA0003052110160000183
The number 1. sup. st column of (a),
Figure BDA0003052110160000184
representation matrix X1Is transferred toY (p) denotes a matrix Y1 TY (1) denotes the matrix Y1 TColumn 1, Y1 TRepresentation matrix Y1Transposing; f. ofp+1The vector of prediction data representing the p +1 th sampling instant, i.e. fp+1＝[y(p+1)T y(p+2)T … y(p+f-1)T]TY (p +1) represents Y1 TColumn p +1, Y (p +2) represents Y1 TIn the p +2 th column, Y (p + f-1) represents Y1 TColumn p + f-1, M denotes the number of vectors of the history data vector and the prediction data vector constructed from n sets of observation data, M is n-p-f +1, n is 500, p is 10, f is 10, and M is 481.
(6) Obtaining a projection matrix Jq
(6.1) obtaining the unitary matrix V according to: h ═ U ∑ VTWherein V represents a unitary matrix constructed when singular value decomposition is performed on the matrix H,
Figure BDA0003052110160000185
is represented by XpAnd XfThe constructed matrix Hankel matrix is used as a matrix,
Figure BDA0003052110160000186
represents the sample covariance of the historical data vector,
Figure BDA0003052110160000187
represents the sample covariance of the prediction data vector,
Figure BDA0003052110160000188
representing the cross covariance of the historical data vector and the predicted data vector.
(6.2) obtaining the projection matrix J according to the following formulaq：
Figure BDA0003052110160000189
Wherein, VqRepresenting a matrix constructed from the first q columns of the unitary matrix V,
Figure BDA00030521101600001810
represents the Hankel matrix XpThe covariance matrix of (2).
(7) Obtaining process-perturbation joint model residuals for the kth window, k ═ 2,3, …, S +1
Acquiring a process-disturbance combined model residual error of the kth window by using the input and output data of the kth window acquired in the step (1): e.g. of the typek'(t)＝yk'(t)-Gm,k'(q-1)uk'(t) wherein ek'(t) represents the process-perturbation joint model residual at the tth sampling instant of the kth' window, yk'(t) denotes the kth' window process output vector at the tth sampling instant uk'(t) represents the process input vector for the kth window at the tth sampling instant, Gm,k'(q-1) Representing the Process model, q, used by the EWMA batch controller in the Current operating regime-1A difference operator is represented.
(8) Construction of a Hankel matrix X according to the following formulap,k'And Xf,k'
(8.1) acquiring a process input dataset according to
Figure BDA00030521101600001811
And outputting the data set
Figure BDA00030521101600001812
Figure BDA00030521101600001813
Wherein the content of the first and second substances,
Figure BDA00030521101600001814
process input vector, m, representing the kth windowuDenotes the number of input variables, mu＝1，
Figure BDA00030521101600001815
The process-perturbation joint model residual vector representing the kth window,
Figure BDA00030521101600001816
process output vector, m, representing the kth windowyDenotes the number of output variables, my＝1。
(8.2) constructing matrix Xk'And Yk'：
Figure BDA00030521101600001817
Wherein u isi,k'The process input vector, u ", representing the k' th windowi,k'The process input vector representing the normalized k' th window, i.e.
Figure BDA00030521101600001818
yj,k'Represents the process output vector, y ″, of the k' th windowj,k'The process output vector representing the normalized k' th window, i.e.
Figure BDA00030521101600001819
ej,k'(t) Process-perturbation Joint model residual error, e, for the kth windowj,k'＝(ej,k'(t))， ej,k'(t)＝yj,k'(t)-(Gm,k′(q-1))juk'(t)，yj,k'(t) denotes the jth process output, u, at the tth sampling instant in the kth windowk'(t) represents the column vector formed by all process inputs at the tth sampling instant in the kth window, (G)m,k′(q-1))jRepresents Gm，k′(q-1) Line j, e ″)j,k'Representing the process-perturbation combined model residual of the normalized k' th window, i.e.
Figure BDA00030521101600001820
(8.3) construction of Hankel matrix X according to the following formulap,k'And Xf,k'：
Figure BDA0003052110160000191
Wherein p isp+1kRepresenting the historical data vector at the p +1 th sampling instant, i.e. pp+1,k'＝[xk'(p)T … xk'(1)Tyk'(p)T … yk'(1)T]T，xk'(p) a representation matrix
Figure BDA0003052110160000192
P column of (2), xk'(1) Representation matrix
Figure BDA0003052110160000193
The number 1. sup. st column of (a),
Figure BDA0003052110160000194
representation matrix Xk'Transpose of (y)k'(p) a representation matrix
Figure BDA0003052110160000195
P column, yk'(1) Representation matrix
Figure BDA0003052110160000196
The number 1. sup. st column of (a),
Figure BDA0003052110160000197
representation matrix Yk'Transposing; f. ofp+1,k'Representing the vector of predicted data at the p +1 th sampling instant in the k' th window, i.e. fp+1,k'＝[y(k')(p+1)Ty(k')(p+2)T … y(k')(p+f-1)T]T，y(k')(p +1) represents
Figure BDA0003052110160000198
P +1 column, y(k')(p +2) represents
Figure BDA0003052110160000199
P +2 column, y(k')(p + f-1) represents
Figure BDA00030521101600001910
Column p + f-1, M indicates the number of vectors of past and future data vectors constructed from n sets of observation data, M is n-p-f +1, n is 500, p is 10, f is 10, and M is 481.
(9) Obtaining a control performance monitoring index eta of an EWMA batch control systemCPIIndex (I)
Acquiring a control performance monitoring index eta of an EWMA batch control system according to the following formulaCPIIndexes are as follows:
Figure BDA00030521101600001911
wherein the content of the first and second substances,
Figure BDA00030521101600001912
represents the amount of control performance monitoring, J, of the EWMA batch-to-batch control system corresponding to the kth windowqRepresenting a projection matrix, pk,k'Representation matrix Xp,k'The k-th column of (1).
(10) Obtaining a reference condition model residual error
Using formulas
Figure BDA00030521101600001913
Obtaining a reference condition model residual error of the kth sampling moment of the kth window, wherein,
Figure BDA00030521101600001914
representing the reference condition disturbance model estimated value y in the step (3.3)k'(t) represents the process output vector of the kth window at the tth sampling instant uk'(t) represents the process input vector for the kth window at the tth sampling instant, Gm,k'(q-1) Representing the process model used by the EWMA batch controller in the k' th window.
(11) Obtaining model mismatch/disturbance dynamic combined monitoring quantity
Figure BDA00030521101600001915
Using formulas
Figure BDA00030521101600001916
Obtaining the model mismatch/disturbance dynamic combined monitoring quantity, wherein,
Figure BDA00030521101600001917
represents the model mismatch/disturbance dynamic joint monitoring quantity corresponding to the k' th window,
Figure BDA00030521101600001918
has a value range of (0, 1)]，εr(t) white noise estimate for the reference condition at the tth sampling time, cov (εr(t)) represents εrCovariance of (t), cov (v)r,k'(t)) represents vr,k'(t) covariance.
(12) Centralizing the closed-loop input and output data of the k' th window
And respectively carrying out centralized processing on the closed-loop input and the process output according to the following formula by utilizing the closed-loop input and output data of the kth window:
Figure BDA00030521101600001919
wherein u'k'(t) represents the process input vector at the t-th sampling instant of the k' -th window after centering, uk'(t) represents the process input vector at the t-th sampling instant of the k' -th window that is not centered,
Figure BDA00030521101600001920
mean vector, m, representing the process input at the t-th sampling instant of the k' -th windowuDenotes the number of input variables, m u1, t denotes the t-th sampling time, n denotes the number of samples of the k '-th window, and n is 500, y'k'(t) represents the process output vector, y, of the k' th window after centeringk'(t) represents the process output vector at the t-th sampling instant of the k' -th window that is not centered,
Figure BDA00030521101600001921
mean vector, m, representing process output at the t-th sampling instantyRepresenting output variablesNumber, my＝1。
(13) Obtaining the estimated value of white noise of actual working condition
Firstly, according to the process output y of the centralization processing obtained in the step (12)k′'(t) constructed data vector: y'P,k'(t)＝[y′k'(t) … y′k'(t-P+1)](ii) a Then, by vector y'P,k'(t) constructed data matrix: y'P,k'(t)＝[y′P,k'(t-1)T … y′P,k'(t-L)T]TWhere P denotes the data window size of the estimated white noise. Obtaining an estimated value epsilon of white noise of the reference working condition according to the following formulah,k'(t)＝y′P,k'(t){I-Y′P,k'(t)T[Y′P,k'(t)Y′P,k'(t)T]-1Y′P,k'(t)}。
(14) Obtaining an estimated value of a disturbance model of an actual working condition
(14.1) obtaining an autoregressive moving average model of the out-of-band source input for the process output of the kth window according to:
Figure BDA0003052110160000201
wherein the content of the first and second substances,
Figure BDA0003052110160000202
which is indicative of the delay factor, is,
Figure BDA0003052110160000203
the coefficients of the auto-regressions are represented,
Figure BDA0003052110160000204
the coefficient representing the input of the external source,
Figure BDA0003052110160000205
coefficient representing moving average, Ik' represents the maximum order of autoregressive, Kk' represents the maximum order of exogenous input; j. the design is a squarek' indicating movementAverage maximum order;
(14.2) according to the self-adaptive minimum absolute value contraction and operator selection method, obtaining the coefficient vector and the order of the autoregressive moving average model input by the out-of-band source in the step (4.3):
Figure BDA0003052110160000206
wherein the content of the first and second substances,
Figure BDA0003052110160000207
Figure BDA0003052110160000208
Wk' denotes the window size of the sampled data,
Figure BDA0003052110160000209
represents WkA data matrix of dimension' x 1,
Figure BDA00030521101600002010
represents Wk'×(Ik'+Kk'+Jk') A dimensional data matrix;
Figure BDA00030521101600002011
the method representing the adaptive minimum absolute value shrinkage and selection operator estimates the coefficient vector of the autoregressive moving average model, i.e.
Figure BDA00030521101600002012
Figure BDA00030521101600002013
Coefficient, phi, representing the estimated autoregressive moving average modelk'A coefficient vector representing a real but unknown type of auto-regressive moving average,
Figure BDA00030521101600002014
Figure BDA00030521101600002015
j-th coefficient phi in coefficients representing autoregressive moving average modelj,k'The weighting factor of (2).
(14.3) obtaining an estimated value of the actual working condition disturbance model according to the autoregressive moving average model of the out-of-band source input obtained in the step (13.2):
Figure BDA00030521101600002016
wherein the content of the first and second substances,
Figure BDA00030521101600002017
representing the coefficients of the auto-regressive part and the coefficients of the moving average part estimated in step (13.2),
Figure BDA00030521101600002018
and representing the estimated value of the actual working condition disturbance model.
(15) Obtaining a model residual error of an actual working condition:
using formulas
Figure BDA00030521101600002019
Wherein the content of the first and second substances,
Figure BDA00030521101600002020
representing the actual condition disturbance model estimated value, G, obtained by the step (14.3)m,k′(q-1) Representing the process model used by the EWMA batch controller in the k' th window.
(16) Obtaining interference characteristic monitoring quantity etaDVI
Using formulas
Figure BDA00030521101600002021
Obtaining interference characteristic monitoring quantity
Figure BDA00030521101600002022
Wherein epsilonh,k'(t) represents the white noise estimate for the actual condition at the tth sampling instant in the kth' window, | cov (εh,k'(t)) | denotes εh,k'(t) of the covariance matrixDeterminant, | cov (εr(t)) | denotes εr(t) determinant of covariance matrix.
(17) Obtaining a monitoring quantity eta of the disturbance modelDMI
Using formulas
Figure BDA00030521101600002023
Obtaining disturbance model monitoring quantity based on KL divergence
Figure BDA00030521101600002024
Figure BDA00030521101600002025
Has a value range of (0, 1)]If, if
Figure BDA00030521101600002026
Near or equal to 1, the perturbation model is unchanged, vr,k'(t) reference condition model residual, v, at the tth sampling timeh,k'(t) represents the actual condition model residual at the t-th sampling instant in the k' -th window, J (f (v)r,k'(t)),g(vh,k'(t))) represents vh,k'(t) and vr,k'(t) KL divergence, expressed as:
Figure BDA00030521101600002027
wherein, f (v)r,k'(t)) represents v constructed using the nuclear density estimation methodr,k'(t) probability density function, f (v)h,k'(t)) represents v acquired by the kernel density estimation methodh,k'(t) probability density function, i ═ 1,2, …, Nf，NfRepresenting the total number of independent identically distributed sample points produced by the kernel density estimate.
(18) Obtaining model mismatch monitoring quantity etaPMI
Using formulas
Figure BDA0003052110160000211
Obtaining model mismatch monitoring quantity
Figure BDA0003052110160000212
Wherein epsilonh,k'(t) represents the white noise estimate of the actual condition at the t-th sampling instant in the k' -th window, vh,k'(t) represents the real-world model residual at the t-th sampling instant in the kth' window, | cov (εh,k'(t)) | denotes εh,k'Determinant of covariance matrix of (t), | cov (v)h,k'(t)) | denotes vh,k'(t) determinant of covariance matrix.
(19) Obtaining control limits of control performance monitoring quantity, model mismatch/disturbance dynamic monitoring quantity, interference characteristic monitoring quantity, disturbance model monitoring quantity and model mismatch monitoring quantity
(19.1) obtaining the control upper limit UP of the control performance monitoring amount according to the following formulaCPI：
Figure BDA0003052110160000213
Wherein the content of the first and second substances,
Figure BDA0003052110160000214
the control performance monitoring amount corresponding to the kth window is K ', K' is 2, 3.., S +1, h represents the kernel bandwidth, and K (·) represents the kernel density estimation with the kernel being a radial basis function.
(19.2) obtaining the control upper limit UP of the model mismatch/disturbance dynamic combined monitoring quantity according to the following formulaPDM：
Figure BDA0003052110160000215
Wherein the content of the first and second substances,
Figure BDA0003052110160000216
is the model mismatch/disturbance dynamic joint monitoring quantity corresponding to the k' th window.
(19.3) obtaining the control upper limit UP of the disturbance characteristic monitoring amount according to the following formulaDVAnd a lower control limit LPDV：
Figure BDA0003052110160000217
Wherein the content of the first and second substances,
Figure BDA0003052110160000218
is the interference characteristic monitoring quantity corresponding to the k' th window.
(19.4) obtaining the lower control limit LP of the disturbance model monitored quantity according to the following formulaDM：
Figure BDA0003052110160000219
Wherein the content of the first and second substances,
Figure BDA00030521101600002110
is the disturbance model monitoring quantity corresponding to the k' th window.
(19.5) obtaining the control upper limit UP of the model mismatch monitoring amount according to the following formulaPMI：
Figure BDA00030521101600002111
Wherein the content of the first and second substances,
Figure BDA00030521101600002112
is the interference characteristic monitoring quantity corresponding to the k' th window.
(20) Construction of Hankel matrix Xp,lAnd Xf,l
(20.1) acquiring N1 groups of input data and output data of the EWMA batch control system of the semiconductor process on line as actual working condition sample data, and dividing the actual working condition sample data into S1An online window, each window containing n groups of samples, wherein S1＝N1-n+1。
(20.2) obtaining an online process-disturbance combined model residual according to the following formula:
el(t)＝yl(t)-Gm,l(q-1)ul(t)，
wherein e isl(t) represents the on-line process-perturbation joint model residual error at the tth sampling instant of the ith window, yl(t) represents the process output vector at the tth sampling instant of the ith window, ul(t) denotes the input vector at the tth sampling instant of the l-th window, t denotes the sampling instant, Gm(q-1) Representing a process model, q, used by the EWMA lot controller in the current operating regime-1Representing a difference operator,/＝1,2,…,S1。
(20.3) acquiring a process input data set of the current working condition by using the process output vector, the process input vector and the process-disturbance combined model residual error which are acquired in real time
Figure BDA00030521101600002113
And output data set Yl 0：
Figure BDA00030521101600002115
Wherein the content of the first and second substances,
Figure BDA00030521101600002116
representing an online process input vector, mu＝1，
Figure BDA00030521101600002117
Representing the on-line process-perturbation joint model residual vector,
Figure BDA00030521101600002118
representing the on-line process output vector, my＝1。
(20.4) constructing matrix X according tolAnd Yl：
Figure BDA00030521101600002119
Wherein u isi,lRepresents the process input vector, u ″, in the ith window of the online datai,lThe process input vector representing the normalized l-th window, i.e.
Figure BDA00030521101600002120
yj,lRepresents the process output vector, y ″, of the ith window of online dataj,lThe process output vector representing the normalized l-th window, i.e.
Figure BDA0003052110160000221
ej,lProcess-perturbation joint model residual vector representing the ith window of online data, ej,l＝(ej,l(t))， ej,l(t)＝yj,l(t)-(Gm,l(q-1))jul(t)，(Gm,l(q-1))jProcess model G representing EWMA inter-batch controller usage in the ith window of online datam，l(q-1) J-th line of (1), yj,l(t) denotes the jth process output at the tth sampling instant of the ith window of online data, ul(t) represents the vector formed by all process inputs at the tth sampling time of the ith window of online data, e ″j,lRepresenting the normalized process-perturbation joint model residual,
Figure BDA0003052110160000222
(20.5) construction of Hankel matrix X according to the following formulap,lAnd Xf,l：
Figure BDA0003052110160000223
Wherein p isp+1,lRepresenting a history data vector at the p +1 th sampling instant in the ith window of line data, i.e. pp+1,l＝[xl(p)T … xl(1)T yl(p)T … yl(1)T]T，xl(p) a representation matrix
Figure BDA0003052110160000224
P column of (2), xl(1) Representation matrix
Figure BDA0003052110160000225
The p-th column of (2),
Figure BDA0003052110160000226
representation matrix XlTranspose of (y)l(p) represents a matrix Yl TP column, yl(1) To representMatrix Yl TColumn 1, Yl TRepresentation matrix YlTransposing; f. ofp+1,lThe vector of predicted data representing the p +1 th sampling instant in the l window of line data, i.e. fp+1,l＝[yl(p+1)T yl(p+2)T … yl(p+f-1)T]T，yl(p +1) represents Yl TP +1 column, yl(p +2) represents Yl TP +2 column, yl(p + f-1) represents Yl TColumn p + f-1, M indicates the number of vectors of past data vectors and future data vectors constructed from n sets of observation data, M is n-p-f +1, n is 500, p is 10, f is 10, and M is 481.
(21) Monitoring control performance of EWMA batch-to-batch control system
Obtaining the control performance monitoring quantity eta of the current working condition according to the following formula by utilizing the online dataCPI：
Figure BDA0003052110160000227
Wherein the content of the first and second substances,
Figure BDA0003052110160000228
indicating the control performance monitoring quantity, p, of the current operating mode corresponding to the ith windowk,lRepresentation matrix Xp,lThe kth column vector of (1).
According to the control upper limit UP of the control performance monitoring quantity obtained in the step (19.1)CPIJudgment of
Figure BDA0003052110160000229
Whether or not less than UPCPIIf, if
Figure BDA00030521101600002210
Less than UPCPIThen the system performance is better, if
Figure BDA00030521101600002211
Greater than UPCPIThe system performance is degraded.
(22) Preliminary diagnosis of factors leading to system performance degradation
(22.1) acquiring a reference working condition model residual error according to the following formula by using the online data:
Figure BDA00030521101600002212
wherein v isr,l(t) represents the reference condition model residual error at the tth sampling time of the ith window of the online data,
Figure BDA00030521101600002213
and representing the reference working condition disturbance model estimation value.
(22.2) obtaining the model mismatch/disturbance dynamic combined monitoring quantity eta according to the following formulaPDMI：
Figure BDA00030521101600002214
Wherein the content of the first and second substances,
Figure BDA00030521101600002215
has a value range of (0, 1)]，vr,l(t) represents the reference condition model residual at the tth sampling instant, | cov (v)r,l(t)) | denotes vr,l(t) determinant of covariance matrix.
(22.3) calculating model mismatching/disturbance dynamic combined monitoring quantity eta by using online collected input data and output dataPDMIUpper control limit UP combined with model mismatch/disturbance dynamic joint monitoring quantityPDMIThe factors for controlling the performance degradation are preliminarily diagnosed if the current working condition is met
Figure BDA00030521101600002216
Greater than UPPDMIIf the model mismatch or the external disturbance dynamics is changed, the system performance is reduced due to the model mismatch or the external disturbance dynamics; if the current operating condition is
Figure BDA00030521101600002217
Less than UPPDMIThen changes in the controller parameters are a factor that causes degradation in system performance.
(23) Further diagnosis of factors leading to system performance degradation
(23.1) acquiring an estimated value epsilon of white noise of the actual working condition according to the steps (13) to (15) by utilizing the online datah,l(t), actual condition disturbance model estimation value
Figure BDA00030521101600002218
Actual condition model residual vh,l(t)。
(23.2) obtaining the interference characteristic monitoring quantity eta of the on-line data according to the following formulaDVI：
Figure BDA00030521101600002219
Wherein the content of the first and second substances,
Figure BDA00030521101600002220
representing the interference characteristic monitoring quantity of the on-line data corresponding to the ith window, representing the change of white noise variance in actual working conditionr(t) represents the white noise estimation value of the reference condition at the t-th sampling moment, epsilonh,l(t) represents the white noise estimate for the actual condition at the tth sampling instant, | cov (εh,l(t)) | denotes εh,l(t) covariance matrix determinant, | cov (εr(t)) | denotes εr(t) covariance matrix determinant.
(23.3) calculating the actual working condition model residual error of the online data:
Figure BDA0003052110160000231
wherein v ish,l(t) represents the model residual error of the online data actual condition corresponding to the ith window,
Figure BDA0003052110160000232
on-line data actual condition disturbance model estimation value y representing ith windowl(t) an online data process output vector, u, at the tth sampling instant of the ith windowl(t) an online data process input vector, G, for the tth sampling instant of the l-th windowm,l(q-1) Represents the process model used by the EWMA batch controller in the ith window of the online data.
(23.4) obtaining the disturbance of the on-line data according to the following formulaModel monitoring quantity etaDMI：
Figure BDA0003052110160000233
Wherein the content of the first and second substances,
Figure BDA0003052110160000234
representing the monitored quantity, v, of the disturbance model corresponding to the l-th windowr,l(t) reference condition model residual, v, at the tth sampling timeh,l(t) represents the real model residual at the tth sampling instant, J (f (v)r,l(t)),g(vh,l(t))) represents vh,l(t) and vr,l(t) KL divergence.
(23.5) obtaining a model mismatch monitoring quantity eta of the on-line data according to the following formulaPMI：
Figure BDA0003052110160000235
Wherein the content of the first and second substances,
Figure BDA0003052110160000236
model mismatch monitoring quantity, epsilon, representing the on-line data corresponding to the ith windowh,l(t) the estimated value of the white noise of the actual working condition at the t-th sampling moment obtained in the step (23.1), vh,l(t) represents the actual condition model residual error, | cov (ε), of the t-th sampling time obtained in step (23.1)h,l(t)) | denotes εh,l(t) determinant of covariance matrix, | cov (v)h,l(t)) | denotes vh,l(t) determinant of covariance matrix.
(23.6) disturbance feature monitoring amount calculated using on-line data
Figure BDA0003052110160000237
Monitoring quantity of disturbance model
Figure BDA0003052110160000238
Amount of model mismatch monitoring
Figure BDA0003052110160000239
Lower control limit LP combined with disturbance characteristic monitoring quantityDVIControl of the amount of monitoring of the interference characteristicUpper limit UPDVILower control limit LP of disturbance model monitoring quantityDMIAnd upper control limit UP of model mismatch monitoring quantityPMIFurther diagnosing and identifying the factors causing the performance degradation of the EWMA batch-to-batch control system, and if the characteristic monitoring quantity is disturbed
Figure BDA00030521101600002310
Greater than UPDVIOr less than LPDVIIf so, the white noise variance changes to cause the performance of the EWMA batch-to-batch control system to be reduced; if the model monitoring quantity is disturbed
Figure BDA00030521101600002311
Less than LPDMIIf so, the performance of the EWMA batch control system is reduced due to the change of the external disturbance model; monitoring quantity if model mismatching
Figure BDA00030521101600002312
Greater than UPPMIThen the model mismatch causes the EWMA inter-batch control system performance to degrade.
In the reference operating mode, the control parameter of the controller is λ ═ 0.35, and the object model used in the controller: gm(q-1)＝159.3q-1I.e., b is 159.3, the variance of white noise is
Figure BDA00030521101600002313
The present simulation will consider six cases shown in table 1:
TABLE 1 parameter settings
Figure BDA00030521101600002314
When a model mismatch occurs, fig. 3 shows the trend of the change in the control performance index and each diagnostic index. As can be seen from FIG. 3(a), the control performance indicator ηCPIExceeding the upper control limit from the 481 st sampling instant indicates a decrease in system performance, and therefore further diagnosis of the cause of system performance degradation is required. As can be seen from FIG. 4(b), η is after the 481-th sampling timingPDMIBeyond its upper control limit UPPDMIModel mismatch or perturbation dynamics changes are accounted for. By observing fig. 3(c) - (e), starting from the 481 th sampling timing, the interference characteristic monitoring amount ηDVIAnd the monitored quantity eta of the disturbance modelDMIAre all kept within the control interval, i.e. LPDVI＜ηDVI＜UPDVIAnd η DMI1, but the model mismatch monitor ηPMIBeyond its control limit, i.e. ηPMI＞UPPMIThis shows that neither the disturbance signature (or white noise variance) nor the disturbance model of the EWMA batch-to-batch control system has changed, but the process model of the batch-to-batch control system does not match the actual process and results in a degradation of system performance. Therefore, the method of the present invention effectively evaluates system performance and identifies model mismatch from factors that control performance degradation.
Fig. 4 shows the variation tendency of the control performance index and each diagnostic index when the white noise variance of the external disturbance changes. As can be seen from FIG. 4(a), the control performance indicator ηCPIExceeding the upper control limit from the 481-th sampling instant indicates a decrease in system performance, and therefore further diagnosis of the cause of system performance degradation is required. As can be seen from FIG. 4(b), η is after the 481-th sampling timingPDMIBeyond its upper control limit UPPDMIModel mismatch or perturbation dynamics changes are illustrated. By observing fig. 4(c) - (e), the model mismatch monitor amount η starts from the 481 th sampling time instantPMIAnd the monitored quantity eta of the disturbance modelDMIAre all kept within a control region, i.e. ηPMI＜UPPMIAnd η DMI1, but the interference characteristic monitoring quantity etaDVIBeyond its control line, i.e. etaDVI＜LPDVIThis indicates that the process model and the actual process matching and disturbance model of the EWMA batch control system have not changed, but the disturbance signature (or white noise variance) has changed and resulted in a degradation of system performance. Thus, the method of the present invention effectively evaluates system performance and identifies changes in the interference signature (or white noise variance) from factors that control performance degradation.
FIG. 5 shows changes in control performance indicators and diagnostic indicators as controller parameters changeTrend. As can be seen from FIG. 5(a), the control performance indicator ηCPIExceeding the upper control limit from the 481-th sampling instant indicates a decrease in system performance, and therefore further diagnosis of the cause of system performance degradation is required. As can be seen from FIG. 5(b), η is after the 481-th sampling timingPDMIWithin its control limits, i.e. etaPDMI＜UPPDMIAccounting for the controller parameter changes and resulting in a degradation of system performance. By observing fig. 5(c) - (e), the model mismatch monitor amount η starts from the 481 th sampling time instantPMIInterference characteristic monitoring quantity etaDVIMonitoring quantity eta of disturbance modelDMIAre all kept within a control region, i.e. ηPMI＜UPPMI、 LPDVI＜ηDVI＜UPDVIAnd ηDMIThis indicates that the process model of the EWMA batch-to-batch control system matches the actual process, and neither the disturbance signature (or white noise variance) nor the disturbance model has changed. Thus, the method of the present invention effectively evaluates system performance and identifies changes in controller parameters from factors that control performance degradation.
Fig. 6 shows the variation tendency of the control performance index and each diagnostic index when the disturbance model is changed. As can be seen from FIG. 6(a), the control performance indicator ηCPIExceeding the upper control limit from the 481 st sampling instant indicates a decrease in system performance, and therefore further diagnosis of the cause of system performance degradation is required. As can be seen from FIG. 6(b), η is after the 481-th sampling timingPDMIBeyond its upper control limit UPPDMIModel mismatch or perturbation dynamics changes are accounted for. By observing fig. 6(c) - (e), the model mismatch monitor amount η starts from the 481 th sampling time instantPMIAnd interference characteristic monitoring quantity etaDVIAre all kept within a control region, i.e. ηPMI＜UPPMIAnd LPDVI＜ηDVI＜UPDVIHowever, the disturbance model monitoring quantity etaDMIBeyond its control line, i.e. etaDMI＜LPDMIThis indicates that the process model and actual process matching and disturbance characteristics (or white noise variance) of the EWMA batch-to-batch control system have not changed, but the disturbance model has changed and resulted in a system performanceAnd (5) reducing. Therefore, the method of the invention effectively evaluates the system performance and identifies the change of the disturbance model from the factors of control performance degradation.
Fig. 7 shows the trend of the change in the control performance index and each diagnostic index when the model mismatch and the white noise variance change occur simultaneously. As can be seen from FIG. 7(a), the control performance index ηCPIExceeding the upper control limit from the 481-th sampling instant indicates that the system performance is degraded, and therefore further diagnosis of the cause of system performance degradation is required. As can be seen from FIG. 7(b), η is after the 481-th sampling timingPDMIBeyond its upper control limit UPPDMIModel mismatch or perturbation dynamics changes are accounted for. By observing fig. 7(c) - (e), starting from the 481 th sampling timing, the interference characteristic monitored quantity ηDVIAnd model mismatch monitor ηPMIAre all changed, i.e. etaDVI＜LPDVIAnd ηPMI＞UPPMIHowever, the disturbance model monitoring quantity etaDMIStill within its control limits, i.e. ηDMIThis indicates that the disturbance model of the EWMA batch-to-batch control system has not changed, but the process model and the actual process mismatch and the disturbance signature (or white noise variance) have changed, and both have resulted in a degradation of the system performance. Thus, the method of the present invention effectively evaluates system performance and identifies interference signature changes and model mismatches from factors that control performance degradation.
When the models are mismatched and the disturbance model is changed, fig. 8 shows the variation tendency of the control performance index and each diagnostic index. As can be seen from FIG. 8(a), the control performance index ηCPIExceeding the upper control limit from the 481 st sampling instant indicates a decrease in system performance, and therefore further diagnosis of the cause of system performance degradation is required. As can be seen from FIG. 8(b), η is after the 481-th sampling timingPDMIBeyond its upper control limit UPPDMIModel mismatch or perturbation dynamics changes are accounted for. By observing fig. 8(c) - (e), starting from the 481 th sampling timing, the interference characteristic monitoring amount ηDVIAnd the monitored quantity eta of the disturbance modelDMIAre all changed, i.e. etaDVI＜LPDVIAnd ηDMI＜LPDMIBut the amount of model mismatch monitoring ηPMIStill within its control limits, i.e. etaPMI＜UPPMIThis shows that the process model of the EWMA batch-to-batch control system matches the actual process, but both the disturbance signature (or white noise variance) and the disturbance model change, and both result in a degradation of the system performance. Thus, the method of the present invention effectively evaluates system performance and identifies changes in disturbance characteristics and disturbance models from factors that control performance degradation.
Through the analysis of the 1-6 and 6 fault conditions of the actual working conditions, the control performance diagnosis method of the semiconductor process batch control system based on the model residual error can be used for accurately diagnosing the factors of the system performance reduction.
The above description is only for the purpose of illustrating the preferred embodiments of the present invention and is not to be construed as limiting the invention, and any modifications, equivalents, improvements and the like that fall within the spirit and principle of the present invention are intended to be included therein.

Claims (10)

1. A control performance diagnosis method of a semiconductor process batch control system is characterized by comprising the following steps:
s1, obtaining N groups of input data and output data of historical normal working conditions of the EWMA batch control system in the semiconductor manufacturing process as sample data, dividing the sample data into S +1 windows, and calculating a white noise estimation value epsilon of the reference working condition by using the input data and the output data of the first windowrAnd reference working condition disturbance model estimation value
Figure FDA0003052110150000011
Obtaining process-disturbance combined model residual e (t), process data set
Figure FDA0003052110150000012
And output data set Y1 0And based on the process data set
Figure FDA0003052110150000013
Obtaining a matrix X of reference conditions1From the output data set Y1 0Obtaining a matrix Y of reference conditions1Construction of Hankel matrix XpAnd XfAnd projection matrix Jq；
S2, constructing a process data set using the input data and the output data of the k' th window
Figure FDA0003052110150000014
And outputting the data set
Figure FDA0003052110150000015
From a process data set
Figure FDA0003052110150000016
Obtaining matrix X of actual working conditionk'From the output data set
Figure FDA0003052110150000017
Obtaining matrix Y of actual working conditionsk'Obtaining the control performance monitoring quantity of the EWMA batch-to-batch control system corresponding to the kth window by adopting a typical variable analysis method
Figure FDA0003052110150000018
Wherein k ═ 2,3, …, S + 1;
s3, according to the reference working condition disturbance model estimated value in the step S1
Figure FDA0003052110150000019
Obtaining reference working condition model residual error v from input data and output data of k' th windowr,k'(t) and model mismatch/perturbation dynamics joint monitoring volume
Figure FDA00030521101500000110
S4, obtaining the actual working condition white noise estimation corresponding to the k 'th window by using the input data and the output data of the k' th windowValue epsilonh,k'(t) and actual condition disturbance model estimation value
Figure FDA00030521101500000111
And according to the actual working condition disturbance model estimated value
Figure FDA00030521101500000112
Obtaining the residual error v of the actual working condition modelh,k'(t)；
S5, according to the reference working condition white noise estimated value epsilon in the step S1rAnd the white noise estimate of the actual condition ε in step S4h,k'(t) calculating an interference characteristic monitoring quantity
Figure FDA00030521101500000113
According to the reference working condition model residual v in the step S3r,k'(t) and the actual condition model residual v in step S4h,k'(t) calculating the amount of monitoring of the disturbance model
Figure FDA00030521101500000114
According to the actual working condition white noise estimated value epsilon in the step S4h,k'(t) and the actual condition model residual vh,k'(t) obtaining a model mismatch monitoring quantity
Figure FDA00030521101500000115
S6, respectively obtaining the control upper limit UP of the control performance monitoring quantity by adopting a kernel density estimation methodCPIUpper control limit UP of combined monitoring quantity of model mismatch/disturbance dynamicsPDMIUpper limit UP of control of the monitored quantity of interference characteristicsDVILower limit of control LP of disturbance characteristic monitoring amountDVILower control limit LP of disturbance model monitoring quantityDMIAnd the upper control limit UP of the model mismatch monitoring quantityPMI；
S7, collecting N1 groups of input data and output data of the EWMA batch control system of the semiconductor process on line, and obtaining the control performance monitoring quantity eta of the current working conditionCPIAccording to the control upper limit UP in step S6CPIJudging whether the control performance of the closed-loop system is reduced or not in real time;
s8, if the control performance monitoring quantity etaCPIControl upper limit UP greater than control performance monitoring amountCPIIf the performance of the EWMA batch control system is reduced, executing the step S9, otherwise, the EWMA batch control system has good performance;
s9, using the reference working condition disturbance model estimated value in the step S1
Figure FDA00030521101500000116
The N1 groups of input data and output data in the step S7 acquire the model mismatch/disturbance dynamic combined monitoring quantity eta of the current working conditionPDMI；
S10, if the model mismatching/disturbance dynamic joint monitoring quantity etaPDMIControl upper limit UP larger than the model mismatch/disturbance dynamic joint monitoring amount in step S6PDMIIf the performance of the EWMA batch control system is reduced due to model mismatch or disturbance dynamics, executing the step S11, otherwise, the performance of the EWMA batch control system is reduced due to the change of the controller parameters;
s11, calculating the disturbance characteristic monitoring quantity eta of the current working conditionDVIMonitoring quantity eta of disturbance modelDMIAnd model mismatch monitor ηPMI(ii) a Disturbance characteristic monitoring quantity eta when current working conditionDVIGreater than the upper control limit UP of the disturbance characteristic monitoring amount in step S6DVIOr less than the lower control limit LP of the disturbance characteristic monitored quantityDVIWhen the interference characteristic changes, the performance of the EWMA batch control system is reduced; when the disturbance model monitors the quantity etaDMILower control limit LP less than disturbance model monitoring amountDMIMeanwhile, the performance of the EWMA batch control system is reduced due to the change of the disturbance model; when model mismatch monitor etaPMIControl upper limit UP greater than model mismatch monitoring amountPMIModel mismatch results in performance degradation of the EWMA batch-to-batch control system.
2. The method of claim 1, wherein the matrix X of baseline operating conditions is a matrix of run-to-run control systems1And Y1Respectively expressed as:
Figure FDA0003052110150000021
Figure FDA0003052110150000022
wherein u isiRepresents the ith process input vector, i 1, …, mu，muDenotes the number of input vectors, t denotes the sampling time, u ″iRepresenting the i-th process input vector after normalization, i.e.
Figure FDA0003052110150000023
Figure FDA0003052110150000024
Means representing the ith process input vector, (σ)u,r)iRepresents the standard deviation of the ith process input vector; y isjRepresents the jth process output vector, j ═ 1, …, my，myIndicates the number of output vectors, y ″)jRepresenting the normalized jth process output vector, i.e.
Figure FDA0003052110150000025
Figure FDA0003052110150000026
Means for representing the jth process output vector, (σ)y,r)jRepresents the standard deviation of the jth process output vector; e.g. of the typej(t) denotes the jth process-perturbation joint model residual vector, where ej＝(ej(t))，ej(t)＝yj(t)-(Gm(q-1))ju(t)，(Gm(q-1))jRepresenting the process model G used in the first windowm(q-1) J-th line of (1), yj(t) representsthe jth process output at the sampling instant t, u (t) representing the column vector, e ″, formed by all the process inputs at the sampling instant t in the first windowjRepresents the normalized jth process-perturbation joint model residual vector,
Figure FDA0003052110150000027
Figure FDA0003052110150000028
represents the mean of the jth process-perturbation combined model residual, (σ)e,r)jRepresenting the standard deviation of the jth process-perturbation combined model residual error;
the projection matrix JqComprises the following steps:
Figure FDA0003052110150000031
wherein, VqRepresenting a matrix constructed from the first q columns of the unitary matrix V,
Figure FDA0003052110150000032
represents the Hankel matrix XpThe covariance matrix of (2).
3. The method of claim 2, wherein the process data set is a process data set
Figure FDA0003052110150000033
And outputting the data set
Figure FDA0003052110150000034
Respectively expressed as:
Figure FDA0003052110150000035
Figure FDA0003052110150000036
wherein the content of the first and second substances,
Figure FDA0003052110150000037
the process input vector representing the kth window,
Figure FDA0003052110150000038
the process-perturbation joint model residual vector representing the kth window,
Figure FDA0003052110150000039
a process output vector representing the kth window;
matrix X of the actual working conditionk'And Yk'Respectively expressed as:
Figure FDA00030521101500000310
Figure FDA00030521101500000311
wherein u isi,k'The process input vector, u ", representing the k' th windowi,k'The process input vector representing the normalized k' th window, i.e.
Figure FDA00030521101500000312
yj,k'Represents the process output vector, y ″, of the k' th windowj,k'The process output vector representing the normalized k' th window, i.e.
Figure FDA00030521101500000313
ej,k'Process-perturbation joint model residual vector representing the kth window, ej,k'＝(ej,k'(t))，ej,k'(t)＝yj,k'(t)-(Gm,k′(q-1))juk'(t)，yj,k'(t) denotes the jth process output, u, at the tth sampling instant in the kth windowk'(t) represents the column vector formed by all process inputs at the tth sampling instant in the kth window, (G)m,k′(q-1))jRepresenting the Process model G used by the EWMA batch-to-batch controller in the k' th windowm，k′(q-1) Line j, e ″)j,k'Representing the process-perturbation combined model residual of the normalized k' th window, i.e.
Figure FDA00030521101500000314
Obtaining control performance monitoring quantity of EWMA batch control system by adopting typical variable analysis method
Figure FDA00030521101500000315
Comprises the following steps:
Figure FDA00030521101500000316
wherein the content of the first and second substances,
Figure FDA00030521101500000317
represents the amount of control performance monitoring, J, of the EWMA batch-to-batch control system corresponding to the kth windowqRepresenting a projection matrix, pk,k'Representation matrix Xp,k'The k-th column of (1).
4. The method as claimed in claim 1 or 3, wherein the model mismatch/disturbance dynamics joint monitoring is
Figure FDA0003052110150000041
Comprises the following steps:
Figure FDA0003052110150000042
wherein the content of the first and second substances,
Figure FDA0003052110150000043
represents the model mismatch/disturbance dynamic joint monitoring quantity corresponding to the k' th window,
Figure FDA0003052110150000044
has a value range of (0, 1)]，εr(t) the white noise estimation value of the reference condition at the t-th sampling moment, vr,k'(t) represents the reference condition model residual at the t-th sampling instant in the kth' window, | cov (εr(t)) | denotes εrDeterminant of covariance matrix of (t), | cov (v)r,k'(t)) | denotes vr,k'(t) determinant of covariance matrix;
the reference working condition model residual vr,k'(t) is:
Figure FDA0003052110150000045
wherein v isr,k'(t) represents the reference condition model residual error corresponding to the k' th window,
Figure FDA0003052110150000046
disturbance model estimation value representing reference working condition
Figure FDA0003052110150000047
Inverse of (y)k'(t) represents the process output vector at the t-th sampling instant of the k' -th window, uk'(t) represents the process input vector at the t-th sampling instant in the k' -th window, Gm,k′(q-1) Representing the process model used by the EWMA batch controller in the k' th window.
5. The method of claim 4, wherein the actual process is performed in a batch control systemCondition model residual vh,k'(t) is:
Figure FDA0003052110150000048
wherein v ish,k'(t) represents the real model residual at the t-th sampling instant of the k' -th window,
Figure FDA0003052110150000049
actual condition disturbance model estimation value representing k' th window
Figure FDA00030521101500000410
Inverse of (y)k'(t) represents the process output vector at the t-th sampling instant of the k' -th window, uk'(t) represents the process input vector at the t-th sampling instant of the k' -th window, Gm,k′(q-1) Representing the process model used by the EWMA batch controller in the k' th window.
6. The method of claim 5, wherein the step S5 includes the steps of:
s5.1, calculating interference characteristic monitoring quantity
Figure FDA00030521101500000411
Figure FDA00030521101500000412
Wherein epsilonr(t) represents the white noise estimation value of the reference condition at the t-th sampling moment, epsilonh,k'(t) represents the white noise estimate for the actual condition at the t-th sampling instant in the kth' window, | cov (εh,k'(t)) | denotes εh,k'Determinant of covariance matrix of (t), | cov (εr(t)) | denotes εr(t) determinant of covariance matrix;
s5.2, calculating the monitoring quantity of the disturbance model
Figure FDA0003052110150000051
Figure FDA0003052110150000052
Wherein the content of the first and second substances,
Figure FDA0003052110150000053
representing the monitored amount of the perturbation model for the k' th window,
Figure FDA0003052110150000054
has a value range of (0, 1)]If, if
Figure FDA0003052110150000055
Near or equal to 1, the perturbation model is unchanged, vr,k'(t) represents the reference condition model residual at the t-th sampling instant in the k' -th window, vh,k'(t) represents the actual condition model residual at the t-th sampling instant in the k' -th window, J (f (v)r,k'(t)),g(vh,k'(t))) represents vh,k'(t) and vr,k'(t) KL divergence, expressed as:
Figure FDA0003052110150000056
wherein, f (v)r,k'(t)) represents v constructed using the nuclear density estimation methodr,k'(t) probability density function, f (v)h,k'(t)) represents v obtained by the nuclear density estimation methodh,k'(t) probability density function, i ═ 1,2, …, Nf，NfRepresenting the total number of independent identically distributed sample points resulting from the kernel density estimation;
s5.3, calculating model mismatch monitoring quantity
Figure FDA0003052110150000057
Figure FDA0003052110150000058
Wherein epsilonh,k'(t) represents the white noise estimate of the actual condition at the t-th sampling instant in the k' -th window, vh,k'(t) represents the real-world model residual at the t-th sampling instant in the kth' window, | cov (εh,k'(t)) | denotes εh,k'Determinant of covariance matrix of (t), | cov (v)h,k'(t)) | denotes vh,k'(t) determinant of covariance matrix.
7. The method of claim 1 or 5, wherein the step S6 includes the steps of:
s6.1, calculating the control upper limit UP of the control performance monitoring quantityCPI：
Figure FDA0003052110150000059
Wherein the content of the first and second substances,
Figure FDA00030521101500000510
the control performance monitoring quantity corresponding to the kth window is K ', K' is 2, 3.., S +1, h represents the kernel bandwidth, and K (·) represents the kernel density estimation with the kernel as a radial basis function;
s6.2, calculating the control upper limit UP of the model mismatch/disturbance dynamic combined monitoring quantityPDMI：
Figure FDA00030521101500000511
Wherein the content of the first and second substances,
Figure FDA0003052110150000061
model mismatch/disturbance dynamic combined monitoring quantity corresponding to the kth' window;
s6.3, calculating the control upper limit UP of the interference characteristic monitoring quantityDVIAnd a lower control limit LPDVI：
Figure FDA0003052110150000062
Wherein the content of the first and second substances,
Figure FDA0003052110150000063
is the interference characteristic monitoring quantity corresponding to the kth window;
s6.4, calculating a control lower limit LP of the monitoring quantity of the disturbance modelDMI：
Figure FDA0003052110150000064
Wherein the content of the first and second substances,
Figure FDA0003052110150000065
is the disturbance model monitoring quantity corresponding to the kth window;
s6.5, calculating the control upper limit UP of the model mismatch monitoring quantityPMI：
Figure FDA0003052110150000066
Wherein the content of the first and second substances,
Figure FDA0003052110150000067
is the interference characteristic monitoring quantity corresponding to the k' th window.
8. The method of claim 2 or 3, wherein the step S7 includes the steps of:
s7.1, collecting the semiconductor system on lineTaking N1 groups of input data and output data of the batch control system of the EWMA as actual working condition sample data, and dividing the actual working condition sample data into S1An online window, each window containing n groups of samples, wherein S1＝N1-n+1；
S7.2, calculating the residual error of the online process-disturbance combined model:
el(t)＝yl(t)-Gm,l(q-1)ul(t)，
wherein e isl(t) represents the on-line process-perturbation joint model residual error at the tth sampling instant of the ith window, yl(t) represents the process output vector at the tth sampling instant of the ith window, ul(t) denotes the input vector at the tth sampling instant of the l-th window, t denotes the sampling instant, Gm,l(q-1) Representing a Process model, q, used by an EWMA batch-to-batch controller in an on-line regime-1Denotes a difference operator, where l is 1,2, …, S1；
S7.3, acquiring a process input data set of the current working condition by using the process output vector, the process input vector and the process-disturbance combined model residual error which are acquired in real time
Figure FDA0003052110150000068
And output data set Y0 L：
Figure FDA00030521101500000610
Figure FDA00030521101500000611
Wherein the content of the first and second substances,
Figure FDA00030521101500000612
representing the on-line process input vector and,
Figure FDA00030521101500000613
representing the on-line process-perturbation joint model residual vector,
Figure FDA0003052110150000071
representing an online process output vector;
s7.4, constructing a matrix XlAnd Yl：
Figure FDA0003052110150000072
Figure FDA0003052110150000073
Wherein u isi,lRepresents the process input vector, u ″, in the ith window of the online datai,lThe process input vector representing the normalized l-th window, i.e.
Figure FDA0003052110150000074
yj,lRepresents the process output vector, y ″, of the ith window of online dataj,lThe process output vector representing the normalized l-th window, i.e.
Figure FDA0003052110150000075
ej,lProcess-perturbation joint model residual vector representing the ith window of online data, ej,l＝(ej,l(t))，ej,l(t)＝yj,l(t)-(Gm,l(q-1))jul(t)，(Gm,l(q-1))jProcess model G representing EWMA inter-batch controller usage in the ith window of online datam，l(q-1) J-th line of (1), yj,l(t) denotes the jth process output at the tth sampling instant of the ith window of online data, ul(t) represents the vector formed by all process inputs at the tth sampling time of the ith window of online data, e ″j,lRepresenting the normalized process-perturbation joint model residual,
Figure FDA0003052110150000076
s7.5, constructing a Hankel matrix Xp,lAnd Xf,l：
Figure FDA0003052110150000077
Figure FDA0003052110150000078
Wherein p isp+1,lRepresenting a history data vector at the p +1 th sampling instant in the ith window of line data, i.e. pp+1,l＝[xl(p)T … xl(1)T yl(p)T … yl(1)T]T，xl(p) a representation matrix
Figure FDA0003052110150000079
P column of (2), xl(1) Representation matrix
Figure FDA00030521101500000710
The p-th column of (2),
Figure FDA00030521101500000711
representation matrix XlTranspose of (y)l(p) represents a matrix Yl TP column, yl(1) Representation matrix Yl TColumn 1, Yl TRepresentation matrix YlTransposing; f. ofp+1,lThe vector of predicted data representing the p +1 th sampling instant in the l window of line data, i.e. fp+1,l＝[yl(p+1)T yl(p+2)T … yl(p+f-1)T]T，yl(p +1) represents Yl TP +1 column, yl(p +2) represents Yl TP +2 column, yl(p + f-1) represents Yl TP + f-1 column (b);
s7.6, acquiring control performance monitoring quantity eta of the current working condition by using online dataCPI：
Figure FDA00030521101500000712
Wherein the content of the first and second substances,
Figure FDA00030521101500000713
indicating the control performance monitoring quantity, p, of the ith window of the on-line datak,lRepresentation matrix Xp,lThe kth column vector of (1);
s7.7, control Upper Limit UP according to control Performance monitoring amountCPIJudgment of
Figure FDA0003052110150000081
Whether or not less than UPCPIIf, if
Figure FDA0003052110150000082
Less than UPCPIThen the system performance is better, if
Figure FDA0003052110150000083
Greater than UPCPIThe system performance is degraded.
9. The method of claim 5, wherein the step S9 includes the steps of:
s9.1, obtaining a reference working condition model residual error by using online data:
Figure FDA0003052110150000084
wherein v isr,l(t) reference for the sampling instant of the ith window of the on-line dataThe residual error of the working condition model,
Figure FDA0003052110150000085
representing the reference condition disturbance model estimate, yl(t) an online data process output vector, u, at the tth sampling instant of the ith windowl(t) an online data process input vector at the tth sampling instant, G, for the ith windowm,l(q-1) Representing a process model used by the EWMA batch controller in the ith window of the online data;
s9.2, calculating model mismatching/disturbance dynamic combined monitoring quantity eta by utilizing online collected input data and output dataPDMI：
Figure FDA0003052110150000086
Wherein the content of the first and second substances,
Figure FDA0003052110150000087
has a value range of (0, 1)]，εr(t) represents the white noise estimate for the reference condition at the tth sampling instant, | cov (εr(t)) | denotes εr(t) determinant of covariance matrix, vr,l(t) represents the on-line data reference case model residual at the tth sampling instant of the ith window, | cov (v)r,l(t)) | denotes vr,l(t) determinant of covariance matrix;
s9.3, combining control upper limit UP of model mismatch/disturbance dynamic joint monitoring quantityPDMIThe factors for controlling the performance degradation are preliminarily diagnosed if the current working condition is met
Figure FDA0003052110150000088
Greater than UPPDMIIf the model mismatch or the external disturbance dynamics is changed, the system performance is reduced due to the model mismatch or the external disturbance dynamics; if the current operating condition is
Figure FDA0003052110150000089
Less than UPPDMIThen changes in the controller tuning parameters are a factor that causes degradation in system performance.
10. The method of claim 8, wherein the step S11 includes the steps of:
s11.1, calculating interference characteristic monitoring quantity eta of online dataDVI：
Figure FDA00030521101500000810
Wherein the content of the first and second substances,
Figure FDA00030521101500000811
representing the interference characteristic monitoring quantity of the online data corresponding to the ith window,
Figure FDA00030521101500000812
for characterizing the change of the white noise variance, ε, in the actual operating conditionsr(t) represents the white noise estimation value of the reference condition at the t-th sampling moment, epsilonh,l(t) represents the white noise estimation value of the actual working condition at the tth sampling moment of the ith window in the online monitoring stage, | cov (epsilon)h,l(t)) | denotes εh,lDeterminant of covariance matrix of (t), | cov (εr(t)) | denotes εr(t) determinant of covariance matrix;
s11.2, calculating the actual working condition model residual error of the online data:
Figure FDA0003052110150000091
wherein v ish,l(t) represents the model residual error of the online data actual condition corresponding to the ith window,
Figure FDA0003052110150000092
on-line data real condition disturbance representing the ith windowDynamic model estimate, yl(t) an online data process output vector, u, at the tth sampling instant of the ith windowl(t) an online data process input vector at the tth sampling instant, G, for the ith windowm,l(q-1) Representing a process model used by the EWMA batch controller in the ith window of the online data;
s11.3, calculating disturbance model monitoring quantity eta of online dataDMI：
Figure FDA0003052110150000093
Wherein the content of the first and second substances,
Figure FDA0003052110150000094
representing the monitored quantity, v, of the disturbance model corresponding to the l-th windowr,l(t) represents the on-line data reference condition model residual at the tth sampling instant of the ith window, vh,l(t) represents the on-line data real-time model residual error at the tth sampling time of the ith window, J (f (v)r,l(t)),g(vh,l(t))) represents vh,l(t) and vr,l(t) KL divergence;
s11.4, calculating model mismatch monitoring quantity eta of online dataPMI：
Figure FDA0003052110150000095
Wherein the content of the first and second substances,
Figure FDA0003052110150000096
model mismatch monitoring quantity, epsilon, of the on-line data corresponding to the ith window of the on-line datah,l(t) represents the white noise estimate of the actual condition at the tth sampling time of the ith window of the on-line data, vh,l(t) represents the real model residual at the tth sampling instant of the l windows of online data, cov (εh,l(t)) represents εh,l(t) determinant of covariance matrix, cov (v)h,l(t)) represents vh,l(t) determinant of covariance matrix;
s11.5, combining the lower control limit LP of the interference characteristic monitoring quantityDVIUpper limit UP of control of the monitored quantity of interference characteristicsDVILower control limit LP of disturbance model monitoring quantityDMIAnd upper control limit UP of model mismatch monitoring quantityPMIFurther diagnosing and identifying the factors causing the performance degradation of the EWMA batch-to-batch control system and a plurality of disturbance characteristic monitoring quantities
Figure FDA0003052110150000097
Greater than UPDVIOr less than LPDVIIf so, the white noise variance changes to cause the performance of the EWMA batch-to-batch control system to be reduced; if the model monitoring quantity is disturbed
Figure FDA0003052110150000098
Less than LPDMIIf so, the performance of the EWMA batch control system is reduced due to the change of the external disturbance model; monitoring quantity if model mismatching
Figure FDA0003052110150000099
Greater than UPPMIThen the model mismatch causes the EWMA inter-batch control system performance to degrade.
CN202110491048.7A 2021-05-06 2021-05-06 Control performance diagnosis method for semiconductor process batch control system Active CN113189967B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110491048.7A CN113189967B (en) 2021-05-06 2021-05-06 Control performance diagnosis method for semiconductor process batch control system

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110491048.7A CN113189967B (en) 2021-05-06 2021-05-06 Control performance diagnosis method for semiconductor process batch control system

Publications (2)

Publication Number Publication Date
CN113189967A true CN113189967A (en) 2021-07-30
CN113189967B CN113189967B (en) 2022-05-27

Family

ID=76983825

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110491048.7A Active CN113189967B (en) 2021-05-06 2021-05-06 Control performance diagnosis method for semiconductor process batch control system

Country Status (1)

Country Link
CN (1) CN113189967B (en)

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2006070689A1 (en) * 2004-12-28 2006-07-06 Tokyo Electron Limited Semiconductor manufacturing apparatus, abnormality detection in such semiconductor manufacturing apparatus, method for specifying abnormality cause or predicting abnormality, and recording medium wherein computer program for executing such method is recorded
TW200745985A (en) * 2006-03-15 2007-12-16 Omron Tateisi Electronics Co Process abnorality analyzing device and method thereof, and memory medium
CN101739472A (en) * 2008-11-18 2010-06-16 上海华虹Nec电子有限公司 Method for constructing and simulating MOSFET mismatch model
CN102361014A (en) * 2011-10-20 2012-02-22 上海大学 State monitoring and fault diagnosis method for large-scale semiconductor manufacture process
CN106971953A (en) * 2015-10-01 2017-07-21 格罗方德半导体公司 Manufacture the error detection method in processing procedure
CN108459579A (en) * 2018-02-02 2018-08-28 郑州轻工业学院 Semiconductor run-to-run process failure diagnosis method based on time series models coefficient
CN108536127A (en) * 2018-04-20 2018-09-14 华中科技大学 A kind of model mismatch diagnostic method of the multivariable control system of data-driven

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2006070689A1 (en) * 2004-12-28 2006-07-06 Tokyo Electron Limited Semiconductor manufacturing apparatus, abnormality detection in such semiconductor manufacturing apparatus, method for specifying abnormality cause or predicting abnormality, and recording medium wherein computer program for executing such method is recorded
TW200745985A (en) * 2006-03-15 2007-12-16 Omron Tateisi Electronics Co Process abnorality analyzing device and method thereof, and memory medium
CN101739472A (en) * 2008-11-18 2010-06-16 上海华虹Nec电子有限公司 Method for constructing and simulating MOSFET mismatch model
CN102361014A (en) * 2011-10-20 2012-02-22 上海大学 State monitoring and fault diagnosis method for large-scale semiconductor manufacture process
CN106971953A (en) * 2015-10-01 2017-07-21 格罗方德半导体公司 Manufacture the error detection method in processing procedure
CN108459579A (en) * 2018-02-02 2018-08-28 郑州轻工业学院 Semiconductor run-to-run process failure diagnosis method based on time series models coefficient
CN108536127A (en) * 2018-04-20 2018-09-14 华中科技大学 A kind of model mismatch diagnostic method of the multivariable control system of data-driven

Also Published As

Publication number Publication date
CN113189967B (en) 2022-05-27

Similar Documents

Publication Publication Date Title
US7603328B2 (en) Dual-phase virtual metrology method
TWI427722B (en) Advanced process control system and method utilizing virtual metrology with reliance index and computer program product thereof
US20150371134A1 (en) Predicting circuit reliability and yield using neural networks
CN109389314B (en) Quality soft measurement and monitoring method based on optimal neighbor component analysis
Zhang et al. Fault detection strategy based on weighted distance of $ k $ nearest neighbors for semiconductor manufacturing processes
CN109409425B (en) Fault type identification method based on neighbor component analysis
CN113094893A (en) Wafer quality virtual measurement method and device, computer equipment and storage medium
CN111860933B (en) Predictive maintenance method for production machine assembly
CN114077876B (en) Strip steel hot continuous rolling multi-mode process monitoring method and device
CN112904810B (en) Process industry nonlinear process monitoring method based on effective feature selection
CN115994337B (en) Method and device for detecting minor faults in non-stationary process of hot continuous rolling of strip steel
CN111638707A (en) Intermittent process fault monitoring method based on SOM clustering and MPCA
CN111709181B (en) Method for predicting fault of polyester filament yarn industrial production process based on principal component analysis
CN111126255A (en) Numerical control machine tool cutter wear value prediction method based on deep learning regression algorithm
CN106845825B (en) Strip steel cold rolling quality problem tracing and control method based on improved PCA
CN111914889A (en) Rectifying tower abnormal state identification method based on brief kernel principal component analysis
CN113189967B (en) Control performance diagnosis method for semiconductor process batch control system
CN110751217A (en) Equipment energy consumption ratio early warning analysis method based on principal component analysis
CN113011102B (en) Multi-time-sequence-based Attention-LSTM penicillin fermentation process fault prediction method
Li et al. Partial Domain Adaptation in Remaining Useful Life Prediction With Incomplete Target Data
CN111914471A (en) Rectification tower fault detection method based on rapid nuclear independent component analysis
Wang et al. Nonlinear fault detection based on an improved kernel approach
CN108268987B (en) Method for estimating quality of various products
CN113420815B (en) Nonlinear PLS intermittent process monitoring method of semi-supervision RSDAE
Mengzhou et al. A modeling method for monitoring tool wear condition based on adaptive dynamic non-bias least square support vector machine

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant