In[]:=
{r=22,l=210^-1,c=110^-4,vi=24,Δvi=0,dt=0.00006};d=DiagonalMatrix[{1,1}];a={{0,-1/l},{1/c,-1/(r*c)}};g=Inverse[a];b={1/l,0};s1=Inverse[d-dt*a+1/2*dt*dt*a.a];s2=Inverse[d-dt*a];star1=s1.g.(b*vi+g.b*Δvi/dt)-g.(b*(vi+Δvi)+g.b*Δvi/dt);star2=s2.g.(b*vi+g.b*Δvi/dt)-g.(b*(vi+Δvi)+g.b*Δvi/dt);il[t_]:=Simplify[(MatrixPower[s1,i]-MatrixPower[s2,i]).{0.,0.}+(d.(d-MatrixPower[s1,i]).Inverse[d-s1]).star1-(d.(d-MatrixPower[s2,i]).Inverse[d-s2]).star2][[1]]/.{i->t/dt}vc[t_]:=Simplify[(MatrixPower[s1,i]-MatrixPower[s2,i]).{0.,0.}+(d.(d-MatrixPower[s1,i]).Inverse[d-s1]).star1-(d.(d-MatrixPower[s2,i]).Inverse[d-s2]).star2][[2]]/.{i->t/dt}Plot[il[t],{t,0,0.08},AxesLabel->{"s","il[t]/A"},PlotLegends->{"0.00006"},PlotStyle->{Red},PlotRange->All]Plot[vc[t],{t,0,0.08},AxesLabel->{"s","vc[t]/V"},PlotLegends->{"0.00006"},PlotStyle->{Blue},PlotRange->All]
Out[]=
Out[]=
