diff --git a/Domain.cpp b/Domain.cpp index 270a474..dc4d484 100644 --- a/Domain.cpp +++ b/Domain.cpp @@ -1,26 +1,25 @@ // // Created by Tristan Liardon on 02.12.2018. // #include "Domain.hpp" // Constructor and destructor Domain::Domain() {} Domain::Domain(double a_, double b_) { a = a_; b = b_; + dimension = 1; } -Domain::Domain(double a_, double b_, double c_, double d_) { +Domain::Domain(vector> boundaries_) { - a = a_; - b = b_; - c = c_; - d = d_; + boundaries = boundaries_; + dimension = (int)boundaries.size(); } Domain::~Domain() {} \ No newline at end of file diff --git a/Domain.hpp b/Domain.hpp index c449d71..880defe 100644 --- a/Domain.hpp +++ b/Domain.hpp @@ -1,30 +1,33 @@ /* * Created by Tristan Liardon on 02.12.2018. * * Domain object is defined by the boundaries (two for each dimension) */ #ifndef PCSC_PROJET_DOMAIN_HPP #define PCSC_PROJET_DOMAIN_HPP +#include + +using namespace std; class Domain { public: // Constructor and destructor Domain(); Domain(double a_, double b_); // 1D - Domain(double a_, double b_, double c_, double d_); // 2D + Domain(vector>); // 2D virtual ~Domain(); // Boundaries double a; double b; - double c; - double d; + vector> boundaries; + int dimension; }; #endif //PCSC_PROJET_DOMAIN_HPP diff --git a/Function.cpp b/Function.cpp index c9c8473..f365326 100644 --- a/Function.cpp +++ b/Function.cpp @@ -1,32 +1,21 @@ // // Created by Tristan Liardon on 02.12.2018. // #include "Function.hpp" -#include -using namespace std; // Constructor and destructor Function::Function() {} -Function::Function(double (*functionPointer_)(double)) { +Function::Function(double (*functionPointer_)(vector)) { functionPointer = functionPointer_; } Function::~Function() {} -// Get result -> WHERE USER HAS TO WRITE HIS FUNCTION -double Function::GetResult(double input) { +// Get result +double Function::GetResult(vector input) { - return functionPointer(input); // R -> R + return functionPointer(input); -} - - double Function::GetResult(double *input) { - - double x = input[0]; - double y = input[1]; - - return x * y; // R2 -> R - -} +} \ No newline at end of file diff --git a/Function.hpp b/Function.hpp index 613fe04..146ff06 100644 --- a/Function.hpp +++ b/Function.hpp @@ -1,32 +1,35 @@ /* * Created by Tristan Liardon on 02.12.2018. * * A function object that needs a pointer to a user defined function. This allows the user defined function to be called * using a simple syntax FunctionObject.GetResult(input) from everywhere without risking mistakes in the complex * writing of the function pointer argument. */ #ifndef PCSC_PROJET_FUNCTION_HPP #define PCSC_PROJET_FUNCTION_HPP +#include + +using namespace std; + class Function { public: // Constructor and destructor Function(); - explicit Function(double (*functionPointer_)(double)); + explicit Function(double (*functionPointer_)(vector)); virtual ~Function(); - // Get result TODO Add template ? - double GetResult(double input); - double GetResult(double *input); + // Get result function + double GetResult(vector input); private: - double(*functionPointer)(double); // Pointer to the user-defined function + double(*functionPointer)(vector); // Pointer to the user-defined function }; #endif //PCSC_PROJET_FUNCTION_HPP diff --git a/IntegralSolver1D.cpp b/IntegralSolver1D.cpp index 67bc59e..56105e0 100644 --- a/IntegralSolver1D.cpp +++ b/IntegralSolver1D.cpp @@ -1,75 +1,87 @@ #include "IntegralSolver1D.hpp" #include "Function.hpp" //constructor, deconstructor IntegralSolver1D::IntegralSolver1D(Function function_, Domain domain, Method method) { a = domain.a; b = domain.b; - n = method.subdivision; + n = method.subdivision1D; type = method.name; function = function_; if (type == "Midpoint") { result = this->IntegrationMidpoint(); } if (type == "Trapz") { result = this->IntegrationTrapeze(); } if (type == "Simpson") { result = this->IntegrationSimpson(); } } IntegralSolver1D::~IntegralSolver1D() {} //Midpoint method: double IntegralSolver1D::IntegrationMidpoint() { double result = 0; double size = (b-a)/n; for (int i = 0; i < n; i++){ - result += size*function.GetResult(size/2+i*size); + vector inputVector; + inputVector.assign(1,size/2+i*size); + result += size*function.GetResult(inputVector); } return result; } //Trapeze method: double IntegralSolver1D::IntegrationTrapeze() { double result = 0; double size = (b-a)/n; for (int i = 0; i < n; i++){ - result += size*( function.GetResult(size*i) + function.GetResult(size*(i+1)) )/2; + vector arg1; + arg1.assign(1, size*i); + vector arg2; + arg2.assign(1, size*(i+1)); + result += size*( function.GetResult(arg1) + function.GetResult(arg2) )/2; } return result; } //Simpson Method: double IntegralSolver1D::IntegrationSimpson() { double result = 0; double size = (b-a)/n; for (int i = 0; i < n; i++){ - result += size/6 * (function.GetResult(size*i) - + 4*function.GetResult(size/2+i*size) - + function.GetResult(size*(i+1)) ); + vector arg1; + vector arg2; + vector arg3; + arg1.assign(1, size*i); + arg2.assign(1, size/2+i*size); + arg3.assign(1, size*(i+1)); + result += size/6 * (function.GetResult(arg1) + + 4*function.GetResult(arg2) + + function.GetResult(arg3) ); } return result; } // Get result method double IntegralSolver1D::GetResult() { return result; } \ No newline at end of file diff --git a/IntegralSolver1D.hpp b/IntegralSolver1D.hpp index 596d116..2a357d4 100644 --- a/IntegralSolver1D.hpp +++ b/IntegralSolver1D.hpp @@ -1,32 +1,33 @@ #ifndef PROJECT_INTEGRALSOLVER1D_HPP #define PROJECT_INTEGRALSOLVER1D_HPP #include "Function.hpp" #include "Domain.hpp" #include "Method.hpp" +#include class IntegralSolver1D { public: //constructor, deconstructor IntegralSolver1D(Function function_, Domain domain, Method method); virtual ~IntegralSolver1D(); //integration methods double IntegrationMidpoint(); double IntegrationTrapeze(); double IntegrationSimpson(); // Get result method double GetResult(); private: double a; double b; double n; std::string type; double result; Function function; }; #endif //PROJECT_INTEGRALSOLVER1D_HPP diff --git a/IntegralSolver2D.cpp b/IntegralSolver2D.cpp index df5542b..5104f57 100644 --- a/IntegralSolver2D.cpp +++ b/IntegralSolver2D.cpp @@ -1,98 +1,145 @@ #include "IntegralSolver2D.hpp" #include -//fonction ici jpp -double f(double x, double y) { return x+y; } - //constructor, deconstructor -IntegralSolver2D::IntegralSolver2D(double a_, double b_, double c_, double d_, int n_, int m_){ //interval [a,b],[c,d] - a = a_; - b = b_; - c = c_; - d = d_; - n = n_; // number of intervals in [a,b] - m = m_; // number of intervals in [c,d] +IntegralSolver2D::IntegralSolver2D(Function function_, Domain domain_, Method method_){ //interval [a,b],[c,d] + + type = method_.name; + a = domain_.boundaries[0][0]; + b = domain_.boundaries[0][1]; + c = domain_.boundaries[1][0]; + d = domain_.boundaries[1][1]; + n = method_.subdivision[0]; // number of intervals in [a,b] + m = method_.subdivision[1]; // number of intervals in [c,d] + function = function_; + if (type == "Midpoint") { + + result = this->IntegrationMidpoint(); + + } + + if (type == "Trapz") { + + result = this->IntegrationTrapeze(); + + } + + if (type == "Simpson") { + + result = this->IntegrationSimpson(); + + } + } IntegralSolver2D::~IntegralSolver2D() {} +// Get result method +double IntegralSolver2D::GetResult() { + + return result; + +} + //Midpoint method: double IntegralSolver2D::IntegrationMidpoint() { double result = 0; double sizeX = (b-a)/n; double sizeY = (d-c)/m; //midpoint value multiplied by the area for the whole grid for (int i = 0; i < n; i++){ for (int j = 0; j < m; j++){ - result += sizeX*sizeY* f( i*sizeX + sizeX/2, j*sizeY + sizeY/2); + + vector arg; + arg.push_back(i*sizeX + sizeX/2); + arg.push_back(j*sizeY + sizeY/2); + result += sizeX * sizeY * function.GetResult(arg); } } return result; } //Trapeze method: double IntegralSolver2D::IntegrationTrapeze() { double result = 0; double sizeX = (b-a)/n; double sizeY = (d-c)/m; // average of the 4 values of the corner multiplied by the area for the whole grid for (int i = 0; i < n; i++){ for (int j = 0; j < m; j++) { - result += sizeX*sizeY* ( f(i*sizeX,j*sizeY) + f((i+1)*sizeX,j*sizeY) + f(i*sizeX,(j+1)*sizeY) + f((i+1)*sizeX,(j+1)*sizeY) )/4; + + vector corner1; + corner1.push_back(i*sizeX); + corner1.push_back(j*sizeY); + + vector corner2; + corner2.push_back((i+1)*sizeX); + corner2.push_back(j*sizeY); + + vector corner3; + corner3.push_back(i*sizeX); + corner3.push_back((j+1)*sizeY); + + vector corner4; + corner4.push_back((i+1)*sizeX); + corner4.push_back((j+1)*sizeY); + + result += sizeX * sizeY / 4 * (function.GetResult(corner1) + function.GetResult(corner2) + + function.GetResult(corner3) + function.GetResult(corner4)); } } return result; } //Simpson Method: double IntegralSolver2D::IntegrationSimpson() { double result = 0; - // algorithm source http://mathfaculty.fullerton.edu/mathews/n2003/SimpsonsRule2DMod.html - double h = (b-a)/(2*n); - double k = (d-c)/(2*m); - - result += f(a,c)+f(a,d)+f(b,c)+f(b,d); - - for (int j = 1; j <= n; j++){ result += 4*f(a , c+(2*j-1)*k ); } - for (int j = 1; j <= n-1; j++){ result += 2*f(a , c+(2*j)*k ); } - for (int j = 1; j <= n; j++){ result += 4*f(b , c+(2*j-1)*k ); } - for (int j = 1; j <= n-1; j++){ result += 2*f(b , c+(2*j)*k ); } - - for (int i = 1; i <= m; i++){ result += 4*f(a+(2*i-1)*k , c); } - for (int i = 1; i <= m; i++){ result += 4*f(a+(2*i-1)*k , d); } - for (int i = 1; i <= m-1; i++){ result += 2*f(a+(2*i)*k , c); } - for (int i = 1; i <= m-1; i++){ result += 2*f(a+(2*i)*k , d); } - - for (int i = 1; i <= n; i++){ - for (int j = 1; j <= m; j++){ - result += 16*f( a+(2*i-1)*h , c+(2*j-1)*k ); - } - } - for (int i = 1; i <= n-1; i++){ - for (int j = 1; j <= m; j++){ - result += 8*f( a+(2*i-1)*h , c+(2*j)*k ); - } - } - for (int i = 1; i <= n; i++){ - for (int j = 1; j <= m-1; j++){ - result += 8*f( a+(2*i)*h , c+(2*j-1)*k ); - } - } - for (int i = 1; i <= n-1; i++){ - for (int j = 1; j <= m-1; j++){ - result += 4*f( a+(2*i)*h , c+(2*j)*k ); - } - } - - result = result*h*k/9; +// // algorithm source http://mathfaculty.fullerton.edu/mathews/n2003/SimpsonsRule2DMod.html +// double h = (b-a)/(2*n); +// double k = (d-c)/(2*m); +// +// result += f(a,c)+f(a,d)+f(b,c)+f(b,d); +// +// for (int j = 1; j <= n; j++){ result += 4*f(a , c+(2*j-1)*k ); } +// for (int j = 1; j <= n-1; j++){ result += 2*f(a , c+(2*j)*k ); } +// for (int j = 1; j <= n; j++){ result += 4*f(b , c+(2*j-1)*k ); } +// for (int j = 1; j <= n-1; j++){ result += 2*f(b , c+(2*j)*k ); } +// +// for (int i = 1; i <= m; i++){ result += 4*f(a+(2*i-1)*k , c); } +// for (int i = 1; i <= m; i++){ result += 4*f(a+(2*i-1)*k , d); } +// for (int i = 1; i <= m-1; i++){ result += 2*f(a+(2*i)*k , c); } +// for (int i = 1; i <= m-1; i++){ result += 2*f(a+(2*i)*k , d); } +// +// for (int i = 1; i <= n; i++){ +// for (int j = 1; j <= m; j++){ +// result += 16*f( a+(2*i-1)*h , c+(2*j-1)*k ); +// } +// } +// for (int i = 1; i <= n-1; i++){ +// for (int j = 1; j <= m; j++){ +// result += 8*f( a+(2*i-1)*h , c+(2*j)*k ); +// } +// } +// for (int i = 1; i <= n; i++){ +// for (int j = 1; j <= m-1; j++){ +// result += 8*f( a+(2*i)*h , c+(2*j-1)*k ); +// } +// } +// for (int i = 1; i <= n-1; i++){ +// for (int j = 1; j <= m-1; j++){ +// result += 4*f( a+(2*i)*h , c+(2*j)*k ); +// } +// } +// +// result = result*h*k/9; return result; } diff --git a/IntegralSolver2D.hpp b/IntegralSolver2D.hpp index 2be8668..25ac99a 100644 --- a/IntegralSolver2D.hpp +++ b/IntegralSolver2D.hpp @@ -1,25 +1,37 @@ #ifndef PROJECT_INTEGRALSOLVER2D_HPP #define PROJECT_INTEGRALSOLVER2D_HPP +#include "Function.hpp" +#include "Method.hpp" +#include "Domain.hpp" +#include + class IntegralSolver2D { public: //constructor, deconstructor - IntegralSolver2D(double a_, double b_, double c_, double d_, int n_, int m_); + IntegralSolver2D(Function function_, Domain domain_, Method method_); virtual ~IntegralSolver2D(); + // Get result method + double GetResult(); + //integration methods double IntegrationMidpoint(); double IntegrationTrapeze(); double IntegrationSimpson(); private: double a; double b; double c; double d; double n; double m; + Function function; + Domain domain; + std::string type; + double result; }; #endif //PROJECT_INTEGRALSOLVER2D_HPP diff --git a/IntegrationSolver.cpp b/IntegrationSolver.cpp index 77b92f8..2e6b9f6 100644 --- a/IntegrationSolver.cpp +++ b/IntegrationSolver.cpp @@ -1,35 +1,41 @@ // // Created by Tristan Liardon on 02.12.2018. // #include "IntegrationSolver.hpp" // Constructor and destructor IntegrationSolver::IntegrationSolver(Function function_, Domain domain_, Method method_) { function = function_; domain = domain_; method = method_; } IntegrationSolver::~IntegrationSolver() {} // Get the result double IntegrationSolver::GetResult() { return result; } // Order of computing the integral double IntegrationSolver::Compute() { - if (1==1) { + if (domain.dimension == 1) { IntegralSolver1D Solver(function, domain, method); result = Solver.GetResult(); } + else if (domain.dimension == 2) { + + IntegralSolver2D Solver(function, domain, method); + result = Solver.GetResult(); + } + } diff --git a/IntegrationSolver.hpp b/IntegrationSolver.hpp index c1fdbd5..c11ac07 100644 --- a/IntegrationSolver.hpp +++ b/IntegrationSolver.hpp @@ -1,40 +1,41 @@ /* * Created by Tristan Liardon on 02.12.2018. * * IntegrationSolver object needs a function object to integrate on a domain (given as an object) using a method (given * as an object). */ #ifndef PCSC_PROJET_INTEGRATIONSOLVER_HPP #define PCSC_PROJET_INTEGRATIONSOLVER_HPP #include "Function.hpp" #include "Domain.hpp" #include "Method.hpp" #include "IntegralSolver1D.hpp" +#include "IntegralSolver2D.hpp" class IntegrationSolver { public: // Constructor and destructor IntegrationSolver(Function function_, Domain domain_, Method method_); virtual ~IntegrationSolver(); // Get the result double GetResult(); // Order the computing double Compute(); private: Function function; Domain domain; Method method; double result; }; #endif //PCSC_PROJET_INTEGRATIONSOLVER_HPP diff --git a/Method.cpp b/Method.cpp index 9bc77f6..e0a4995 100644 --- a/Method.cpp +++ b/Method.cpp @@ -1,16 +1,23 @@ // // Created by Tristan Liardon on 02.12.2018. // #include "Method.hpp" // Constructor and destructor Method::Method() {} Method::Method(std::string name_, int subdivision_) { + name = name_; + subdivision1D = subdivision_; + +} + +Method::Method(std::string name_, vector subdivision_) { + name = name_; subdivision = subdivision_; } Method::~Method() {} \ No newline at end of file diff --git a/Method.hpp b/Method.hpp index 77174bb..b464353 100644 --- a/Method.hpp +++ b/Method.hpp @@ -1,28 +1,33 @@ /* * Created by Tristan Liardon on 02.12.2018. * * Method object is defined by the method's name and the specific parameters needed by the method (mainly the number of * subdivision for each dimension) */ #ifndef PCSC_PROJET_METHOD_HPP #define PCSC_PROJET_METHOD_HPP #include +#include + +using namespace std; class Method { public: // Constructor and destructor Method(); Method(std::string name_, int subdivision_); + Method(std::string name_, vector subdivision_); virtual ~Method(); std::string name; - int subdivision; + int subdivision1D; + vector subdivision; }; #endif //PCSC_PROJET_METHOD_HPP diff --git a/main.cpp b/main.cpp index 95efcc8..d78276b 100644 --- a/main.cpp +++ b/main.cpp @@ -1,45 +1,82 @@ #include #include #include "IntegralSolver1D.hpp" #include "IntegralSolver2D.hpp" #include "IntegrationSolver.hpp" #include "Function.hpp" #include "Domain.hpp" #include "Method.hpp" // User defined function -double myFunction(double input) { +double myFunction1D(vector x) { - return input * input; + return x[0] * x[0]; } +double myFunction2D(vector x) { + + return x[0]*x[1]; +} + int main(int argc, char *argv[]) { + ////////// 2D test ////////// + // Creation of the pointer to the function to integrate + double (*functionPointer2D)(vector); + functionPointer2D = &myFunction2D; + + // Object creation + Function function2D(functionPointer2D); + vector> boundaries2D; + vector dim1; + vector dim2; + dim1.push_back(0); + dim1.push_back(5); + dim2.push_back(0); + dim2.push_back(5); + boundaries2D.push_back(dim1); + boundaries2D.push_back(dim2); + Domain domain2D(boundaries2D); + vector subdivision; + subdivision.push_back(100); + subdivision.push_back(100); + Method method2D("Trapz", subdivision); + IntegrationSolver Integration4(function2D, domain2D, method2D); + + // Asking to compute the result + double result4 = Integration4.Compute(); + std::cout << "Result is : " << result4 << " using " << method2D.name << "'s method" << std::endl; + + ////////// 1D test ////////// // Creation of the pointer to the function to integrate - double (*functionPointer)(double); - functionPointer = &myFunction; - - Function function(functionPointer); - Domain domain(0, 5); - Method method("Simpson", 100); - IntegrationSolver Object = IntegrationSolver(function, domain, method); - double result = Object.Compute(); - std::cout << "Result is : " << result << " using " << method.name << "'s method" << std::endl; - - // interval of the domain [a,b],[c,d] - double a = 0; - double b = 5; - double c = 0; - double d = 10; - int n = 100; // number of intervals in [a,b] - int m = 100; // number of intervals in [c,d] - - IntegralSolver2D function2D(a,b,c,d,n,m); // TODO la fonction est définie dans integralsolver2d.cpp... à changer! - std::cout << "Midpoint integration result: " << function2D.IntegrationMidpoint() << "\n"; - std::cout << "Trapeze integration result: " << function2D.IntegrationTrapeze() << "\n"; - std::cout << "Simpson integration result: " << function2D.IntegrationSimpson() << "\n"; + double (*functionPointer1D)(vector); + functionPointer1D = &myFunction1D; + + // Object creation + Function function1D(functionPointer1D); + Domain domain1D(0, 5); + Method method1("Simpson", 100); + IntegrationSolver Integration1 = IntegrationSolver(function1D, domain1D, method1); + + // Asking to compute the result + double result1 = Integration1.Compute(); + std::cout << "Result is : " << result1 << " using " << method1.name << "'s method" << std::endl; + + Method method2("Midpoint", 100); + IntegrationSolver Integration2 = IntegrationSolver(function1D, domain1D, method2); + + // Asking to compute the result + double result2 = Integration2.Compute(); + std::cout << "Result is : " << result2 << " using " << method2.name << "'s method" << std::endl; + + Method method3("Trapz", 100); + IntegrationSolver Integration3 = IntegrationSolver(function1D, domain1D, method3); + + // Asking to compute the result + double result3 = Integration3.Compute(); + std::cout << "Result is : " << result3 << " using " << method3.name << "'s method" << std::endl; return 0; }