XProfan

Italia

p.specht

Für eine Funktion 3. Grades necessario das Verfahren nach R. Brent [1973] nur 8 Durchgänge, um auf 16-stellige Genauigkeit zu kommen. Allerdings sind Voraussetzungen zu erfüllen: Die Anfangswerte sollten unterschiedliche Vorzeichen haben. Dazu muss man schon etwas Ahnung haben, wie der Funktionsgraph ungefähr aussieht. Da das nachfolgende Programm hier nur eine Fehlermeldung abgibt, sich also nicht (wie die Profi-Programme) selber auf die Cerca nach geeigneten Werten begibt, handelt es sich erst einmal nur um eine Prinzipstudie. Immerhin funktioniert das Prinzip aber sehr gut!

Kompilieren Markieren Separieren

WindowTitle "Prinzipstudie: Nichtlineare Gleichung numerisch Lösen nach Brent [1973]"
WindowStyle 24:CLS:font 2
declare c0!,c1!,c2!,c3!,c4!,c5!,x!,a!,b!,c!,d!,e!,tol!,fa!,fb!,fc!
declare iter&,maxiter&,eps!,m!,r!,s!,t!,p!,q!,xs!
'Funktionstyp der zu lösenden Gleichung (Hier:Polynom max. 5.Grades)

proc f :parameters x!

    return c0!+x!*c1!+x!*x!*c2!+x!*x!*x!*c3!+x!*x!*x!*x!*c4!+x!*x!*x!*x!*x!*c5!

endproc

'Funktionsparameter: Testfunktion ist hier vom Grad 3 = Kubische Parabel
c0!=-5
c1!=8.17
c2!=-1.894
c3!=.114
c4!=0
c5!=0
'Gewünschte Genauigkeit (max. 15 signifikannte Stellen):
eps!=val("1e-16")
maxiter&=1000
t!=0'kleine additive Konstante (nnur bei sehr schwierigen Funktionen)
'Suchbereich 1:
print "\n Test-Sollwerte: ",6.208317, 9.67556
a!=6
b!=7
BRENT :Waitinput 1000
a!=7
b!=10
BRENT :beep:Waitinput
END

Proc BRENT

    Declare x!, c!,d!,e!,tol!,fa!,fb!,fc!
    Declare iter&,m!,r!,s!,p!,q!,xs!
    print "\n\n In den Grenzen ";format$("%g",a!),"bis",format$("%g",b!);" wurden gefunden:"
    fa!=f(a!)
    fb!=f(b!)

    if (fa!*fb!)>0

        print "\n ERROR: f(a) und f(b) sollten unterschiedliche Vorzeichen haben"
        beep:waitinput:end

    endif

    c!=a!:fc!=fa!' Zu Beginn ist c = a
    d!=b!-a!
    e!=d!
    iter&=0

    while iter&<maxiter&

        iter&=iter&+1

        if (fb!*fc!)>0

            c!=a!:fc!=fa!:d!=b!-a!:e!=d!

        endif

        if abs(fc!)<abs(fb!)

            a!=b!:b!=c!:c!=a!
            fa!=fb!: fb!=fc!: fc!=fa!

        endif

        tol!=2*eps!*abs(b!)+t!'Toleranz
        m!=(c!-b!)/2

        if (abs(m!)>tol!) AND (abs(fb!)>0)' Verfahren muss noch durchgeführt werden

            if (abs(e!)<tol!) OR (abs(fa!)<=abs(fb!))

                d!=m!:e!=m!

            else

                s!=fb!/fa!

                if a!=c!

                    p!=2*m!*s!
                    q!=1-s!

                else

                    q!=fa!/fc!
                    r!=fb!/fc!
                    p!=s!*(2*m!*q!*(q!-r!)-(b!-a!)*(r!-1))
                    q!=(q!-1)*(r!-1)*(s!-1)

                endif

                if p!>0

                    q!=-q!

                else

                    p!=-p!

                endif

                s!=e!
                e!=d!

                if ( (2*p!)<(3*m!*q!-abs(tol!*q!)) ) AND (p!<abs(s!*q!/2))

                    d!=p!/q!

                else

                    d!=m!
                    e!=m!

                endif

            endif

            a!=b!
            fa!=fb!

            if abs(d!)>tol!

                b!=b!+d!

            else

                if m!>0

                    b!=b!+tol!

                else

                    b!=b!-tol!

                endif

            endif

        else

            break

        endif

        fb!=f(b!)

    endwhile

    xs!=b!
    print "\n Nullwurzel: ";format$("%g",xs!)," mit "
    print "\n Restfehler: ";format$("%g",fb!)," gefunden in ";iter&;" Durchgängen."

Endproc