c++ Newtonverfahren mit dynamischen Arrays

floriandarter

Hallo,

für die Uni soll ich ein c++ Programm schreiben(mit Visual Studio 2010), mit dem die Nullstellen einer Polynomfunktion berechnet werden können. Der benutzer soll den Grad n des Polynoms eingeben und dann soll ein Array mit n+1 Feldern angelegt werden. Mein Problem ist, dass ich das Array nicht an die Funktion evalPolynom übergeben kann und auch nich weiß, ob die funktion das Polynom auswertet. Weiterhin sollen auch die ersten beiden Ableitungen in einem dynamischen Array abgespeichert werden. Hier der bisherige Code:

#include <iostream>
#include <cmath>
using namespace std;

double evalPolynom(double x, const double koeff[], int n)
{
int grad;
double fx=0;
for(int i=0; i<n; i++) {
grad = n-i-1;
fx = fx + koeff[i]*pow(x,grad);
}
return fx;
}

int main()
{
cout<<"Newton-Verfahren zur Nullstellenbestimmung "<<endl;
unsigned short n;
cout << "Grad der Funktion: ";
cin >> n;

double *koeff=new double[n+1];
double *koeffabl1=new double[n];
double *koeffabl2=new double[n-1];

for(int i=n; i>=0; i--)
{
cout<<"\nBitte geben sie a" << i << " ein: ";
cin>>koeff[i];

}
double x0,x=0,xn=0,xn1,fx,fx1=0,epsilon;
cout<<"Bitte geben Sie den Startwert ein:";
cin>>x0;

epsilon=10e-12;

xn1=x0;

while (abs(xn1-xn) > epsilon)
{
xn=xn1;
fx=evalPolynom(x,koeff,n);

xn1=xn-(fx/fx1);
}

//koeffabl1[n]=koeff[n]*n;
//n=n-1;

//...

delete [] koeff;

}

Ich freue mich über ein paar tipps und anregungen.
Vielen Dank.

SeppJ

Dynamisches Array heißt auf C++ std::vector. Damit dürften auch alle anderen beschriebenen Probleme entfallen, denn std::vector verhält sich, ganz im Gegensatz zu Arrays, ganz intuitiv bei der Parameterübergabe an Funktionen. Damit wären deine technischen Schwierigkeiten beseitigt und du kannst dich auf die eigentliche Aufgabe konzentrieren.

Oder ist das so einer dieser Kurse, wo alles genau so gemacht werden muss, wie der Prof es sagt, der C++ vor 30 Jahren anhand einer Übersetzungstabelle C <-> C++ gelernt hat?

floriandarter

Es soll schon so gemacht werden, wie ich angefangen habe, also mit dem new-Operator.

SeppJ

Hier hast du new[]. Gefahr von Speicherlöchern: Unmöglich. Standardbibliothek: Nicht benutzt. Pointergefrickel: Nicht mehr nötig.

namespace seppj
{
  typedef long unsigned int size_t;

  template<class InputIterator, class OutputIterator>
  OutputIterator copy ( InputIterator first, InputIterator last, OutputIterator result )
  {
    while (first!=last) *result++ = *first++;
    return result;
  }

  template <class T> void swap ( T& a, T& b )
  {
    T c(a); a=b; b=c;
  }

  template<class ElementType> class vector
  {
    size_t size_;
    ElementType * data_;
  public:
    vector(size_t s): size_(s), data_(new ElementType[size_]) {}
    size_t size() { return size_; }
    ElementType& operator[](size_t x) { return data_[x]; }
    ElementType const& operator[](size_t x) const { return data_[x]; }
    ~vector() { delete[] data_; }
    vector const& operator=(vector other) 
    {
      swap(size_, other.size_); 
      swap(data_, other.data_);
      return *this;
    }
    vector(vector const& other): size_(other.size_), data_(new ElementType[size_])
    {
      copy(other.data_, other.data_ + size_, data_);
    }
  };
}

//////////////////////////////////
//            DEMO              //
//////////////////////////////////

#include <iostream>
using namespace std;

void print(seppj::vector<int> arg)
{
  cout << "Dies ist eine Kopie: ";
  for (seppj::size_t i = 0; i < arg.size(); ++i) cout << arg[i] << ' ';
  cout << '\n';
}

void bar(seppj::vector<int> &arg)
{
  for (seppj::size_t i = 0; i < arg.size(); ++i) ++arg[i];
}

seppj::vector<int> foo()
{
  seppj::vector<int> daten(3);
  daten[0] = 99;
  daten[1] = 98;
  daten[2] = 1;
  return daten;
}

int main()
{
  seppj::vector<int> daten(4);
  daten[0] = 1;
  daten[1] = 5;
  daten[2] = 10;
  daten[3] = -963;

  cout << "Einfache Übergabe per Kopie. ";
  print(daten);

  bar(daten);
  cout << "Daten in Funktion geändert, dann kopiert. ";
  print(daten);

  cout << "Rückgabe als Kopie aus Funktion. ";
  print(foo());
}

Werner Salomon

floriandarter schrieb:

für die Uni soll ich ein c++ Programm schreiben(mit Visual Studio 2010), mit dem die Nullstellen einer Polynomfunktion berechnet werden können. Der benutzer soll den Grad n des Polynoms eingeben und dann soll ein Array mit n+1 Feldern angelegt werden. Mein Problem ist, dass ich das Array nicht an die Funktion evalPolynom übergeben kann und auch nich weiß, ob die funktion das Polynom auswertet.

Hallo floriandarter,

solche Art Aufgabenstellungen lesen zu müssen, ist ein Trauerspiel. Da wird gefordert, new zu verwenden, anscheinend bevor überhaupt was über std::vector end/oder std::string erzählt wird.

Die Aufgabestellung 'Nullstellen eines Polynoms mit Newton-Verfahren' schreit doch geradezu nach einem Funktionstemplate für das Newtonverfahren und einer Klasse für das Polynom. In unten stehenden Code habe ich das mal skizziert.

Deine Funktion 'evalPolynom' habe ich durch eine operator() ersetzt (Zeile 20-28). Dies hat den Vorteil, dass man ein Objekt der Klasse Polynom wie eine Funktion nutzen kann (siehe Zeile 68 und 70). Weiter benutze ich zur Berechnung des Polynoms das Horner-Schema - ich hoffe für Dich, dass das was Du in 'evalPolynom' codiert hast, nicht an Deiner Uni gelehrt wird.
Die Klasse Polynom hat noch eine Methode 'ableitung', mit der die Ableitung des Polynoms bestimmt und als Returnwert zurück gegeben wird (Zeile 35-44).
Und außerdem habe ich die Methode 'sub' hinzugefügt (Zeile 46-59). Diese reduziert das Polynom um eine Nullstelle, indem es das betreffende Objekt in ein Polynom umwandelt, welches das Ergebnis der Division Polynom/(x - x0) ist. Damit wird sicher gestellt, dass nicht eine Nullstelle mehrfach gefunden wird. Wird eine Nullstelle gefunden, so wird sie mittels 'sub' einfach aus dem Polynom entfernt (Zeile 96).

#include <iostream>
#include <vector>
#include <cassert>
#include <cmath> // std::abs
#include <stdexcept> // std::invalid_argument
#include <memory> // std::unique_ptr

class Polynom
{
public:
    Polynom() : m_koeff( 1 ) {}

    template< typename I >
    Polynom( I von, I bis )
        : m_koeff( von, bis )
    {
        assert( !m_koeff.empty() );
    }

    double operator()( double x ) const
    {
        // Berechnung des Funktionwertes an der Stelle 'x' nach dem
        //  http://de.wikipedia.org/wiki/Horner-Schema
        double ergebnis = 0;
        for( auto i = m_koeff.rbegin(); i != m_koeff.rend(); ++i )
            ergebnis = ergebnis * x + *i;
        return ergebnis;
    }

    int grad() const 
    {
        return int(m_koeff.size()) - 1;
    }

    Polynom ableitung() const
    {
        if( grad() <= 0 )
            return Polynom();
        std::vector< double > koeff;
        int faktor = 1;
        for( auto i = m_koeff.begin() + 1; i != m_koeff.end(); ++i, ++faktor )
            koeff.push_back( faktor * *i );
        return Polynom( koeff.begin(), koeff.end() );
    }

    double sub( double x0 )
    {   // --   dividiert das Polynom durch (x - x0) und git den Rest zurück
        assert( grad() > 0 );
        std::vector< double > koeff;
        auto i = m_koeff.rbegin();
        double rest = *i++;
        for( ; i != m_koeff.rend(); ++i )
        {
            koeff.push_back( rest );
            rest = x0 * rest + *i;
        }
        m_koeff.assign( koeff.rbegin(), koeff.rend() );
        return rest;
    }
private:
    std::vector< double > m_koeff;
};

template< typename Fkt, typename DFkt >
double newton( Fkt func, DFkt dfunc, double x, double eps )
{
    int check = 0;
    for( double y; std::abs( y = func( x ) ) > eps; ++check )
    {
        x -= y / dfunc( x );
        if( check > 100 ) // verhindert ggf. die Endlosschleife
        {   // Bem.: Exception ist hier vielleich übertrieben ...
            throw std::invalid_argument( "Newton Verfahren: x konvertiert nicht" );
        }
    }
    return x;
}

int main()
{
    using namespace std;
    cout<<"Newton-Verfahren zur Nullstellenbestimmung "<<endl;
    int n;
    cout << "Grad der Funktion: ";
    if( cin >> n && n > 0 ) // nur ein Grad des Polynoms >0 ist zulässig
    {
        unique_ptr< double[] > koeff( new double[n+1] ); // wenn schon 'new', dann mit smart-Pointer
        int i = 0;
        for( double a_i; i <= n && (cout<<"\nBitte geben sie a" << i << " ein: ", cin >> a_i); ++i )
            koeff[i] = a_i;
        if( cin ) // Eingabe war ok?
        {
            try
            {
                Polynom f( koeff.get(), koeff.get()+n+1 );
                for( double x0; f.grad() > 0; f.sub( x0 ) )
                    cout << "\nNullstelle: " << (x0 = newton( f, f.ableitung(), 0.0, 1.E-8 )) << endl;
            }
            catch( const std::exception& ex )
            {
                cerr << "\n### Exception: " << ex.what() << endl;
            }
        }
    }
    return 0;
}

.. an welcher Uni studierst Du, und was genau?

Gruß
Werner