Commit e04ee33d authored by PIACIBELLO Cyrille's avatar PIACIBELLO Cyrille

First release ready

parent 39e47a4b
......@@ -127,7 +127,7 @@ int main(int ac , char ** av){
//Fill with the diag of Mat
for(int j=0 ; j<mat->dim1 ; ++j){
my_rpc->datas[j] = mat->datas[j*(mat->dim1)+j];
my_rpc->datas[j] = 1/(mat->datas[j*(mat->dim1)+j]);
}
//Give to the library the PreConditionner and the method to call
ibgmresdr_set_rightprecond(RightPreCondxBlockVect,my_rpc,handle);
......
......@@ -3,10 +3,11 @@
#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
* Matrix Class : need to implement MatBlockVect(), MatBaseproduct(),
* useRightPreCond(), preCondBlockVect(),preCondBaseProduct() and
* size() member function in order to works with BlockGMRes algorithm
*
* An abstract class coulbe be used to enforce that.
* An abstract class could be used to enforce that.
*/
template<class Scalar>
class InputMatrix{
......@@ -26,11 +27,21 @@ public:
}
void initMatrix(std::random_device& rd){
std::mt19937 gen(rd());
std::uniform_real_distribution<> dis(0, 2);
// std::mt19937 gen(rd());
// std::uniform_real_distribution<> dis(0, 2);
// for(int i=0 ; i<dim ; ++i){
// for(int j=0 ; j<dim ; ++j){
// data[dim*i+j] = dis(gen);
// }
// }
Scalar diag = dim;
Scalar no_diag = 1;
for(int i=0 ; i<dim ; ++i){
for(int j=0 ; j<dim ; ++j){
data[dim*i+j] = dis(gen);
if(i==j)
data[i*dim + j] = diag+j;
else
data[i*dim + j] = no_diag;
}
}
}
......@@ -50,6 +61,24 @@ public:
out.getPtr(),out.getLeadingDim(),
Scalar(1),Scalar(0));
}
bool useRightPreCond(){
return true;
}
void preCondBlockVect(Block<Scalar>& input, Block<Scalar>& output){
for(int j=0 ; j<output.getSizeBlock() ; ++j){
for(int i=0 ; i<dim ; ++i){
output.getPtr(j)[i] = input.getPtr(j)[i]*(1/data[i*dim+i]);
}
}
}
void preCondBaseProduct(Scalar * ptr,Block<Scalar>& output){
for(int j=0 ; j<output.getSizeBlock() ; ++j){
for(int i=0 ; i<dim ; ++i){
output.getPtr(j)[i] = ptr[j*dim+i]*1/(data[i*dim+i]);
}
}
}
};
......@@ -58,7 +87,7 @@ int main(int ac, char ** av){
for(auto ag : args)
std::cout << ag << "\n";
int SizeBlock=5,nbBlock=10;
int SizeBlock=1,nbBlock=4;
if (args.size()>1){
SizeBlock = std::stoi(args[1]);
......@@ -108,9 +137,8 @@ int main(int ac, char ** av){
Logger<Scalar> log;
int nb_iter = BlockGMRes<InputMatrix<Scalar>,Scalar>(mat,RHS,
X0,nbBlock,sol,nbBlock,log);
X0,nbBlock,sol,nbBlock,log,
0.0000000001);
//Check solution
{
//Compute what i need for my metrics
......
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