Commit a79903b3 authored by antoine lucas's avatar antoine lucas
Browse files

refactor bigvec and bigvec_q

parent e21517b8
......@@ -42,15 +42,10 @@ bigmod & bigmod::operator= (const bigmod& rhs)
return(*this);
}
bigmod & bigmod::inv ()
bigmod bigmod::inv () const
{
if (inverse != NULL){
inverse = NULL;
delete inverse;
}
if(value.isNA() || modulus.isNA()) {
inverse = new DefaultBigMod();
return *inverse;
return bigmod();
}
mpz_t val;
......@@ -61,11 +56,10 @@ bigmod & bigmod::inv ()
if(wOpt != R_NilValue && Rf_asInteger(wOpt))
warning(_("inv(x) returning NA as x has no inverse"));
inverse = new DefaultBigMod();
return *inverse; // return NA; was
return bigmod(); // return NA; was
}
inverse = new DefaultBigMod(val, modulus );
return *inverse;
return bigmod(val, modulus);
}
......
......@@ -36,14 +36,10 @@ extern "C" {
*/
class bigmod {
private:
bigmod * inverse;
/** \brief copy operator */
bigmod(const bigmod & rhs) :
value(((bigmod)rhs).getValue()),modulus(((bigmod)rhs).getModulus()) {
}
/** optional source */
biginteger * value_ptr;
biginteger * modulus_ptr;
protected:
/** \brief Value of our bigmod -- only references*/
biginteger & value;
......@@ -51,13 +47,55 @@ class bigmod {
biginteger & modulus;
public:
/** keep both references value / modulus
*/
bigmod(biginteger& value_,
biginteger& modulus_) :
inverse(NULL),
value_ptr(NULL),
modulus_ptr(NULL),
value(value_),modulus(modulus_) {};
/** keep references value / modulus is new object.
*/
bigmod(biginteger& value_) :
value_ptr(NULL),
modulus_ptr(new biginteger()),
value(value_),modulus(*modulus_ptr) {};
/**
* create 2 new objects valus / modulus.
*/
bigmod(const biginteger& value_,
const biginteger& modulus_) :
value_ptr(new biginteger(value_)),
modulus_ptr(new biginteger(modulus_)),
value(*value_ptr),modulus(*modulus_ptr) {};
bigmod(const biginteger& value_) :
value_ptr(new biginteger(value_)),
modulus_ptr(new biginteger()),
value(*value_ptr),modulus(*modulus_ptr) {};
bigmod() :
value_ptr(new biginteger()),
modulus_ptr(new biginteger()),
value(*value_ptr),modulus(*modulus_ptr) {};
/** \brief copy operator */
bigmod(const bigmod & rhs) :
value_ptr(new biginteger(rhs.getValue())),
modulus_ptr(new biginteger(rhs.getModulus())),
value(*value_ptr),modulus(*modulus_ptr) {
};
virtual ~bigmod(){
if(inverse != NULL) delete inverse;
if(value_ptr != NULL) delete value_ptr;
if(modulus_ptr != NULL) delete modulus_ptr;
};
/**
......@@ -75,7 +113,7 @@ class bigmod {
return(mpz_sgn(getValue().getValueTemp()));
}
bigmod & inv () ;
bigmod inv () const;
biginteger & getValue() {
......
......@@ -300,53 +300,6 @@ bool operator!=(const bigvec & rhs, const bigvec& lhs)
//
// \brief substract lambda[0] * line j to line i
//
void bigvec::subLine(unsigned int i,unsigned int j,const bigvec & lambda)
{
if(nrow <= 0)
error(_("Need matrix with at least one row to do this operation"));
unsigned int k, n = value.size() / (unsigned int) nrow;
if(modulus.size() != 1)
{
for(k=0; k < n; ++k)
value[i + k*nrow] = value[i + k*nrow] - ( value[j + k*nrow] * lambda.value[0] ) ;
}
else
for(k=0; k < n ; ++k)
{
value[i + k*nrow] = value[i + k*nrow] - ( value[j + k*nrow] * lambda.value[0] ) ;
value[i + k*nrow] = value[i + k*nrow] % modulus[0];
}
}
/*
* \brief multiply line i by lambda
*/
void bigvec::mulLine(unsigned int i,const bigvec & lambda)
{
if(nrow <= 0)
error(_("Need matrix with at least one row to do this operation"));
unsigned int k;
// n number of columns
unsigned int n = value.size() / (unsigned int) nrow;
if(modulus.size() != 1)
{
for(k=0; k < n; ++k)
value[i + k*nrow] = value[i + k*nrow] * lambda.value[0] ;
}
else
for(k=0; k < n ; ++k)
{
value[i + k*nrow] = value[i + k*nrow] * lambda.value [0] ;
value[i + k*nrow] = value[i + k*nrow] % modulus[0];
}
}
// never used
void bigvec::print()
{
......
......@@ -119,17 +119,6 @@ class bigvec : public matrix::Matrix<bigmod> {
*/
bigvec & operator= (const bigvec & rhs);
/**
* \brief substract lambda * line j to line i
*/
void subLine(unsigned int i,unsigned int j,const bigvec & lambda);
/**
* \brief multiply line i by lambda
*/
void mulLine(unsigned int i,const bigvec & lambda);
/** \brief print matrix to stdout
*
* use for debug purpose
......
......@@ -41,11 +41,27 @@ bigvec_q & bigvec_q::operator= (const bigvec_q & rhs)
}
bigrational bigvec_q::operator[] (unsigned int i) const
const bigrational & bigvec_q::operator[] (unsigned int i) const
{
return(value[i]);
}
bigrational & bigvec_q::operator[] (unsigned int i)
{
return(value[i]);
}
bigrational & bigvec_q::get(unsigned int row, unsigned int col) {
return (*this)[row + col*nrow];
}
void bigvec_q::set(unsigned int row, unsigned int col, const bigrational & val) {
set( row + col*nrow,val);
}
void bigvec_q::set(unsigned int i,const bigrational & val)
{
//DEBUG !!
......@@ -79,6 +95,10 @@ unsigned int bigvec_q::size() const
return(value.size());
}
unsigned int bigvec_q::nRows() const {
return abs(nrow);
}
void bigvec_q::resize(unsigned int n)
{
value.resize(n);
......@@ -92,32 +112,6 @@ void bigvec_q::clear()
//
// \brief substract lambda[0] * line j to line i
//
void bigvec_q::subLine(unsigned int i,unsigned int j,bigvec_q lambda)
{
if(nrow <= 0)
error(_("Need matrix with at least one row to do this operation"));
unsigned int k, n = value.size() / nrow;
for(k=0; k < n; ++k)
value[i + k*nrow] = value[i + k*nrow] - ( value[j + k*nrow] * lambda.value[0] ) ;
}
/*
* \brief multiply line i by lambda[0]
*/
void bigvec_q::mulLine(unsigned int i, bigvec_q lambda)
{
if(nrow <= 0)
error(_("Need matrix with at least one row to do this operation"));
unsigned int k, n = value.size() / nrow;
for(k=0; k < n; ++k)
value[i + k*nrow] = value[i + k*nrow] * lambda.value[0];
}
// never used
void bigvec_q::print()
......
......@@ -22,7 +22,7 @@
* It is a class composed of a vector of bigrational
* and an nrow parameter.
*/
class bigvec_q {
class bigvec_q : public matrix::Matrix<bigrational> {
public:
/** \brief The real value */
......@@ -62,9 +62,18 @@ class bigvec_q {
/**
* \brief extract a value.
* \note it returns a copy of the value.
* \note it returns a reference on the value.
*/
bigrational operator[] (unsigned int i) const;
const bigrational & operator[] (unsigned int i) const;
/**
* \brief extract a value.
* \note it returns a reference on the value.
*/
bigrational & operator[] (unsigned int i) ;
bigrational & get(unsigned int row, unsigned int col) ;
void set(unsigned int row, unsigned int col, const bigrational & val);
/**
* \brief Set a value
......@@ -92,22 +101,14 @@ class bigvec_q {
* \brief resize value
*/
void resize(unsigned int n);
unsigned int nRows() const;
/**
* \brief clear all.
*/
void clear();
/**
* \brief substract lambda[0] * line j to line i
*/
void subLine(unsigned int i,unsigned int j,bigvec_q lambda);
/**
* \brief multiply line i by lambda[0]
*/
void mulLine(unsigned int i, bigvec_q lambda);
/** \brief print matrix to stdout
*
* use for debug purpose
......
......@@ -52,31 +52,26 @@ namespace solve_gmp_R
* A is of dimension nxn X nxm and B nxm (X will be return a B address)
* We use the Gauss algorithm
*/
template< class T> void solve (T & A , T & B)
template< class T> void solve (matrix::Matrix<T> & A , matrix::Matrix<T> & B)
{
// n = A.nrow
int i;
// a temp value
T tmpValeur(1);
// A [ i ,j] = A[ i + j * A.nrow]
for(int k = 0 ; k < A.nrow; ++k)
for(int k = 0 ; k < A.nRows(); ++k)
{
if(A[k + k* A.nrow ].sgn() == 0 )
if(A.get(k, k).sgn() == 0 )
Rf_error("System is singular");
// l_k <- (1/akk) l_k
tmpValeur.set(0, A[k + k*A.nrow].inv() );
T tmpValeur =A.get(k , k).inv() ;
A.mulLine(k,tmpValeur);
B.mulLine(k,tmpValeur);
for(i = 0; i< A.nrow; ++i)
for(int i = 0; i< A.nRows(); ++i)
{
if(i == k)
continue;
// l_i <- l_i - a_ik l_k
tmpValeur.set(0, A[i + k*A.nrow]) ;
tmpValeur= A.get(i, k) ;
A.subLine(i,k, tmpValeur) ;
B.subLine(i,k, tmpValeur) ;
}
......
......@@ -76,7 +76,7 @@ namespace matrix{
*/
void mulLine(unsigned int i,const T & lambda){
for(int k=0; k < nCols(); ++k)
set(i , k) = get(i,k) * lambda ;
set(i , k, get(i,k) * lambda );
}
Matrix & transpose();
......@@ -93,15 +93,15 @@ namespace matrix{
: source(source_p) {
};
unsigned int size(){
unsigned int size() const{
return source.size();
}
unsigned int nRows(){
unsigned int nRows() const {
return source.nCols();
}
T & get(unsigned int i, unsigned int j) const{
T & get(unsigned int i, unsigned int j) {
return source.get(j,i);
}
......
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