Recoding Matlab code of MR-GRAS to R
REFERENCES : Temursho , U., Oosterhaven , J. and M.A. Cardenete (2019) , A multi - regional generalized RAS updating technique , IOpedia Research Paper No. 2, September 2019 , www. IOpedia .eu
#REFERENCES article : Temursho , U., Oosterhaven , J. and M.A. Cardenete (2019) , A multi - regional generalized RAS updating technique , IOpedia Research Paper No. 2, September 2019 , www. IOpedia .eu
#M 1 function [X,r,s,T] = mrgras (X0 ,u,v,G,Q,W, eps )
mrgras <- function(X0 ,u,v,G,Q,W, eps=0.00001,theshld=5000) {
#M 2 % PURPOSE : estimate a new multi - regional (or any partitioned ) matrix X as
#M 3 % close as possible to a given matrix X0 subject to the row sums , column
#M 4 % sums and non - overlapping aggregation constraints , using MR - GRAS approach
#M 5 % -------------------------------------------------------------------------
#M 6 % USAGE :
#M 7 % Write X = mrgras (X0 ,u,v,G,Q,W) OR [X,r,s,T] = mrgras (X0 ,u,v,G,Q,W) with
#M 8 % or without eps as the seventh optional input argument , where
#M 9 %
#M 10 % INPUT :
#M 11 % -> X0 = benchmark ( original ) matrix , not necessarily square
#M 12 % -> u = column vector of row totals (new row sums )
#M 13 % -> v = column vector of column totals ( new column sums )
#M 14 % -> W = non - overlapping aggregation constraints matrix
#M 15 % -> G & Q = the row and column aggregator matrices such that G*X0*Q = W;
#M 16 % non - overlapping aggregation necessarily implies that the
#M 17 % column sums of G and the row sums of Q must be all unity
#M 18 % -> eps = convergence tolerance level ; if empty , the default threshold
#M 19 % level is 0.1e -5 (=0.000001)
#M 20 % -> In case of missing w_IJ , set the corresponding missing number to
#M 21 % w_IJ = 1010101 ( assuming that 1010101 is not among the w_IJ 's)
#M 22 %
#M 23 % OUTPUT (input - output analysis - related interpretation ):
#M 24 % -> X = updated / balanced / adjusted / projected matrix
#M 25 % -> r = substitution effects ( row multipliers )
#M 26 % -> s = fabrication effects ( column multipliers )
#M 27 % -> T = technology or regional effects ( aggregation multipliers )
#M 28 % -------------------------------------------------------------------------
#M 29 % REFERENCES :
#M 30 % Temursho , U., Oosterhaven , J. and M.A. Cardenete (2019) , A multi - regional
#M 31 % generalized RAS updating technique , IOpedia Research Paper No. 2,
#M 32 % September 2019 , www. IOpedia .eu
#M 33 % -------------------------------------------------------------------------
#M 34 % NOTE FROM THE AUTHOR : Using this program and publishing the results in
#M 35 % the form of a report , working / discussion paper , journal article , etc .
#M 36 % requires citation of the above reference paper .
#M 37 % -------------------------------------------------------------------------
#M 38 % Written by: Umed Temursho , May 2019
#M 39 % E- mail : utemursho@gmail . com
#M 40
#M 41 [m,n] = size (X0);
m<-nrow(X0)
n<-ncol(X0)
#M 42 [h,k] = size (W);
h<-nrow(W)
k<-ncol(W)
#M 43 N = zeros (m,n);
N<-matrix(0,m,n)
#M 44 N(X0 <0) = -X0(X0 <0) ;
N[X0 <0]<--X0[X0<0]
#M 45 N = sparse (N); % could save memory with large - scale matrices
#M 46 P = X0+N;
P<-X0+N
#M 47 P = sparse (P);
#M 48 %
#M 49 r0 = ones (m ,1) ; % initial guess for r in step 0
r0<-matrix(1,m,1)
r0<-as.vector(r0)
#M 50 T0 = ones (h,k); % initial guess for T in step 0
T0<-matrix(1,h,k)
#M 51 Te = G '* T0*Q '; % T expanded to fit the dimention of X0
Te<-t(G)%*%T0%*%t(Q)
#M 52 p_rt = (P.* Te) '*r0;
p_rt<-t(P*Te)%*%r0
p_rt<-as.vector(p_rt)
#M 53 n_rt = (N.* invM (Te)) '* invd (r0)* ones (m ,1);
OneM<-matrix(1,m,1)
OneM<-as.vector(OneM)
n_rt <-t(N*invM(Te))%*%invd(r0)%*%OneM
n_rt<-as.vector(n_rt)
#M 54 s1 = invd (2* p_rt )*(v+ sqrt (v .^2+4* p_rt .* n_rt )); % first step s
s1<-invd(2*p_rt)%*%(v+sqrt(v^2+4*p_rt*n_rt))
s1<-as.vector(s1)
#M 55 ss = -invd (v)* n_rt ;
ss<--invd(v)%*%n_rt
ss<-as.vector(ss)
#M 56 s1( p_rt ==0) = ss( p_rt ==0) ;
s1[p_rt==0]<-ss[p_rt==0]
#M 57 %
#M 58 p_st = (P.* Te)*s1;
p_st<-(P*Te)%*%s1
p_st<-as.vector(p_st)
#M 59 n_st = (N.* invM (Te))* invd (s1)* ones (n ,1) ;
OneN<-matrix(1,n,1)
OneN<-as.vector(OneN)
# n_st<-(N*invM(Te))*invd(s1)*vect_1_n Version origine
n_st<-(N*invM(Te))%*%invd(s1)%*%OneN
n_st<-as.vector(n_st)
#M 60 r1 = invd (2* p_st )*(u+ sqrt (u .^2+4* p_st .* n_st )); % first step r
r1<-invd(2*p_st)%*%(u+sqrt(u^2+4*p_st*n_st ))
r1<-as.vector(r1)
#M 61 rr = -invd (u)* n_st ;
rr<--invd(u)%*%n_st
rr<-as.vector(rr)
#M 62 r1( p_st ==0) = rr( p_st ==0) ;
r1[p_st==0]<-rr[p_st==0]
#M 63 %
#M 64 P_rs = G*( diag (r1)*P* diag (s1))*Q;
P_rs<-G%*%(diag(r1)%*%P%*%diag(s1))%*%Q
#M 65 N_rs = G*( invd (r1)*N* invd (s1))*Q;
N_rs<-G%*%(invd(r1)%*%N%*%invd(s1))%*%Q
#M 66 T1 = invM (2* P_rs ) .*( W+ sqrt (W .^2+4* P_rs .* N_rs )); % first step T
T1<-invM(2*P_rs)*(W+sqrt(W^2+4*P_rs*N_rs))
#M 67 TT = -invM (W).* N_rs ;
TT<--invM(W)*N_rs
#M 68 T1( P_rs ==0) = TT( P_rs ==0) ;
T1[P_rs==0]<-TT[P_rs==0]
#M 69 T1(W ==1010101) = 1; % set t_IJ =1 for missing aggregation total w_IJ
T1[W ==1010101]<-1
#M 70 Te = G '* T1*Q ';
T1<-as.matrix(T1)
Te<-t(G)%*%T1%*%t(Q)
#M 71 %
#M 72 p_rt = (P.* Te) '*r1;
p_rt<-t(P*Te)%*%r1
p_rt<-as.vector(p_rt)
#M 73 n_rt = (N.* invM (Te)) '* invd (r1)* ones (m ,1);
n_rt<-t(N*invM(Te))%*%invd(r1)%*%OneM
n_rt<-as.vector(n_rt)
#M 74 s2 = invd (2* p_rt )*(v+ sqrt (v .^2+4* p_rt .* n_rt )); % second step s
s2<-invd(2*p_rt)%*%(v+sqrt(v^2+4*p_rt*n_rt ))
s2<-as.vector(s2)
#M 75 ss = -invd (v)* n_rt ;
ss<--invd(v)%*%n_rt
ss<-as.vector(ss)
#M 76 s2( p_rt ==0) = ss( p_rt ==0) ;
s2[p_rt==0]<-ss[p_rt==0]
#M 77 %
#M 78 p_st = (P.* Te)*s2;
p_st<-(P*Te)%*%s2
p_st<-as.vector(p_st)
#M 79 n_st = (N.* invM (Te))* invd (s2)* ones (n ,1) ;
n_st<-(N*invM(Te))%*%invd(s2)%*%OneN
n_st<-as.vector(n_st)
#M 80 r2 = invd (2* p_st )*(u+ sqrt (u .^2+4* p_st .* n_st )); % second step r
r2<-invd(2*p_st)%*%(u+sqrt(u^2+4*p_st*n_st ))
r2<-as.vector(r2)
#M 81 rr = -invd (u)* n_st ;
rr<--invd(u)%*%n_st
rr<-as.vector(rr)
#M 82 r2( p_st ==0) = rr( p_st ==0) ;
r2[p_st==0]<-rr[p_st==0]
#M 83 %
#M 84 P_rs = G*( diag (r2)*P* diag (s2))*Q;
P_rs<-G%*%(diag(r2)%*%P%*%diag(s2))%*%Q
#M 85 N_rs = G*( invd (r2)*N* invd (s2))*Q;
N_rs<-G%*%(invd(r2)%*%N%*%invd(s2))%*%Q
#M 86 T2 = invM (2* P_rs ) .*( W+ sqrt (W .^2+4* P_rs .* N_rs )); % second step T
T2<-invM(2*P_rs)*(W+sqrt(W^2+4*P_rs*N_rs ))
#M 87 TT = -invM (W).* N_rs ;
TT<--invM(W)*N_rs
#M 88 T2( P_rs ==0) = TT( P_rs ==0) ;
T2[P_rs==0]<-TT[P_rs==0]
#M 89 T2(W ==1010101) = 1; % set t_IJ =1 for missing w_IJ
T2[W==1010101]<-1
#M 90 %
#M 91 tmax = max ( max ( abs (T2 -T1)));
tmax<-max(abs(T2-T1))
#M 92 dif = [s2 -s1;r2 -r1; tmax ];
dif<-rbind(c(s2-s1),c(r2-r1),tmax)
#M 93 iter = 1 % first iteration
iter<-1
#M 94 if nargin < 7 || isempty ( eps)
#M 95 eps = 0.1e -5; % default tolerance level
#M 96 end
if(nargs()<7 || is.null(eps)){eps<-0.000001}
#M 97 M = max ( abs ( dif ));
M<-max(abs(dif))
#M 98 while (M > eps )
iter=0
while(M>eps & iter < theshld){ # Default conditions : M>eps or 5000 iterations
#M 99 s1 = s2;
s1<-s2
#M 100 r1 = r2;
r1<-r2
#M 101 T1 = T2;
T1<-T2
#M 102 Te = G '* T1*Q ';
T1<-as.matrix(T1)
Te<-t(G)%*%T1%*%t(Q)
#M 103 %
#M 104 p_rt = (P.* Te) '*r1;
p_rt<-t(P*Te)%*%r1
p_rt<-as.vector(p_rt)
#M 105 n_rt = (N.* invM (Te)) '* invd (r1)* ones (m ,1);
n_rt<-t(N*invM(Te))%*%invd(r1)%*%OneM
n_rt<-as.vector(n_rt)
#M 106 s2 = invd (2* p_rt )*(v+ sqrt (v .^2+4* p_rt .* n_rt )); % next step s
s2<-invd(2*p_rt)%*%(v+sqrt(v^2+4*p_rt*n_rt))
s2<-as.vector(s2)
#M 107 ss = -invd (v)* n_rt ;
ss<--invd(v)%*%n_rt
ss<-as.vector(ss)
#M 108 s2( p_rt ==0) = ss( p_rt ==0) ;
s2[p_rt==0]<-ss[p_rt==0]
#M 109 %
#M 110 p_st = (P.* Te)*s2;
p_st<-(P*Te)%*%s2
p_st<-as.vector(p_st)
#M 111 n_st = (N.* invM (Te))* invd (s2)* ones (n ,1) ;
n_st<-(N*invM(Te))%*%invd(s2)%*%OneN
n_st<-as.vector(n_st)
#M 112 r2 = invd (2* p_st )*(u+ sqrt (u .^2+4* p_st .* n_st )); % next step r
r2<-invd(2*p_st)%*%(u+sqrt(u^2+4*p_st* n_st))
r2<-as.vector(r2)
#M 113 rr = -invd (u)* n_st ;
rr<--invd(u)%*%n_st
rr<-as.vector(rr)
#M 114 r2( p_st ==0) = rr( p_st ==0) ;
r2[p_st==0]<-rr[p_st==0]
#M 115 %
#M 116 P_rs = G*( diag (r2)*P* diag (s2))*Q;
P_rs<-G%*%(diag(r2)%*%P%*%diag(s2))%*%Q
#M 117 N_rs = G*( invd (r2)*N* invd (s2))*Q;
N_rs<-G%*%(invd(r2)%*%N%*%invd(s2))%*%Q
#M 118 T2 = invM (2* P_rs ) .*( W+ sqrt (W .^2+4* P_rs .* N_rs )); % next step T
T2<-invM(2*P_rs)*(W+sqrt(W^2+4*P_rs*N_rs))
#M 119 TT = -invM (W).* N_rs ;
TT<--invM(W)*N_rs
#M 120 T2( P_rs ==0) = TT( P_rs ==0) ;
T2[P_rs==0]<-TT[P_rs==0]
#M 121 T2(W ==1010101) = 1; % set t_IJ =1 for missing w_IJ
T2[W==1010101]<-1
#M 122 %
#M 123 tmax = max ( max ( abs (T2 -T1)));
tmax<-max(abs(T2-T1))
#M 124 dif = [s2 -s1;r2 -r1; tmax ];
dif<-rbind(c(s2-s1),c(r2-r1),tmax)
#M 125 iter = iter +1
iter <- iter +1
#M 126 M = max ( abs ( dif ));
M<-max(abs(dif))
#M 127 end
}
print(paste0(iter," iterations required to ensure convergence (if equal to threshold : check convergence)"))
#M 128 s = s2; % final step s
s <- s2
#M 129 r = r2; % final step r
r <- r2
#M 130 T = T2; % final step T
Tfin <- T2
#M 131 Te = G '*T*Q ';
Tfin<-as.matrix(Tfin)
Te <- t(G)%*%Tfin%*%t(Q)
#M 132 %
#M 133 X = Te .*( diag (r)*P* diag (s))-invM (Te) .*( invd (r)*N* invd (s)); % updated matrix
X <- Te*(diag(r)%*%P%*%diag(s))-invM(Te)*(invd(r)%*%N%*%invd(s))
#M 134 end
return(list(X,r,s,Tfin))
}
#M 135
#M 136 function invd = invd (x) % auxiliary function used above
invd <- function(x) {
#M 137 invd = 1./ x;
sortie<-1/x
#M 138 invd (x ==0) = 1;
sortie[x==0]<-1
#M 139 invd = diag ( invd );
#M 140 end
return(diag(sortie))
}
#M 141
#M 142 function invM = invM (X) % auxiliary function used above
invM <- function(X) {
#M 143 invM = 1./ X;
sortie<-1/X
#M 144 invM (X ==0) = 1;
sortie[X==0]<-1
#M 145 end
return(sortie)
}
#M 146 % --------------------------------------------------------------------------
#M 147 % END OF THE CODE
#M 148 % -------------------------------------------------------------------------
Test of the code
We use the example given in the article (Updating with exhaustive constraints p.23)
datamatTEI<-c(63,9,14,9,-18,75,-14,53,-10,66,69,66,16,56,-21,9,93,-25,53,16,74,72,-1,80,4,-48,14,64,51,99,61,-1,84,6,16,27)
datamatUcible<-c(160,194,145,320,134,151)
datamatVcible<-c(197,71,151,242,178,265)
Wcible<-c(230,0,250,123,75,130,86,174,36)
datamatTEIupdatedTarget<-c(74.2,8.2,16.4,10.6,-21.5,72.1,-13.4,44.4,-10.4,68.5,52.8,52.2,18.8,64.8,-19.3,10.5,98.3,-28,61.7,14.5,85.5,83.5,-1.2,76,4,-59.6,12.9,63.9,37.5,75.3,51.7,-1.2,65.9,5.1,12.2,17.4)
Gagg<-cbind(diag(3),diag(3))
Qagg<-rbind(diag(3),diag(3))
Mcible<-matrix(Wcible,nrow=3,ncol=3,byrow=TRUE)
Xtarget<-matrix(datamatTEIupdatedTarget,nrow=6,ncol=6,byrow=TRUE)
print(Xtarget)
[,1] [,2] [,3] [,4] [,5] [,6]
[1,] 74.2 8.2 16.4 10.6 -21.5 72.1
[2,] -13.4 44.4 -10.4 68.5 52.8 52.2
[3,] 18.8 64.8 -19.3 10.5 98.3 -28.0
[4,] 61.7 14.5 85.5 83.5 -1.2 76.0
[5,] 4.0 -59.6 12.9 63.9 37.5 75.3
[6,] 51.7 -1.2 65.9 5.1 12.2 17.4
# With paper's names
X0<-matrix(datamatTEI,nrow=6,ncol=6,byrow=TRUE)
U<-as.vector(datamatUcible)
V<-as.vector(datamatVcible)
W<-matrix(Wcible,nrow=3,ncol=3,byrow=TRUE)
G<-Gagg
Q<-Qagg
WO<-G %*% X0 %*% Q
resultat<-mrgras(X0 ,U,V,G,Q,W, eps=0.0000001,theshld=1000)
[1] "12 iterations required to ensure convergence (if equal to threshold : check convergence)"
print(resultat[[1]])
[,1] [,2] [,3] [,4] [,5]
[1,] 74.245480 8.241533 16.36994 10.556140 -21.513217
[2,] -13.422031 44.359339 -10.39911 68.515197 52.766704
[3,] 18.777638 64.750609 -19.31884 10.512274 98.251786
[4,] 61.732934 14.480949 85.51896 83.465446 -1.209264
[5,] 4.012449 -59.633774 12.94707 63.894384 37.507731
[6,] 51.653535 -1.198657 65.88199 5.056554 12.196261
[,6]
[1,] 72.10012
[2,] 52.17990
[3,] -27.97346
[4,] 76.01098
[5,] 75.27214
[6,] 17.41032
The algorithm converge before the threshold, and we will check consistency and precision.
If the algorithm had’nt converge, you should investigate the convergence and check source data (is there any NA ou NaN or NULL or something like that ?).
Deviations at fine levels :
print(resultat[[1]]-Xtarget)
[,1] [,2] [,3] [,4]
[1,] 0.04547997 0.041532526 -0.0300577161 -0.043859806
[2,] -0.02203075 -0.040660515 0.0008870869 0.015197177
[3,] -0.02236207 -0.049390616 -0.0188439845 0.012273654
[4,] 0.03293418 -0.019050787 0.0189603394 -0.034554344
[5,] 0.01244911 -0.033774104 0.0470695859 -0.005615539
[6,] -0.04646537 0.001342976 -0.0180146936 -0.043446215
[,5] [,6]
[1,] -0.013217309 0.0001221746
[2,] -0.033296244 -0.0200971144
[3,] -0.048213688 0.0265375885
[4,] -0.009264431 0.0109752022
[5,] 0.007730863 -0.0278595585
[6,] -0.003738672 0.0103210897
The convergence is good at 1 decimal, wich is consistent with the target.
The algorithm reaches convergence at 10^-6 in 11 iterations, and in R it is near to be the same.
Deviations at agregated level :
Wtarget<-G %*% resultat[[1]] %*% Q
print(Wtarget)
[,1] [,2] [,3]
[1,] 230 -7.105427e-15 250
[2,] 123 7.500000e+01 130
[3,] 86 1.740000e+02 36
print(Wtarget-Mcible)
[,1] [,2] [,3]
[1,] 0.000000e+00 -7.105427e-15 0.000000e+00
[2,] 2.842171e-14 -1.421085e-14 2.842171e-14
[3,] 0.000000e+00 2.842171e-14 7.105427e-15
We find a very good precision on the 3rd dimension.
Consistency of the algorithm :
Mresult<-resultat[[1]]
print(rowSums(Mresult)-U)
[1] -1.556154e-07 -3.618795e-07 8.867970e-07 1.556154e-07
[5] 3.618795e-07 -8.867970e-07
print(colSums(Mresult)-V)
[1] 5.073102e-06 -5.193373e-07 6.179617e-07 -5.073102e-06
[5] 5.193373e-07 -6.179616e-07
print(sum(rowSums(Mresult)))
[1] 1104
print(sum(colSums(Mresult)))
[1] 1104
print(Mresult[1,1]+Mresult[4,1]+Mresult[1,4]+Mresult[4,4]) # we should find 230
[1] 230
Discrepancies un row sums and col sums are under the threshold, and even better on the 3rd dimension.
LS0tDQp0aXRsZTogIlJlY29kaW5nIE1hdGxhYiBjb2RlIG9mIE1SLUdSQVMgdG8gUiAtIEFsZXhhbmRyZSBCb3VyZ2VvaXMgMDEvMDgvMjAyMSINCm91dHB1dDoNCiAgcGRmX2RvY3VtZW50OiBkZWZhdWx0DQogIGh0bWxfbm90ZWJvb2s6IGRlZmF1bHQNCi0tLQ0KDQoqKlJlY29kaW5nIE1hdGxhYiBjb2RlIG9mIE1SLUdSQVMgdG8gUioqDQoNClJFRkVSRU5DRVMgOiBUZW11cnNobyAsIFUuLCBPb3N0ZXJoYXZlbiAsIEouIGFuZCBNLkEuIENhcmRlbmV0ZSAoMjAxOSkgLCBBIG11bHRpIC0gcmVnaW9uYWwgZ2VuZXJhbGl6ZWQgUkFTIHVwZGF0aW5nIHRlY2huaXF1ZSAsIElPcGVkaWEgUmVzZWFyY2ggUGFwZXIgTm8uIDIsIFNlcHRlbWJlciAyMDE5ICwgd3d3LiBJT3BlZGlhIC5ldQ0KDQpgYGB7cn0NCiNSRUZFUkVOQ0VTIGFydGljbGUgOiBUZW11cnNobyAsIFUuLCBPb3N0ZXJoYXZlbiAsIEouIGFuZCBNLkEuIENhcmRlbmV0ZSAoMjAxOSkgLCBBIG11bHRpIC0gcmVnaW9uYWwgZ2VuZXJhbGl6ZWQgUkFTIHVwZGF0aW5nIHRlY2huaXF1ZSAsIElPcGVkaWEgUmVzZWFyY2ggUGFwZXIgTm8uIDIsIFNlcHRlbWJlciAyMDE5ICwgd3d3LiBJT3BlZGlhIC5ldQ0KDQojTSAxIGZ1bmN0aW9uIFtYLHIscyxUXSA9IG1yZ3JhcyAoWDAgLHUsdixHLFEsVywgZXBzICkNCm1yZ3JhcyA8LSBmdW5jdGlvbihYMCAsdSx2LEcsUSxXLCBlcHM9MC4wMDAwMSx0aGVzaGxkPTUwMDApIHsNCiAgI00gMiAlIFBVUlBPU0UgOiBlc3RpbWF0ZSBhIG5ldyBtdWx0aSAtIHJlZ2lvbmFsIChvciBhbnkgcGFydGl0aW9uZWQgKSBtYXRyaXggWCBhcw0KICAjTSAzICUgY2xvc2UgYXMgcG9zc2libGUgdG8gYSBnaXZlbiBtYXRyaXggWDAgc3ViamVjdCB0byB0aGUgcm93IHN1bXMgLCBjb2x1bW4NCiAgI00gNCAlIHN1bXMgYW5kIG5vbiAtIG92ZXJsYXBwaW5nIGFnZ3JlZ2F0aW9uIGNvbnN0cmFpbnRzICwgdXNpbmcgTVIgLSBHUkFTIGFwcHJvYWNoDQogICNNIDUgJSAtLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tDQogICNNIDYgJSBVU0FHRSA6DQogICNNIDcgJSBXcml0ZSBYID0gbXJncmFzIChYMCAsdSx2LEcsUSxXKSBPUiBbWCxyLHMsVF0gPSBtcmdyYXMgKFgwICx1LHYsRyxRLFcpIHdpdGgNCiAgI00gOCAlIG9yIHdpdGhvdXQgZXBzIGFzIHRoZSBzZXZlbnRoIG9wdGlvbmFsIGlucHV0IGFyZ3VtZW50ICwgd2hlcmUNCiAgI00gOSAlDQogICNNIDEwICUgSU5QVVQgOg0KICAjTSAxMSAlIC0+IFgwID0gYmVuY2htYXJrICggb3JpZ2luYWwgKSBtYXRyaXggLCBub3QgbmVjZXNzYXJpbHkgc3F1YXJlDQogICNNIDEyICUgLT4gdSA9IGNvbHVtbiB2ZWN0b3Igb2Ygcm93IHRvdGFscyAobmV3IHJvdyBzdW1zICkNCiAgI00gMTMgJSAtPiB2ID0gY29sdW1uIHZlY3RvciBvZiBjb2x1bW4gdG90YWxzICggbmV3IGNvbHVtbiBzdW1zICkNCiAgI00gMTQgJSAtPiBXID0gbm9uIC0gb3ZlcmxhcHBpbmcgYWdncmVnYXRpb24gY29uc3RyYWludHMgbWF0cml4DQogICNNIDE1ICUgLT4gRyAmIFEgPSB0aGUgcm93IGFuZCBjb2x1bW4gYWdncmVnYXRvciBtYXRyaWNlcyBzdWNoIHRoYXQgRypYMCpRID0gVzsNCiAgI00gMTYgJSBub24gLSBvdmVybGFwcGluZyBhZ2dyZWdhdGlvbiBuZWNlc3NhcmlseSBpbXBsaWVzIHRoYXQgdGhlDQogICNNIDE3ICUgY29sdW1uIHN1bXMgb2YgRyBhbmQgdGhlIHJvdyBzdW1zIG9mIFEgbXVzdCBiZSBhbGwgdW5pdHkNCiAgI00gMTggJSAtPiBlcHMgPSBjb252ZXJnZW5jZSB0b2xlcmFuY2UgbGV2ZWwgOyBpZiBlbXB0eSAsIHRoZSBkZWZhdWx0IHRocmVzaG9sZA0KICAjTSAxOSAlIGxldmVsIGlzIDAuMWUgLTUgKD0wLjAwMDAwMSkNCiAgI00gMjAgJSAtPiBJbiBjYXNlIG9mIG1pc3Npbmcgd19JSiAsIHNldCB0aGUgY29ycmVzcG9uZGluZyBtaXNzaW5nIG51bWJlciB0bw0KICAjTSAyMSAlIHdfSUogPSAxMDEwMTAxICggYXNzdW1pbmcgdGhhdCAxMDEwMTAxIGlzIG5vdCBhbW9uZyB0aGUgd19JSiAncykNCiAgI00gMjIgJQ0KICAjTSAyMyAlIE9VVFBVVCAoaW5wdXQgLSBvdXRwdXQgYW5hbHlzaXMgLSByZWxhdGVkIGludGVycHJldGF0aW9uICk6DQogICNNIDI0ICUgLT4gWCA9IHVwZGF0ZWQgLyBiYWxhbmNlZCAvIGFkanVzdGVkIC8gcHJvamVjdGVkIG1hdHJpeA0KICAjTSAyNSAlIC0+IHIgPSBzdWJzdGl0dXRpb24gZWZmZWN0cyAoIHJvdyBtdWx0aXBsaWVycyApDQogICNNIDI2ICUgLT4gcyA9IGZhYnJpY2F0aW9uIGVmZmVjdHMgKCBjb2x1bW4gbXVsdGlwbGllcnMgKQ0KICAjTSAyNyAlIC0+IFQgPSB0ZWNobm9sb2d5IG9yIHJlZ2lvbmFsIGVmZmVjdHMgKCBhZ2dyZWdhdGlvbiBtdWx0aXBsaWVycyApDQogICNNIDI4ICUgLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLQ0KICAjTSAyOSAlIFJFRkVSRU5DRVMgOg0KICAjTSAzMCAlIFRlbXVyc2hvICwgVS4sIE9vc3RlcmhhdmVuICwgSi4gYW5kIE0uQS4gQ2FyZGVuZXRlICgyMDE5KSAsIEEgbXVsdGkgLSByZWdpb25hbA0KICAjTSAzMSAlIGdlbmVyYWxpemVkIFJBUyB1cGRhdGluZyB0ZWNobmlxdWUgLCBJT3BlZGlhIFJlc2VhcmNoIFBhcGVyIE5vLiAyLA0KICAjTSAzMiAlIFNlcHRlbWJlciAyMDE5ICwgd3d3LiBJT3BlZGlhIC5ldQ0KICAjTSAzMyAlIC0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0NCiAgI00gMzQgJSBOT1RFIEZST00gVEhFIEFVVEhPUiA6IFVzaW5nIHRoaXMgcHJvZ3JhbSBhbmQgcHVibGlzaGluZyB0aGUgcmVzdWx0cyBpbg0KICAjTSAzNSAlIHRoZSBmb3JtIG9mIGEgcmVwb3J0ICwgd29ya2luZyAvIGRpc2N1c3Npb24gcGFwZXIgLCBqb3VybmFsIGFydGljbGUgLCBldGMgLg0KICAjTSAzNiAlIHJlcXVpcmVzIGNpdGF0aW9uIG9mIHRoZSBhYm92ZSByZWZlcmVuY2UgcGFwZXIgLg0KICAjTSAzNyAlIC0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0NCiAgI00gMzggJSBXcml0dGVuIGJ5OiBVbWVkIFRlbXVyc2hvICwgTWF5IDIwMTkNCiAgI00gMzkgJSBFLSBtYWlsIDogdXRlbXVyc2hvQGdtYWlsIC4gY29tDQogICNNIDQwDQogICNNIDQxIFttLG5dID0gc2l6ZSAoWDApOw0KICBtPC1ucm93KFgwKQ0KICBuPC1uY29sKFgwKQ0KICAjTSA0MiBbaCxrXSA9IHNpemUgKFcpOw0KICBoPC1ucm93KFcpDQogIGs8LW5jb2woVykgIA0KICAjTSA0MyBOID0gemVyb3MgKG0sbik7DQogIE48LW1hdHJpeCgwLG0sbikNCiAgI00gNDQgTihYMCA8MCkgPSAtWDAoWDAgPDApIDsNCiAgTltYMCA8MF08LS1YMFtYMDwwXQ0KICAjTSA0NSBOID0gc3BhcnNlIChOKTsgJSBjb3VsZCBzYXZlIG1lbW9yeSB3aXRoIGxhcmdlIC0gc2NhbGUgbWF0cmljZXMNCiAgI00gNDYgUCA9IFgwK047DQogIFA8LVgwK04NCiAgI00gNDcgUCA9IHNwYXJzZSAoUCk7DQogICNNIDQ4ICUNCiAgI00gNDkgcjAgPSBvbmVzIChtICwxKSA7ICUgaW5pdGlhbCBndWVzcyBmb3IgciBpbiBzdGVwIDANCiAgcjA8LW1hdHJpeCgxLG0sMSkNCiAgcjA8LWFzLnZlY3RvcihyMCkNCiAgI00gNTAgVDAgPSBvbmVzIChoLGspOyAlIGluaXRpYWwgZ3Vlc3MgZm9yIFQgaW4gc3RlcCAwDQogIFQwPC1tYXRyaXgoMSxoLGspDQogICNNIDUxIFRlID0gRyAnKiBUMCpRICc7ICUgVCBleHBhbmRlZCB0byBmaXQgdGhlIGRpbWVudGlvbiBvZiBYMA0KICBUZTwtdChHKSUqJVQwJSoldChRKQ0KICAjTSA1MiBwX3J0ID0gKFAuKiBUZSkgJypyMDsNCiAgcF9ydDwtdChQKlRlKSUqJXIwDQogIHBfcnQ8LWFzLnZlY3RvcihwX3J0KQ0KDQogICNNIDUzIG5fcnQgPSAoTi4qIGludk0gKFRlKSkgJyogaW52ZCAocjApKiBvbmVzIChtICwxKTsNCiAgT25lTTwtbWF0cml4KDEsbSwxKQ0KICBPbmVNPC1hcy52ZWN0b3IoT25lTSkNCiAgbl9ydCA8LXQoTippbnZNKFRlKSklKiVpbnZkKHIwKSUqJU9uZU0NCiAgbl9ydDwtYXMudmVjdG9yKG5fcnQpDQogIA0KICAjTSA1NCBzMSA9IGludmQgKDIqIHBfcnQgKSoodisgc3FydCAodiAuXjIrNCogcF9ydCAuKiBuX3J0ICkpOyAlIGZpcnN0IHN0ZXAgcw0KICBzMTwtaW52ZCgyKnBfcnQpJSolKHYrc3FydCh2XjIrNCpwX3J0Km5fcnQpKSANCiAgczE8LWFzLnZlY3RvcihzMSkNCg0KICAjTSA1NSBzcyA9IC1pbnZkICh2KSogbl9ydCA7DQogIHNzPC0taW52ZCh2KSUqJW5fcnQNCiAgc3M8LWFzLnZlY3RvcihzcykNCiAgDQogICNNIDU2IHMxKCBwX3J0ID09MCkgPSBzcyggcF9ydCA9PTApIDsNCiAgczFbcF9ydD09MF08LXNzW3BfcnQ9PTBdDQogIA0KICAjTSA1NyAlDQogICNNIDU4IHBfc3QgPSAoUC4qIFRlKSpzMTsNCiAgcF9zdDwtKFAqVGUpJSolczENCiAgcF9zdDwtYXMudmVjdG9yKHBfc3QpDQogICNNIDU5IG5fc3QgPSAoTi4qIGludk0gKFRlKSkqIGludmQgKHMxKSogb25lcyAobiAsMSkgOw0KICBPbmVOPC1tYXRyaXgoMSxuLDEpDQogIE9uZU48LWFzLnZlY3RvcihPbmVOKQ0KICAjIG5fc3Q8LShOKmludk0oVGUpKSppbnZkKHMxKSp2ZWN0XzFfbiAgICBWZXJzaW9uIG9yaWdpbmUNCiAgbl9zdDwtKE4qaW52TShUZSkpJSolaW52ZChzMSklKiVPbmVOICANCiAgbl9zdDwtYXMudmVjdG9yKG5fc3QpDQogICNNIDYwIHIxID0gaW52ZCAoMiogcF9zdCApKih1KyBzcXJ0ICh1IC5eMis0KiBwX3N0IC4qIG5fc3QgKSk7ICUgZmlyc3Qgc3RlcCByDQogIHIxPC1pbnZkKDIqcF9zdCklKiUodStzcXJ0KHVeMis0KnBfc3Qqbl9zdCApKQ0KICByMTwtYXMudmVjdG9yKHIxKQ0KICAjTSA2MSByciA9IC1pbnZkICh1KSogbl9zdCA7DQogIHJyPC0taW52ZCh1KSUqJW5fc3QNCiAgcnI8LWFzLnZlY3RvcihycikNCiAgI00gNjIgcjEoIHBfc3QgPT0wKSA9IHJyKCBwX3N0ID09MCkgOw0KICByMVtwX3N0PT0wXTwtcnJbcF9zdD09MF0NCiAgDQogICNNIDYzICUNCiAgI00gNjQgUF9ycyA9IEcqKCBkaWFnIChyMSkqUCogZGlhZyAoczEpKSpROw0KICBQX3JzPC1HJSolKGRpYWcocjEpJSolUCUqJWRpYWcoczEpKSUqJVENCiAgI00gNjUgTl9ycyA9IEcqKCBpbnZkIChyMSkqTiogaW52ZCAoczEpKSpROw0KICBOX3JzPC1HJSolKGludmQocjEpJSolTiUqJWludmQoczEpKSUqJVENCiAgI00gNjYgVDEgPSBpbnZNICgyKiBQX3JzICkgLiooIFcrIHNxcnQgKFcgLl4yKzQqIFBfcnMgLiogTl9ycyApKTsgJSBmaXJzdCBzdGVwIFQNCiAgVDE8LWludk0oMipQX3JzKSooVytzcXJ0KFdeMis0KlBfcnMqTl9ycykpDQogICNNIDY3IFRUID0gLWludk0gKFcpLiogTl9ycyA7DQogIFRUPC0taW52TShXKSpOX3JzDQogICNNIDY4IFQxKCBQX3JzID09MCkgPSBUVCggUF9ycyA9PTApIDsNCiAgVDFbUF9ycz09MF08LVRUW1BfcnM9PTBdDQogICNNIDY5IFQxKFcgPT0xMDEwMTAxKSA9IDE7ICUgc2V0IHRfSUogPTEgZm9yIG1pc3NpbmcgYWdncmVnYXRpb24gdG90YWwgd19JSg0KICBUMVtXID09MTAxMDEwMV08LTENCiAgI00gNzAgVGUgPSBHICcqIFQxKlEgJzsNCiAgVDE8LWFzLm1hdHJpeChUMSkNCiAgVGU8LXQoRyklKiVUMSUqJXQoUSkgDQogICNNIDcxICUNCiAgI00gNzIgcF9ydCA9IChQLiogVGUpICcqcjE7DQogIHBfcnQ8LXQoUCpUZSklKiVyMQ0KICBwX3J0PC1hcy52ZWN0b3IocF9ydCkNCiAgI00gNzMgbl9ydCA9IChOLiogaW52TSAoVGUpKSAnKiBpbnZkIChyMSkqIG9uZXMgKG0gLDEpOw0KICBuX3J0PC10KE4qaW52TShUZSkpJSolaW52ZChyMSklKiVPbmVNIA0KICBuX3J0PC1hcy52ZWN0b3Iobl9ydCkNCiAgI00gNzQgczIgPSBpbnZkICgyKiBwX3J0ICkqKHYrIHNxcnQgKHYgLl4yKzQqIHBfcnQgLiogbl9ydCApKTsgJSBzZWNvbmQgc3RlcCBzDQogIHMyPC1pbnZkKDIqcF9ydCklKiUoditzcXJ0KHZeMis0KnBfcnQqbl9ydCApKQ0KICBzMjwtYXMudmVjdG9yKHMyKQ0KICAjTSA3NSBzcyA9IC1pbnZkICh2KSogbl9ydCA7DQogIHNzPC0taW52ZCh2KSUqJW5fcnQNCiAgc3M8LWFzLnZlY3RvcihzcykNCiAgI00gNzYgczIoIHBfcnQgPT0wKSA9IHNzKCBwX3J0ID09MCkgOw0KICBzMltwX3J0PT0wXTwtc3NbcF9ydD09MF0NCiAgDQogICNNIDc3ICUNCiAgI00gNzggcF9zdCA9IChQLiogVGUpKnMyOw0KICBwX3N0PC0oUCpUZSklKiVzMg0KICBwX3N0PC1hcy52ZWN0b3IocF9zdCkNCiAgI00gNzkgbl9zdCA9IChOLiogaW52TSAoVGUpKSogaW52ZCAoczIpKiBvbmVzIChuICwxKSA7DQogIG5fc3Q8LShOKmludk0oVGUpKSUqJWludmQoczIpJSolT25lTg0KICBuX3N0PC1hcy52ZWN0b3Iobl9zdCkNCiAgI00gODAgcjIgPSBpbnZkICgyKiBwX3N0ICkqKHUrIHNxcnQgKHUgLl4yKzQqIHBfc3QgLiogbl9zdCApKTsgJSBzZWNvbmQgc3RlcCByDQogIHIyPC1pbnZkKDIqcF9zdCklKiUodStzcXJ0KHVeMis0KnBfc3Qqbl9zdCApKQ0KICByMjwtYXMudmVjdG9yKHIyKQ0KICAjTSA4MSByciA9IC1pbnZkICh1KSogbl9zdCA7DQogIHJyPC0taW52ZCh1KSUqJW5fc3QNCiAgcnI8LWFzLnZlY3RvcihycikNCiAgI00gODIgcjIoIHBfc3QgPT0wKSA9IHJyKCBwX3N0ID09MCkgOw0KICByMltwX3N0PT0wXTwtcnJbcF9zdD09MF0NCiAgDQogICNNIDgzICUNCiAgI00gODQgUF9ycyA9IEcqKCBkaWFnIChyMikqUCogZGlhZyAoczIpKSpROw0KICBQX3JzPC1HJSolKGRpYWcocjIpJSolUCUqJWRpYWcoczIpKSUqJVENCiAgI00gODUgTl9ycyA9IEcqKCBpbnZkIChyMikqTiogaW52ZCAoczIpKSpROw0KICBOX3JzPC1HJSolKGludmQocjIpJSolTiUqJWludmQoczIpKSUqJVENCiAgI00gODYgVDIgPSBpbnZNICgyKiBQX3JzICkgLiooIFcrIHNxcnQgKFcgLl4yKzQqIFBfcnMgLiogTl9ycyApKTsgJSBzZWNvbmQgc3RlcCBUDQogIFQyPC1pbnZNKDIqUF9ycykqKFcrc3FydChXXjIrNCpQX3JzKk5fcnMgKSkNCiAgI00gODcgVFQgPSAtaW52TSAoVykuKiBOX3JzIDsNCiAgVFQ8LS1pbnZNKFcpKk5fcnMNCiAgI00gODggVDIoIFBfcnMgPT0wKSA9IFRUKCBQX3JzID09MCkgOw0KICBUMltQX3JzPT0wXTwtVFRbUF9ycz09MF0NCiAgI00gODkgVDIoVyA9PTEwMTAxMDEpID0gMTsgJSBzZXQgdF9JSiA9MSBmb3IgbWlzc2luZyB3X0lKDQogIFQyW1c9PTEwMTAxMDFdPC0xDQogIA0KICAjTSA5MCAlDQogICNNIDkxIHRtYXggPSBtYXggKCBtYXggKCBhYnMgKFQyIC1UMSkpKTsNCiAgdG1heDwtbWF4KGFicyhUMi1UMSkpDQogICNNIDkyIGRpZiA9IFtzMiAtczE7cjIgLXIxOyB0bWF4IF07DQogIGRpZjwtcmJpbmQoYyhzMi1zMSksYyhyMi1yMSksdG1heCkNCiAgI00gOTMgaXRlciA9IDEgJSBmaXJzdCBpdGVyYXRpb24NCiAgaXRlcjwtMSANCiAgI00gOTQgaWYgbmFyZ2luIDwgNyB8fCBpc2VtcHR5ICggZXBzKQ0KICAjTSA5NSBlcHMgPSAwLjFlIC01OyAlIGRlZmF1bHQgdG9sZXJhbmNlIGxldmVsDQogICNNIDk2IGVuZA0KICBpZihuYXJncygpPDcgfHwgaXMubnVsbChlcHMpKXtlcHM8LTAuMDAwMDAxfQ0KICAjTSA5NyBNID0gbWF4ICggYWJzICggZGlmICkpOw0KICBNPC1tYXgoYWJzKGRpZikpDQogIA0KICAjTSA5OCB3aGlsZSAoTSA+IGVwcyApDQogIGl0ZXI9MA0KICB3aGlsZShNPmVwcyAmIGl0ZXIgPCB0aGVzaGxkKXsgICMgRGVmYXVsdCBjb25kaXRpb25zIDogTT5lcHMgb3IgNTAwMCBpdGVyYXRpb25zDQogICAgI00gOTkgczEgPSBzMjsNCiAgICBzMTwtczINCiAgICAjTSAxMDAgcjEgPSByMjsNCiAgICByMTwtcjINCiAgICAjTSAxMDEgVDEgPSBUMjsNCiAgICBUMTwtVDINCiAgICAjTSAxMDIgVGUgPSBHICcqIFQxKlEgJzsNCiAgICBUMTwtYXMubWF0cml4KFQxKQ0KICAgIFRlPC10KEcpJSolVDElKiV0KFEpIA0KICAgICNNIDEwMyAlDQogICAgI00gMTA0IHBfcnQgPSAoUC4qIFRlKSAnKnIxOw0KICAgIHBfcnQ8LXQoUCpUZSklKiVyMQ0KICAgIHBfcnQ8LWFzLnZlY3RvcihwX3J0KQ0KICAgICNNIDEwNSBuX3J0ID0gKE4uKiBpbnZNIChUZSkpICcqIGludmQgKHIxKSogb25lcyAobSAsMSk7DQogICAgbl9ydDwtdChOKmludk0oVGUpKSUqJWludmQocjEpJSolT25lTQ0KICAgIG5fcnQ8LWFzLnZlY3RvcihuX3J0KQ0KICAgICNNIDEwNiBzMiA9IGludmQgKDIqIHBfcnQgKSoodisgc3FydCAodiAuXjIrNCogcF9ydCAuKiBuX3J0ICkpOyAlIG5leHQgc3RlcCBzDQogICAgczI8LWludmQoMipwX3J0KSUqJSh2K3NxcnQodl4yKzQqcF9ydCpuX3J0KSkNCiAgICBzMjwtYXMudmVjdG9yKHMyKQ0KICAgICNNIDEwNyBzcyA9IC1pbnZkICh2KSogbl9ydCA7DQogICAgc3M8LS1pbnZkKHYpJSolbl9ydA0KICAgIHNzPC1hcy52ZWN0b3Ioc3MpDQogICAgI00gMTA4IHMyKCBwX3J0ID09MCkgPSBzcyggcF9ydCA9PTApIDsNCiAgICBzMltwX3J0PT0wXTwtc3NbcF9ydD09MF0NCiAgICANCiAgICAjTSAxMDkgJQ0KICAgICNNIDExMCBwX3N0ID0gKFAuKiBUZSkqczI7DQogICAgcF9zdDwtKFAqVGUpJSolczINCiAgICBwX3N0PC1hcy52ZWN0b3IocF9zdCkNCiAgICAjTSAxMTEgbl9zdCA9IChOLiogaW52TSAoVGUpKSogaW52ZCAoczIpKiBvbmVzIChuICwxKSA7DQogICAgbl9zdDwtKE4qaW52TShUZSkpJSolaW52ZChzMiklKiVPbmVODQogICAgbl9zdDwtYXMudmVjdG9yKG5fc3QpDQogICAgI00gMTEyIHIyID0gaW52ZCAoMiogcF9zdCApKih1KyBzcXJ0ICh1IC5eMis0KiBwX3N0IC4qIG5fc3QgKSk7ICUgbmV4dCBzdGVwIHINCiAgICByMjwtaW52ZCgyKnBfc3QpJSolKHUrc3FydCh1XjIrNCpwX3N0KiBuX3N0KSkNCiAgICByMjwtYXMudmVjdG9yKHIyKQ0KICAgICNNIDExMyByciA9IC1pbnZkICh1KSogbl9zdCA7DQogICAgcnI8LS1pbnZkKHUpJSolbl9zdCANCiAgICBycjwtYXMudmVjdG9yKHJyKQ0KICAgICNNIDExNCByMiggcF9zdCA9PTApID0gcnIoIHBfc3QgPT0wKSA7DQogICAgcjJbcF9zdD09MF08LXJyW3Bfc3Q9PTBdDQogICAgDQogICAgI00gMTE1ICUNCiAgICAjTSAxMTYgUF9ycyA9IEcqKCBkaWFnIChyMikqUCogZGlhZyAoczIpKSpROw0KICAgIFBfcnM8LUclKiUoZGlhZyhyMiklKiVQJSolZGlhZyhzMikpJSolUQ0KICAgICNNIDExNyBOX3JzID0gRyooIGludmQgKHIyKSpOKiBpbnZkIChzMikpKlE7DQogICAgTl9yczwtRyUqJShpbnZkKHIyKSUqJU4lKiVpbnZkKHMyKSklKiVRDQogICAgI00gMTE4IFQyID0gaW52TSAoMiogUF9ycyApIC4qKCBXKyBzcXJ0IChXIC5eMis0KiBQX3JzIC4qIE5fcnMgKSk7ICUgbmV4dCBzdGVwIFQNCiAgICBUMjwtaW52TSgyKlBfcnMpKihXK3NxcnQoV14yKzQqUF9ycypOX3JzKSkNCiAgICAjTSAxMTkgVFQgPSAtaW52TSAoVykuKiBOX3JzIDsNCiAgICBUVDwtLWludk0oVykqTl9ycw0KICAgICNNIDEyMCBUMiggUF9ycyA9PTApID0gVFQoIFBfcnMgPT0wKSA7DQogICAgVDJbUF9ycz09MF08LVRUW1BfcnM9PTBdDQogICAgI00gMTIxIFQyKFcgPT0xMDEwMTAxKSA9IDE7ICUgc2V0IHRfSUogPTEgZm9yIG1pc3Npbmcgd19JSg0KICAgIFQyW1c9PTEwMTAxMDFdPC0xDQogICAgDQogICAgI00gMTIyICUNCiAgICAjTSAxMjMgdG1heCA9IG1heCAoIG1heCAoIGFicyAoVDIgLVQxKSkpOw0KICAgIHRtYXg8LW1heChhYnMoVDItVDEpKQ0KICAgICNNIDEyNCBkaWYgPSBbczIgLXMxO3IyIC1yMTsgdG1heCBdOw0KICAgIGRpZjwtcmJpbmQoYyhzMi1zMSksYyhyMi1yMSksdG1heCkNCiAgICAjTSAxMjUgaXRlciA9IGl0ZXIgKzENCiAgICBpdGVyIDwtIGl0ZXIgKzENCiAgICAjTSAxMjYgTSA9IG1heCAoIGFicyAoIGRpZiApKTsNCiAgICBNPC1tYXgoYWJzKGRpZikpDQogICAgI00gMTI3IGVuZA0KICB9DQogIHByaW50KHBhc3RlMChpdGVyLCIgaXRlcmF0aW9ucyByZXF1aXJlZCB0byBlbnN1cmUgY29udmVyZ2VuY2UgKGlmIGVxdWFsIHRvIHRocmVzaG9sZCA6IGNoZWNrIGNvbnZlcmdlbmNlKSIpKQ0KICAjTSAxMjggcyA9IHMyOyAlIGZpbmFsIHN0ZXAgcw0KICBzIDwtIHMyDQogICNNIDEyOSByID0gcjI7ICUgZmluYWwgc3RlcCByDQogIHIgPC0gcjINCiAgI00gMTMwIFQgPSBUMjsgJSBmaW5hbCBzdGVwIFQNCiAgVGZpbiA8LSBUMg0KICAjTSAxMzEgVGUgPSBHICcqVCpRICc7DQogIFRmaW48LWFzLm1hdHJpeChUZmluKQ0KICBUZSA8LSB0KEcpJSolVGZpbiUqJXQoUSkNCiAgI00gMTMyICUNCiAgI00gMTMzIFggPSBUZSAuKiggZGlhZyAocikqUCogZGlhZyAocykpLWludk0gKFRlKSAuKiggaW52ZCAocikqTiogaW52ZCAocykpOyAlIHVwZGF0ZWQgbWF0cml4DQogIFggPC0gVGUqKGRpYWcociklKiVQJSolZGlhZyhzKSktaW52TShUZSkqKGludmQociklKiVOJSolaW52ZChzKSkNCiAgI00gMTM0IGVuZA0KICByZXR1cm4obGlzdChYLHIscyxUZmluKSkNCn0NCiNNIDEzNQ0KI00gMTM2IGZ1bmN0aW9uIGludmQgPSBpbnZkICh4KSAlIGF1eGlsaWFyeSBmdW5jdGlvbiB1c2VkIGFib3ZlDQppbnZkIDwtIGZ1bmN0aW9uKHgpIHsNCiAgI00gMTM3IGludmQgPSAxLi8geDsNCiAgc29ydGllPC0xL3gNCiAgI00gMTM4IGludmQgKHggPT0wKSA9IDE7DQogIHNvcnRpZVt4PT0wXTwtMQ0KICAjTSAxMzkgaW52ZCA9IGRpYWcgKCBpbnZkICk7DQogICNNIDE0MCBlbmQNCiAgcmV0dXJuKGRpYWcoc29ydGllKSkNCn0NCiNNIDE0MQ0KI00gMTQyIGZ1bmN0aW9uIGludk0gPSBpbnZNIChYKSAlIGF1eGlsaWFyeSBmdW5jdGlvbiB1c2VkIGFib3ZlDQppbnZNIDwtIGZ1bmN0aW9uKFgpIHsNCiAgI00gMTQzIGludk0gPSAxLi8gWDsNCiAgc29ydGllPC0xL1gNCiAgI00gMTQ0IGludk0gKFggPT0wKSA9IDE7DQogIHNvcnRpZVtYPT0wXTwtMQ0KICAjTSAxNDUgZW5kDQogIHJldHVybihzb3J0aWUpDQp9DQojTSAxNDYgJSAtLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLQ0KI00gMTQ3ICUgRU5EIE9GIFRIRSBDT0RFDQojTSAxNDggJSAtLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tDQpgYGANCg0KKipUZXN0IG9mIHRoZSBjb2RlKioNCg0KV2UgdXNlIHRoZSBleGFtcGxlIGdpdmVuIGluIHRoZSBhcnRpY2xlIChVcGRhdGluZyB3aXRoIGV4aGF1c3RpdmUgY29uc3RyYWludHMgcC4yMykNCg0KYGBge3J9DQpkYXRhbWF0VEVJPC1jKDYzLDksMTQsOSwtMTgsNzUsLTE0LDUzLC0xMCw2Niw2OSw2NiwxNiw1NiwtMjEsOSw5MywtMjUsNTMsMTYsNzQsNzIsLTEsODAsNCwtNDgsMTQsNjQsNTEsOTksNjEsLTEsODQsNiwxNiwyNykNCmRhdGFtYXRVY2libGU8LWMoMTYwLDE5NCwxNDUsMzIwLDEzNCwxNTEpDQpkYXRhbWF0VmNpYmxlPC1jKDE5Nyw3MSwxNTEsMjQyLDE3OCwyNjUpDQpXY2libGU8LWMoMjMwLDAsMjUwLDEyMyw3NSwxMzAsODYsMTc0LDM2KQ0KZGF0YW1hdFRFSXVwZGF0ZWRUYXJnZXQ8LWMoNzQuMiw4LjIsMTYuNCwxMC42LC0yMS41LDcyLjEsLTEzLjQsNDQuNCwtMTAuNCw2OC41LDUyLjgsNTIuMiwxOC44LDY0LjgsLTE5LjMsMTAuNSw5OC4zLC0yOCw2MS43LDE0LjUsODUuNSw4My41LC0xLjIsNzYsNCwtNTkuNiwxMi45LDYzLjksMzcuNSw3NS4zLDUxLjcsLTEuMiw2NS45LDUuMSwxMi4yLDE3LjQpDQpHYWdnPC1jYmluZChkaWFnKDMpLGRpYWcoMykpDQpRYWdnPC1yYmluZChkaWFnKDMpLGRpYWcoMykpDQoNCk1jaWJsZTwtbWF0cml4KFdjaWJsZSxucm93PTMsbmNvbD0zLGJ5cm93PVRSVUUpDQpYdGFyZ2V0PC1tYXRyaXgoZGF0YW1hdFRFSXVwZGF0ZWRUYXJnZXQsbnJvdz02LG5jb2w9NixieXJvdz1UUlVFKQ0KcHJpbnQoWHRhcmdldCkNCg0KIyBXaXRoIHBhcGVyJ3MgbmFtZXMNClgwPC1tYXRyaXgoZGF0YW1hdFRFSSxucm93PTYsbmNvbD02LGJ5cm93PVRSVUUpDQpVPC1hcy52ZWN0b3IoZGF0YW1hdFVjaWJsZSkNClY8LWFzLnZlY3RvcihkYXRhbWF0VmNpYmxlKQ0KVzwtbWF0cml4KFdjaWJsZSxucm93PTMsbmNvbD0zLGJ5cm93PVRSVUUpDQpHPC1HYWdnDQpRPC1RYWdnDQpXTzwtRyAlKiUgWDAgJSolIFENCg0KcmVzdWx0YXQ8LW1yZ3JhcyhYMCAsVSxWLEcsUSxXLCBlcHM9MC4wMDAwMDAxLHRoZXNobGQ9MTAwMCkNCnByaW50KHJlc3VsdGF0W1sxXV0pDQoNCmBgYA0KDQpUaGUgYWxnb3JpdGhtIGNvbnZlcmdlIGJlZm9yZSB0aGUgdGhyZXNob2xkLCBhbmQgd2Ugd2lsbCBjaGVjayBjb25zaXN0ZW5jeSBhbmQgcHJlY2lzaW9uLg0KDQpJZiB0aGUgYWxnb3JpdGhtIGhhZCdudCBjb252ZXJnZSwgeW91IHNob3VsZCBpbnZlc3RpZ2F0ZSB0aGUgY29udmVyZ2VuY2UgYW5kIGNoZWNrIHNvdXJjZSBkYXRhIChpcyB0aGVyZSBhbnkgTkEgb3UgTmFOIG9yIE5VTEwgb3Igc29tZXRoaW5nIGxpa2UgdGhhdCA/KS4NCg0KRGV2aWF0aW9ucyBhdCBmaW5lIGxldmVscyA6DQpgYGB7cn0NCnByaW50KHJlc3VsdGF0W1sxXV0tWHRhcmdldCkNCmBgYA0KDQpUaGUgY29udmVyZ2VuY2UgaXMgZ29vZCBhdCAxIGRlY2ltYWwsIHdpY2ggaXMgY29uc2lzdGVudCB3aXRoIHRoZSB0YXJnZXQuIA0KDQpUaGUgYWxnb3JpdGhtIHJlYWNoZXMgY29udmVyZ2VuY2UgYXQgMTBeLTYgaW4gMTEgaXRlcmF0aW9ucywgYW5kIGluIFIgaXQgaXMgbmVhciB0byBiZSB0aGUgc2FtZS4NCg0KRGV2aWF0aW9ucyBhdCBhZ3JlZ2F0ZWQgbGV2ZWwgOg0KYGBge3J9DQpXdGFyZ2V0PC1HICUqJSByZXN1bHRhdFtbMV1dICUqJSBRDQpwcmludChXdGFyZ2V0KQ0KcHJpbnQoV3RhcmdldC1NY2libGUpDQpgYGANCg0KV2UgZmluZCBhIHZlcnkgZ29vZCBwcmVjaXNpb24gb24gdGhlIDNyZCBkaW1lbnNpb24uDQoNCkNvbnNpc3RlbmN5IG9mIHRoZSBhbGdvcml0aG0gOg0KYGBge3J9DQpNcmVzdWx0PC1yZXN1bHRhdFtbMV1dDQpwcmludChyb3dTdW1zKE1yZXN1bHQpLVUpDQpwcmludChjb2xTdW1zKE1yZXN1bHQpLVYpDQpwcmludChzdW0ocm93U3VtcyhNcmVzdWx0KSkpDQpwcmludChzdW0oY29sU3VtcyhNcmVzdWx0KSkpDQpwcmludChNcmVzdWx0WzEsMV0rTXJlc3VsdFs0LDFdK01yZXN1bHRbMSw0XStNcmVzdWx0WzQsNF0pICMgd2Ugc2hvdWxkIGZpbmQgMjMwDQpgYGANCg0KRGlzY3JlcGFuY2llcyB1biByb3cgc3VtcyBhbmQgY29sIHN1bXMgYXJlIHVuZGVyIHRoZSB0aHJlc2hvbGQsIGFuZCBldmVuIGJldHRlciBvbiB0aGUgM3JkIGRpbWVuc2lvbi4gIA0K