Commit 8c81a7dc authored by PIACIBELLO Cyrille's avatar PIACIBELLO Cyrille

API with Example from G Sylvand working

parent 990dd741
......@@ -12,6 +12,7 @@ include_directories(src)
#define the test of the API
set(IBGMRES_API_TEST
tests/testApiGMres.c
tests/testApiGMresUserCase.c
)
add_library(ibgmresdr ${IBGMRES_API_SRC})
......
......@@ -126,10 +126,30 @@ void ibgmresdr_set_parameters(int max_ite, int restart, void * tolerance,
*
* @param nbRHS : number of right hand side (will be used as
* SizeBlock)
* @param RHS : Right Hand Side
* @param RHS : Right Hand Side (Colomn major mode)
* @param X0 : Starting block of vector
*/
int ibgmresdr_solve(int nbRHS,void * RHS,void *X0,
IbGMResDr_handle handle);
/**
* @brief Retreive the block of solution computed. The solution is
* stored the same way RHS have been given in the first place.
*
* @attention The User need to copy data if he wants to read it after
* deallocating ibgmresdr
*/
void * ibgmresdr_get_results(IbGMResDr_handle handle);
/**
* @brief Get the history of convergence as an array containing 3
* double for each steps : Min over the vectors, max over the vectors,
* and timestamp.
*
* @attention The User need to copy data if he wants to read it after
* deallocating ibgmresdr
*/
double * ibgmresdr_get_logs(IbGMResDr_handle handle);
#endif //IBGMRESDR_H
......@@ -22,6 +22,10 @@ public:
virtual void set_parameters(int inMaxIte,int restart, void * inTolerance){};
virtual int solve(int, void * , void *){};
virtual void * get_results(){};
virtual double * get_logs(){};
};
......@@ -38,19 +42,28 @@ class IbGMResDrEngine : public EngineDispatcher{
S tolerance;
int restart;
//Solution storage
Block<T>* sol;
Logger<S>* log;
public:
IbGMResDrEngine(void * userMatrix, int dim, void * inUserEnv) :
matrix(nullptr),
userEnvPtr(inUserEnv),
max_iter(0),
tolerance(0),
restart(0)
restart(0),
sol(nullptr),
log(nullptr)
{
matrix = new UserMatrix<T>(userEnvPtr,userMatrix,dim);
}
~IbGMResDrEngine(){
delete matrix;
if(sol)
delete sol;
if(log)
delete log;
}
/**
......@@ -71,18 +84,37 @@ public:
}
//Solve method : RHS and X0 have the same size
int solve(int nbRHS, void * RHS , void * X0){
int solve(int nbRHS, void * RHS , void * X0)override{
//Convert raw data to Block
Block<T> X_init(nbRHS,matrix->size());
Block<T> B(nbRHS,matrix->size());
X_init.InitBlock(reinterpret_cast<T*>(X0));
B.InitBlock(reinterpret_cast<T*>(RHS));
//Init an empty block to store solution
Block<T> Sol(nbRHS,matrix->size());
sol = new Block<T>(nbRHS,matrix->size());
log = new Logger<S>();
//Call to BlockGMRes with parameters
int res = BlockGMRes<UserMatrix<T>,T,S>(*matrix,B,X_init,max_iter,Sol,max_iter);
int res = BlockGMRes<UserMatrix<T>,T,S>(*matrix,B,X_init,max_iter,*sol,
max_iter,*log,tolerance);
return res;
}
void* get_results()override{
if(sol){
return sol->getPtr();
}else{
return nullptr;
}
}
double* get_logs()override{
if(log){
double * res = log->writeDownArray();
return res;
}else{
return nullptr;
}
}
};
......@@ -111,4 +143,12 @@ extern "C" void ibgmresdr_dealloc(IbGMResDr_handle handle){
delete (reinterpret_cast<EngineDispatcher*>(handle));
}
extern "C" void * ibgmresdr_get_results(IbGMResDr_handle handle){
return (reinterpret_cast<EngineDispatcher*>(handle))->get_results();
}
extern "C" double * ibgmresdr_get_logs(IbGMResDr_handle handle){
return (reinterpret_cast<EngineDispatcher*>(handle))->get_logs();
}
#endif //IBGMRESENGINE_HPP
......@@ -103,6 +103,8 @@ int main(int ac, char ** av){
int nb_iter = ibgmresdr_solve(userEnv->nbRHS,RHS,X0,handle);
double * res = ibgmresdr_get_results(handle);
ibgmresdr_dealloc(handle);
free(userEnv);
......
......@@ -10,6 +10,7 @@ project(ib-gmres-dr C CXX)
#Trouver un moyen de tester que le compilateur support la norme 2011
#(regarder éventuellement Simgrid/boost etc)
set(CMAKE_CXX_FLAGS "-std=c++11")
set(CMAKE_C_FLAGS "-std=c99")
#Ajout à la liste des repertoires dans lesquels CMAKE cherche ses modules.
list(APPEND CMAKE_MODULE_PATH "${CMAKE_SOURCE_DIR}/cmake_modules/morse/find")
......
......@@ -121,7 +121,7 @@ public:
}
/**
* @brief Init block with rax datas
* @brief Init block with raw datas
*/
void InitBlock(Scalar * rawDatas){
memcpy(getPtr(),rawDatas, sizeof(Scalar)*ldb*SizeBlock);
......
......@@ -6,7 +6,7 @@
#include "Hess.hpp"
#include "Base.hpp"
#include "Block.hpp"
#include "Logger.hpp"
/**
* @brief Computation of the GMRES Block.
......@@ -18,6 +18,7 @@ template<class Matrix,typename Scalar,typename Primary=Scalar>
int BlockGMRes(Matrix& A, Block<Scalar>& B, Block<Scalar>& X0,
//Matrix preCond Deal with that later
int max_iter, Block<Scalar>& Sol, int nb_iter_restart,
Logger<Primary>& log,
Primary epsilon = 0.00001){
// if(preCond.size() != 0){
......@@ -40,6 +41,7 @@ int BlockGMRes(Matrix& A, Block<Scalar>& B, Block<Scalar>& X0,
Block<Scalar> Yj(SizeBlock,dim);
//Add a loop for restart
int k=0,nbRestart=0, globalStep=0;
log.init_log();
while(k<max_iter){
std::cout<<"\t\t\tRestarting number "<<nbRestart<<" ("<<nb_iter_restart<<")\n";
......@@ -93,6 +95,7 @@ int BlockGMRes(Matrix& A, Block<Scalar>& B, Block<Scalar>& X0,
{//Display residual to see convergence
std::pair<Primary,Primary> MinMax
= Hess.template displayBlockResidual<Primary>(Yj,R1,normB);
log.add_ite(j,MinMax.first,MinMax.second);
std::cout<<"["<<j<<"]\t"<<"Min "<<MinMax.first
<<"\tMax "<<MinMax.second<<"\n";
if(MinMax.second <= epsilon){
......
#ifndef LOGGER_HPP
#define LOGGER_HPP
#include <time.h>
#include <sys/time.h>
#include <functional>
template<class Primary>//Primary will be double or float
class Logger{
typedef struct ite{
int curr_ite;
Primary min;
Primary max;
double timestamp;
}Ite;
std::vector<Ite> iterations;
double last_tic;
double * array;
public:
Logger():
iterations(0),
last_tic(0),array(nullptr)
{
}
~Logger(){
if(array){
delete [] array;
}
}
void init_log(){
last_tic = GetTime();
}
void add_ite(int curr,Primary min, Primary max){
double a = GetTime();
iterations.push_back(Ite{curr,min,max,a-last_tic});
last_tic = a;
}
static double GetTime(){
timeval t;
gettimeofday(&t, NULL);
return double(t.tv_sec) + (double(t.tv_usec)/1000000.0);
}
void forEachIte(std::function<void(int,Primary,Primary,double)> func){
for(auto & ite : iterations){
func(ite.curr_ite,ite.min,ite.max,ite.timestamp);
}
}
double * writeDownArray(){
array = new double[3 * iterations.size()];
for(int i =0 ; i<iterations.size() ; ++i){
array[3*i] = iterations[i].min;
array[3*i+1] = iterations[i].max;
array[3*i+2] = iterations[i].timestamp;
}
return array;
}
};
#endif //LOGGER_HPP
#include "../src/Block.hpp"
#include "../src/BlockGMRes.hpp"
#include "../src/LapackInterface.hpp"
#include "../src/Logger.hpp"
/**
* Matrix Class : need to implement MatBlockVect, MatBaseproduct and
* size member function in order to works with BlockGMRes algorithm
......@@ -105,8 +105,11 @@ int main(int ac, char ** av){
Block<Scalar> sol(SizeBlock,dim);
Logger<Scalar> log;
int nb_iter = BlockGMRes<InputMatrix<Scalar>,Scalar>(mat,RHS,
X0,nbBlock,sol,nbBlock);
X0,nbBlock,sol,nbBlock,log);
//Check solution
{
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment