Actual code:
Code:
//Newton-Raphson Method
//http://en.wikipedia.org/wiki/Newton's_method
#include <cstdlib>
#include <iostream>
#include <math.h>
using namespace std;
//------------------------------------------------------------------------------
class point
{
public:
double x1, x2;
double length()
{
return sqrt(pow(x1,2)+pow(x2,2)); //sqrt(x1^2+x2^2)
}
//I write that overloaded operator is a friend of this class.
friend ostream & operator<<(ostream & theStream, point & thePoint);
};
//overloading operators point*double, double*point, point+point, point-point
point operator*(point x, double a)
{
point result;
result.x1 = x.x1 * a;
result.x2 = x.x2 * a;
return result;
}
point operator*(double a, point x)
{
point result;
result.x1 = x.x1 * a;
result.x2 = x.x2 * a;
return result;
}
point operator+(point & thePoint1, point & thePoint2)
{
point sum;
sum.x1 = thePoint1.x1 + thePoint2.x1;
sum.x2 = thePoint1.x2 + thePoint2.x2;
return sum;
}
point operator-(point & thePoint1, point & thePoint2)
{
point sum;
sum.x1 = thePoint1.x1 - thePoint2.x1;
sum.x2 = thePoint1.x2 - thePoint2.x2;
return sum;
}
//overloaded operator <<
ostream & operator<<(ostream & theStream, point & thePoint)
{
theStream<<"["<<thePoint.x1<<", "<<thePoint.x2<<"]";
return theStream;
}
//------------------------------------------------------------------------------
class hessian
{
public: //I've got a problem with determinant (in main function there is
//a problem with compiling if I use 'private' here!!!!!!!!
double x00, x01, x10, x11;
//|x00 x01|
//|x10 x11|
double determinant;
public:
friend point operator*(hessian theHessian, point thePoint);
friend hessian operator*(hessian theHes1, hessian theHes2);
void det();
void invertHessian();
void hessianAtPoint(point thePoint);
friend ostream&operator<<(ostream & theStream, hessian & theHessian);
};
//internal functions - determinant, invert hessian and hessian at point
void hessian::det()
{ determinant = x00 * x11 - x10 * x01; }
void hessian::invertHessian()
{
//inverted matrix: A^(-1) = A^D / det(A)
//A = [(a, b), (c, d)], A^D = [(d, -b), (-c, a)]
double _x00, _x01, _x10, _x11;
_x00 = x11/determinant;
_x01 = -x01/determinant;
_x10 = x10/determinant;
_x11 = x00/determinant;
x00 = _x00; x01 = _x01; x10 = _x10; x11 = _x11;
}
void hessian::hessianAtPoint(point thePoint)
{
x00 = x11 = 2;
x01 = x10 = 0;
}
//overloading operators << and hessian*hessian
ostream&operator<<(ostream & theStream, hessian & theHessian)
{
theStream<<"|"<<theHessian.x00<<" "<<theHessian.x01<<"|]n";
theStream<<"|"<<theHessian.x10<<" "<<theHessian.x11<<"|]n";
}
hessian operator*(hessian theHes1, hessian theHes2)
{
hessian result;
//simple multiplication matrix by matrix
//[a b][e f] = [ae+bg af+bh]
//[c d][g h] [ce+dg cf+dh]
result.x00 = theHes1.x00 * theHes2.x00 + theHes1.x01 * theHes2.x10;
result.x01 = theHes1.x00 * theHes2.x01 + theHes1.x01 * theHes2.x11;
result.x10 = theHes1.x10 * theHes2.x00 + theHes1.x11 * theHes2.x10;
result.x11 = theHes1.x10 * theHes2.x01 + theHes1.x11 * theHes2.x11;
}
//------------------------------------------------------------------------------
class gradient
{
public: //???????????????????
point itsPoint;
public:
void gradientAtPoint(point thePoint);
friend ostream&operator<<(ostream & theStream, gradient & theGradient);
};
void gradient::gradientAtPoint(point thePoint)
{
itsPoint.x1 = 2 * thePoint.x1 - 1;
itsPoint.x2 = 2 * thePoint.x2 - 1;
}
ostream&operator<<(ostream & theStream, gradient & theGradient)
{
theStream<<theGradient.itsPoint;
return theStream;
}
//------------------------------------------------------------------------------
//overloading other operator: hessian*point
point operator*(hessian theHessian, point thePoint)
{
point result;
//[x00 x01] * [x1 x2] = [x00*x1 + x01*x2]
//[x10 x11] [x10*x1 + x11*x2]
result.x1 = theHessian.x00 * thePoint.x1 + theHessian.x01 * thePoint.x2;
result.x2 = theHessian.x10 * thePoint.x1 + theHessian.x11 * thePoint.x2;
}
//main function
double myFunction(point thePoint)
{ return (pow(thePoint.x1, 2)+pow(thePoint.x2, 2)); };
//------------------------------------------------------------------------------
int main(int argc, char *argv[])
{
//declaring variables
gradient myGradient;
hessian myHessian;
point myPoint;
point tempPoint;
double accuracy;
int n = 0;
//stating input variables
cout<<"Write the input value x1 = ";
cin>>myPoint.x1;
cout<<"Write the input value x2 = ";
cin>>myPoint.x2;
cout<<"Write the accuracy = ";
cin>>accuracy;
//main algorithm
myGradient.gradientAtPoint(myPoint);
while (myGradient.itsPoint.length() > accuracy)
{
n++;
myGradient.gradientAtPoint(myPoint);
myHessian.hessianAtPoint(myPoint);
myHessian.det();
myHessian.invertHessian();
//myPoint = myPoint - (myHessian * (myGradient.itsPoint)); //[1]
tempPoint = myHessian * myGradient.itsPoint; //[2]
myPoint = myPoint - tempPoint; //[3]
//Here I changed line [1] into lines [2] and [3].
//I think both [1] and [2][3] are the same but [1] doesn't compile.
}
//writing output
cout<<"\n[x1, x2] = ["<<myPoint.x1<<", "<<myPoint.x2<<"]";
cout<<"\nminimum value = "<<myFunction(myPoint);
cout<<"\naccuracy = "<<accuracy;
cout<<"\nnumber of iterations = "<<n<<"\n";
system("PAUSE");
return EXIT_SUCCESS;
}
//And by the way - I've got other question. If I'd like the user to be able to
//change the function, is there any way to allow this application to read the
//fuction (and gradient, hessian) from the file and put inside those lines which
//are responsible for calculations?
This code compiles but still there may be some logical errors. I changed 198 into 199 and 200.
Code:
//myPoint = myPoint - (myHessian * (myGradient.itsPoint)); //[1]
tempPoint = myHessian * myGradient.itsPoint; //[2]
myPoint = myPoint - tempPoint; //[3]
//Here I changed line [1] into lines [2] and [3].
//I think both [1] and [2][3] are the same but [1] doesn't compile.
myPoint is declared in line 176.
point has got operator-
point operator-(point & thePoint1, point & thePoint2) //line 51