SVF filter saturation in feedback path / diode

[andy-cytomic]

Actually, that's not quite what I had in mind, rather I was suggesting you could only tabulate atanh(BP), because as soon as you substitute the result into the equation for HP the system becomes linear and can be solved in the usual way.

If you only tabulate ArcTanh(BP) then you would still need to compute the Tanh of that to get your BP which is needed for the other terms, and also state updates.

On second thought though, rather than tabulating atanh(BP) directly, you could tabulate a fully converged linearisation. The potential advantage is that if the tabulated result is slightly inaccurate (ie. we don't quite hit the linearisation point, but fairly close to it), then this might allow for a reduced error.

In fact, this leads to the approach I've been suggesting (in PMs at least) in the past as well. While I still haven't tried this in practice yet, in theory you could always add extra dimensions for each of the non-linearities (which you probably should be doing with NR as well), factor the system (eg. using LU) until you only have the non-linear sub-block left (same as with Newton), but then rather than iterating until convergence, fetch either the direct result or the fully converged linearisation for the whole non-linear block from a pre-computed MIMO table.

Additionally, if the table lookup is used simply as a short-cut for NR, then it becomes straight-forward to throw in an extra NR iteration or two to further refine the results if desired. This could allow for a lower resolution table, while still avoiding the pathological cases with NR by always providing a good initial guess that (hopefully) guarantees quadratic converge straight away.

In any case, most of this is just random thoughts really.

I have tried similar things in the past, but they didn't work that well for me. It was always the "bendy bits" of the approximations which needed the most help, and in that case unless you've got a really good approximation the first NR refinement could possibly go wrong and you're in an as bad if not worse position since you've just blown you cache doing the table lookup! I should probably look at it again though just to make sure.

[mystran]

Actually, that's not quite what I had in mind, rather I was suggesting you could only tabulate atanh(BP), because as soon as you substitute the result into the equation for HP the system becomes linear and can be solved in the usual way.

If you only tabulate ArcTanh(BP) then you would still need to compute the Tanh of that to get your BP which is needed for the other terms, and also state updates.

In the equation as given, simply looking up atanh(bp) and substituting it into the HP equation as constant makes the whole system linear, which means you can solve for HP and then BP=g*HP+s1.

[andy-cytomic]

In the equation as given, simply looking up atanh(bp) and substituting it into the HP equation as constant makes the whole system linear, which means you can solve for HP and then BP=g*HP+s1.

Gotcha. Then you're trading the interpolation cost for the cost of another division, which is probably better since as you point out your non-linearity can be a little off and you still solve the rest of the system properly. This is actually exactly what the DK-Method does, you solve for non-linear component contributions only (usually as currents), then with those in place solve the rest of the now linear circuit. The MIMO part comes in with multiple non-linear contributions.

[Marvinh]

Thank you guys for these post. Sorry if old.

I had to do a negative resonance 0 to negative 1

I saturated a+b*BP in a newton raphson loop

Then

BP became

{{bp->(a gk+gkeq)/(1-b gk),lp->-((g (a gk+gkeq))/(-1+b gk))+s2}}

HP i dont think should use arctanh(BP)

HP=IN-LP is safer for high gains but im not sure what type of HP that is or if its correct.

HP=IN-LP (is probably a 1 pole not 2)

Code: Select all

    template <typename FloatType>
    class NLSVF {
    public:
                
        NLSVF() {
            
        }
        void init(FloatType srIn){
            s1=0;
            s2=0;
            sr = srIn;
            isr = 1.0/sr;
        }
        double R = 0.0;
        void setParams(FloatType cutoffIn, FloatType resIn) {
            cutoffIn=std::clamp(cutoffIn, static_cast<FloatType>(10.0), static_cast<FloatType>(22049.0));
            FloatType w = juce::MathConstants<FloatType>::pi*cutoffIn*isr;
            FloatType co = sin(w*0.5);
            co = 1-2*(co*co);
            g = sin(w)/co;
            k = 2-2*resIn;
            R = -resIn;
        }
        
        FloatType process(FloatType IN) {
            

           
            FloatType HP = (IN-(g+k)*s1-s2)/(1+g*(g+k));
            FloatType BP = g*HP+s1;
            FloatType LP = g*BP+s2;

            
            FloatType a = (s1+g*(IN-s2))/g;
            FloatType b = (-1-g*(g+R))/g;
            
            for(int n = 0; n < 50; n++) {
                FloatType x = a+b*BP;
                FloatType fk = tanh(x);
                FloatType gk = 1-fk*fk;
                FloatType gkeq = fk-gk*x;
                
                FloatType nx = (a*gk+gkeq)/(1-b*gk);
                
                if(abs(nx-BP)<1e-8){
                    BP=nx;
                    break;
                }
                BP=nx;
                
            }
            LP = g*BP+s2;
            //this doesnt work well
//            HP = IN - k*BP - atanh(BP) - LP;
            //HP = IN-LP is better im not sure if correct;
            s1 = 2*BP-s1;
            s2 = 2*LP-s2;
            
            return LP;
            
        }
    private:
        FloatType s1{0};
        FloatType s2{0};
        FloatType sr{44100.0};
        FloatType isr{1.0/44100.0};
        FloatType g{0.0};
        FloatType k{0.0};
    };
    

[Urs]

Uhm, you're basically solving linear, then replacing BP with a possibly correctly solved solution, but then you integrate over the linear LP and non-linear BP. Hence, the result of the the non-linearity enters s2 at the wrong point.

(or do I need more coffee?)

Addendum: You just need to calculate LP after the loop and before the integration step (if the loop correctly calculates BP)

[andy-cytomic]

In any equations you want to add an implicit non-linearity you just replace f(x) with m*x + b, where m = f'(x), b = f(x) - m*x, then solve as normal. For a full reference on this look at my ADC2020 talk here: https://cytomic.com/technical-papers

[Marvinh]

Uhm, you're basically solving linear, then replacing BP with a possibly correctly solved solution, but then you integrate over the linear LP and non-linear BP. Hence, the result of the the non-linearity enters s2 at the wrong point.

(or do I need more coffee?)

Addendum: You just need to calculate LP after the loop and before the integration step (if the loop correctly calculates BP)

Thank you I fixed it but it still had the same result. very weird. Linear solution was for initial guess.

The state update of s1 might have kept LP from not saturating.

[Marvinh]

In any equations you want to add an implicit non-linearity you just replace f(x) with m*x + b, where m = f'(x), b = f(x) - m*x, then solve as normal. For a full reference on this look at my ADC2020 talk here: https://cytomic.com/technical-papers

Make sense where the a+b*bp comes from?