I have a Delay Differential Equation as given below

b = 1;
V0 = 300;
I0 = 0.001;
N1 = 7/4;
p1 = (\[Pi]*b*N1)/V0;
a1 = 1.001; a2 = 0.123; a3 = -3.622*10^-3; b1 = 0.001959; b2 = 0.031; \
b3 = 0.003241; G = 0.5*10^-5; \[CapitalOmega] = 0; C1 = p1/2000;
f = 1*10^3;
Is = 0.001;
\[Tau] = 10;
NL = 1000;

sol1[t_] = 
 NDSolve[{x'[t] - (I0*p1)/(
      2*C1)*((a1*x[t - \[Tau]] - G/(p1*I0)*x[t]) - 
        3/4*x[t - \[Tau]]^3*a2 - 5/8*a3*x[t - \[Tau]]^5) - 
     Is*p1*Cos[y[t]] == 0, 
   y'[t] - \[CapitalOmega] - (I0*p1)/(
      2*C1*x[t])*(b1*x[t - \[Tau]] + 3/4*x[t - \[Tau]]^3*b2 - 
        5/8*b3*x[t - \[Tau]]^5) + (Is*p1)/(2*C1*x[t])*Sin[y[t]] == 0, 
   x[0] == 0.003, y[0] == 0.001}, {x, y}, {t, 0, NL}]

I want to compute the Lyapunov exponent for this system by varying Tau from 3 to 22. For This I have used the following program.

b = 1;
V0 = 300;
I0 = 0.001;
N1 = 7/4;
p1 = (\[Pi]*b*N1)/V0;
a1 = 1.001; a2 = 0.123; a3 = -3.622*10^-3; b1 = 0.001959; b2 = 0.031; \
b3 = 0.003241; G = 0.5*10^-5; \[CapitalOmega] = 0; C1 = p1/2000;
f = 1*10^3;
Is = 0.001;
\[Tau] = 10;
NL = 1;

deq1 = (I0*p1)/(
   2*C1)*((a1*x1[t - \[Tau]] - G/(p1*I0)*x1[t]) - 
     3/4*x1[t - \[Tau]]^3*a2 - 5/8*a3*x1[t - \[Tau]]^5) + 
  Is*p1*Cos[y1[t]]; deq2 = \[CapitalOmega] + (I0*p1)/(
   2*C1*x1[t])*(b1*x1[t - \[Tau]] + 3/4*x1[t - \[Tau]]^3*b2 - 
     5/8*b3*x1[t - \[Tau]]^5) - (Is*p1)/(2*C1*x1[t])*
   Sin[y1[t]]; x10 = 0.003; y10 = 0.001; dx0 = 
 10^-8; tin = 0; tfin = 201; tstep = NL; acc = 12; lcedata = {};
sum = 0; d0 = Sqrt[(x10)^2 + (y10)^2]; For[i = 1, i < tfin/tstep, i++,
  sdeq = {x1'[t] == deq1, y1'[t] == deq2, x1[0] == x10, y1[0] == y10};
  sol = NDSolve[sdeq, {x1[t], y1[t]}, {t, 0, tstep}, 
   MaxSteps -> Infinity, Method -> "StiffnessSwitching", 
   PrecisionGoal -> acc, AccuracyGoal -> acc]; 
 xx1[t_] = x1[t] /. sol[[1]]; yy1[t_] = y1[t] /. sol[[1]]; 
 d1 = Sqrt[(xx1[tstep])^2 + (yy1[tstep])^2]; sum += Log10[d1/d0]; 
 dlce = sum/(tstep*i); AppendTo[lcedata, {tstep*i, Log10[dlce]}]; 
 w1 = (xx1[tstep])*(d0/d1); w2 = (yy1[tstep])*(d0/d1); 
 x10 = xx1[tstep]; y10 = yy1[tstep]; x20 = x10 + w1; y20 = y10 + w2; 
 i = i++;
 If[Mod[tstep*i, 1] == 0, 
  Print[" For t = ", tstep*i, " , ", " LCE = ", Log10[dlce]]]]
 S0 = ListLinePlot[{lcedata}, PlotRange -> Automatic, 
  AxesLabel -> {"t", "log10(LCE)"}, 
  AxesStyle -> Directive["Black", 13], GridLines -> Automatic]

This program gives me some values of the Lyapunov exponent for a given Tau (say Tau=10), but if I change the value of Tau the Lyapunov exponent values are not changing. Please help me to modify the program so that one can get different Lyapunov exponent values for different Tau.