Re: [Gnucap-devel] discontinuities

[al davis]

On Thursday 16 January 2014, Felix Salfelder wrote:
> this is surprising. particulatly because the pulse patch was
> about zero rise time..

Zero rise time was the catylist and initial test case, but 
that's not all.  It has to do with timing of events around break 
points in general.  It still isn't final.  Every break point 
needs an explicit before and after step.  The original code put 
a step *approximately at* a break point, which could be exactly 
at, a little before, or a little after.  The modified code put a 
step a little before.  It also needs a step a little after, and 
must explicitly avoid exactly at.  More work is needed.

break point definition:  a signal or its first or second 
derivative with respect to either time or its input is 
discontinuous.

Note that I am using a different definition of break point than 
Spice does.

> . anyway i(r1) looks pretty wrong to
> me and the ringing in v(n2) is clearly a trap artifact.

clearly a trap artifact.

Try Gear.  The Gear artifacts are pretty wierd on this one too!

why is there a problem here?

Here's a stab at the answer.  

Trap is valid across break points.  It is common to believe that 
it is not valid, but it is valid.  I am not sure about Gear..

The analysis if Gear is more difficult.  It is a weighted 
average of 2 steps.  Higher order Gear uses more steps so it is 
even harder to properly analyze.  All of the books show a simple 
analysis where all steps are the same size.  In practice, often 
all steps are not the same size.  Gnucap tries to keep steps the 
same size when it makes sense to do so, but around break points 
steps are not the same size.  In the gnucap variable step Gear 
method, I did at one time prove that the method is still 
"consistent" and "convergent" but did not rigorously analyze its 
"stability".  (Quoted terms .. see a numerical analysis text for 
definitions.)  Based on this, Gear fallback is probably not a 
good idea.

The issue here is in time step control.  Although the method is 
valid through break points, the step control is not.   Say that 
again ....  Truncation error and curve fit time step control is 
not valid through break points.

For now, let's ignore the capacitors and just look at the step 
control when there are no storage devices.  Steps are chosen 
such that the result is a smooth curve, and that interpolating 
between actual calculation points gives accuracy within the 
specified tolerance.

The method is to try to fit a curve of order 1, 2, or 3 
depending on the option "trsteporder", which defaults to 3.

To estimate the error, take one more derivative.  That would be 
the 4th derivative for a 3rd order fit, or the second derivative 
for a first order fit.  Think of a Taylor series, estimating the 
error you get by truncating the series.  Truncation error.

In case it isn't obvious, a first order fit is line segments, a 
third order fit is cubic splines.  A higher order here allows 
larger time steps for the same interpolated accuracy.

The estimated value of the nth derivative here is calculated by 
divided differences.  3 points are needed to calculate the 
second derivative, used when trsteporder==1.  5 points are 
needed to calculate the 4th derivative, used when trsteporder=3.

If that derivative of a signal is discontinuous the 
estimation of truncation error is not valid.  The points before 
the discontinuity are not valid for this purpose.

In most cases, it is not a serious problem because the 
consequence is non-optimal time-stepping, not incorrect results.

In this case, large steps are used up to the break then suddenly 
very small steps are needed, but not demanded, hence the trap 
ringing.  Trap uses order 2 step control.  That requires 4 
samples to calculate the 3rd derivative.  Immediately after a 
break, 3 of the samples are before the break, when it was flat, 
so it still asks for large steps.

In this case, setting trsteporder to 1 eliminates the trap 
overshoot, because something else (the voltage source) now 
demands a small enough step.

So what is needed here is not method fallback (fallback to 
euler) but rather step control fallback.  If the step control 
falls back to 1 after a break, that should fix it.

Still working on it.  more to come.


Also ... while working on this, I think I fixed the "zero-time-
step" and "very backward time step" bug.  It has to do with 
rejecting a step then the quantizer unrejects it scheduling a 
new step at exactly the same time as the previous one.  More to 
follow on this too.