On 23.06.2014, at 06:37, robert bristow-johnson <r...@audioimagination.com>
wrote:
> the other thing Urs brought up for discussion is an iterative and recursive
> process that converges on a result value, given an input. i am saying that
> this can be rolled out into a non-recursive equivalent, if the number of
> iterations needed for convergence is finite. this is a totally different
> issue.
Yes, it's a totally different issue, but I arrived there from the same
assumptions Andy did (actually, I think it was Andy who brought me there):
1. The filter algorithm can be derived from circuit analysis directly, without
taking a detour over an abstraction layer
2. If done so, non-linear elements of the circuit can be modeled at the right
place
Regarding the iterative method, unrolling like you did
> y0 = y[n-1]
> y1 = g * ( x[n] - tanh( y0 ) ) + s
> y2 = g * ( x[n] - tanh( y1 ) ) + s
> y3 = g * ( x[n] - tanh( y2 ) ) + s
> y[n] = y3
is *not* what I described in general. It's a subset won't ever converge in 3
iterations :^)
Thing is, while starting with y[n-1] is a possibility, it's not the only one
and in my experience it's hardly ever a good one.
A better solution is to start calculating the boundaries of the solution. In
this case it can be done easily:
y_max = g * ( x - (+1) ) + s;
y_min = g * ( x - (-1) + s;
... because +1 and -1 is the highest absolute value we'll ever get for tanh( y
).
Then you can immediately find the root of the linear equivalent by taking (
y_max + y_min ) / 2
This would be a better starting point, and it would *not* be derived from any
history.
The actual iteration is then something like this, which is an example of the
bisection method
while ( ( y_max - y_min ) > 0.000000001 )
{
y_test = ( y_max + y_min ) / 2
y_result = g * ( x + tanh( y_test ) ) + s;
if( y_result > y_test ) y_min = y_test;
else y_max = y_test;
}
http://en.wikipedia.org/wiki/Bisection_method
(not sure about getting the compare right, I need more coffee...)
To make this a filter of course we have to then also update s, which in case of
Euler or the trapezoidal rule can be done afterwards.
Cheers,
- Urs
--
dupswapdrop -- the music-dsp mailing list and website:
subscription info, FAQ, source code archive, list archive, book reviews, dsp
links
http://music.columbia.edu/cmc/music-dsp
http://music.columbia.edu/mailman/listinfo/music-dsp
