Exponential curve fitting

We use linear least squares fitting to observe the exponential behaviour in the voltage of a capacitor as it charges or discharges in an RC circuit as a function of time.

The fitting of an exponential function can be accomplished through linear least squares (LLSQ) after taking the logarithm of the data.

> with(plots):

We define a linear least squares procedure: (developed in the worksheet DataAnalysis.mws )

> addme:=arg->evalf(Sum(arg[_i],_i=1..nops(arg)));

> addpr:=(arg1,arg2)->evalf(Sum(arg1[_i]*arg2[_i],_i=1..nops(arg1)));

> LLSQ:=proc(xv,yv) local _i,N,a,b,c,d,Delta,A,B,sq,sigsq,dA,dB;

> N:=nops(xv);

> if nops(yv)<N then RETURN(`Mismatch in sizes of xv and yv`); fi;

> c:=addpr(xv,xv);

> d:=addme(xv);

> a:=addpr(xv,yv);

> b:=addme(yv);

> Delta:=evalf(N*c-d^2);

> A:=(c*b-d*a)/Delta;

> B:=(N*a-d*b)/Delta;

> sq:=[]; _i:='_i': for _i from 1 to N do:

> sq:=[op(sq),evalf(yv[_i]-A-B*xv[_i])]: od:

> sigsq:=addpr(sq,sq)/(N-2);

> dA:=sigsq*c/Delta;

> dB:=evalf(N*sigsq)/Delta;

> print(`[A,B,sigma,dA,dB]`);

> [A,B,evalf(sqrt(sigsq),4),evalf(sqrt(dA),4),evalf(sqrt(dB),4)];

> end:

>

The time at which the capacitor voltage was measured (in seconds):

> Times:=[]:

> for i from 1 to 21 do:

> Time:=(i-1)*10: Times:=[op(Times),Time]: od:

> Times;

The voltage across a big electrolytic capacitor while being connected to 8.1 Volts through a 1000 Ohm resistor at these times was measured to be:

> V_Cchg:=[0.025,0.555,1.04,1.47,1.85,2.16,2.43,2.66,2.83,3.05,3.22,3.38,3.52,3.66,3.79,3.9,3.97,4.15,4.25,4.34,4.43];

The battery voltage was 8.1 Volts. The capacitor charged maximally to 7.4 Volts when fed through a 1000 Ohm resistor! A leaking current of 0.7 mA accounting for the voltage drop across the resistor in the steady-state regime was observed for the electrolytic capacitor. A non-ideal world after all...

We generate a graph of the data points. First for the charging regime, then for the discharge.

> i:='i': PL_c:=[seq([Times[i],V_Cchg[i]],i=1..nops(Times))]:

> P_c:=plot(PL_c,style=point,symbol=cross,color=red): display(P_c);

The capacitor was discharged through the same resistor. We reset time to zero and use the same time interval as before.

> V_Cdis:=[4.29,3.96,3.67,3.42,3.19,2.98,2.79,2.62,2.45,2.3,2.15,2.02,1.9,1.78,1.67,1.57,1.47,1.39,1.3,1.23,1.15];

> i:='i': PL_d:=[seq([Times[i],V_Cdis[i]],i=1..nops(Times))]:

> P_d:=plot(PL_d,style=point,symbol=cross,color=green):

> display(P_d);

The discharge regime is easy to analyze, we have an exponential decay. Taking the log of the data allows to extract the decay constant using LLSQ fitting:

> ln_of_V_Cdis:=map(ln,V_Cdis);

> i:='i': PL_dl:=[seq([Times[i],ln_of_V_Cdis[i]],i=1..nops(Times))]:

> P_dl:=plot(PL_dl,style=point,symbol=cross,color=blue):

> display(P_dl);

This looks perfectly linear, and thus we have confidence in the fit.

> i:='i': sol_d:=LLSQ(Times,ln_of_V_Cdis);

> tau:=-1/sol_d[2];

We discharged the capacitor with a 1 kOhm resistor.

> evalf(tau/1000);

> evalf(%*10^6);

The capacitor was actually rated at 110 000 microFarads. We graph the fit together with the original voltage measurements:

> y_d:=exp(sol_d[1]+sol_d[2]*t);

> P_dfit:=plot(y_d,t=0..200,color=black):

> display({P_d,P_dfit},scaling=unconstrained);

The charging regime is described by the expression

To analyze the charging data we need a map that first divides out the final voltage (the one reached for asymptotic times, which in our case is 7.4V and not the battery voltage due to the leakage) from the datapoints, and computes the difference to 1. Then we take the logarithm of the voltages. The real part in the map below is taken to avoid a wrong complex-valued choice of branch.

> myfun:=x->Re(ln(1-x/7.4));

> ln_of_V_Cchg:=map(myfun,V_Cchg);

> i:='i': PL_cl:=[seq([Times[i],ln_of_V_Cchg[i]],i=1..nops(Times))]:

> P_cl:=plot(PL_cl,style=point,symbol=cross,color=blue):

> display(P_cl,scaling=unconstrained);

We observe that the data points behave linearly for the first third, and then turn over to a different slope. Thus, we pick out a part of the data set for the fitting.

> nops(Times);

> Times1:=[seq(Times[i],i=1..7)]: ln_of_V_Cchg1:=[seq(ln_of_V_Cchg[i],i=1..7)]:

> sol_c:=LLSQ(Times1,ln_of_V_Cchg1);

> tau_c:=-1/sol_c[2];

> C_c:=evalf(tau_c/1000);

> evalf(C_c*10^6*mu*F);

The value of the capacitance is apparently quite close to the one obtained from the discharge cycle.

Exercise:

Use the uncertainty on the slopes for the charge and discharge cycles to determine whether the two determinations of the capacitance are consistent with each other.

Exercise:

Use Maple's built-in LLSQ procedure (use ?fit to determine its usage) to perform a fit for which the intercept is fixed as A =0 (only the slope is determined from the fit). Compare your results for the capacitance.

Systematic Error

If we had used all the data points from the charging cycle, we would have obtained a smaller slope in the linear fit, and as a result a substantially larger capacitance

(about 30% higher, but also a much larger uncertainty). To use such a fit blindly (without checking whether the data to be fitted look like a straight line with statistical scatter) would result in the introduction of a systematic error. The physical reason for this systematic error is the strange behaviour in the charge cycle, whereby the capacictor never charges up to the full voltage, and a leakage current persists in the asymptotic time regime. The capacitor we used was an old specimen (garage sale purchase), and is likely to be faulty, i.e., probably performed better when it was new.

We compare the fit result with the data over the full charge time that was observed.

> y_c:=7.4*(1-exp(-t/tau_c));

> P_cfit:=plot(y_c,t=0..200,color=black):

> display({P_cfit,P_c});

>