Skip to content

Commit

Permalink
Merge branch 'master' of github.com:linbox-team/linbox
Browse files Browse the repository at this point in the history
  • Loading branch information
@jgdumas
jgdumas committed Jul 10, 2019
2 parents 520c8c6 + e1cfc07 commit f9ef2c7
Show file tree
Hide file tree
Showing 10 changed files with 174 additions and 198 deletions.
5 changes: 0 additions & 5 deletions examples/mats.C
View file
Original file line number Diff line number Diff line change
Expand Up @@ -62,11 +62,6 @@ using namespace std;
#include <linbox/util/timer.h>
#include <linbox/matrix/dense-matrix.h>

// place A: Edit here and at place B for ring change
//#include <linbox/ring/pir-modular-int32.h>
//#include <linbox/ring/pir-ntl-zz_p.h>
#include <linbox/ring/pir-ntl-zz_p.h>

using namespace LinBox;

template <class PIR>
Expand Down
64 changes: 19 additions & 45 deletions linbox/algorithms/smith-form-adaptive.inl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -35,6 +35,7 @@
#include "linbox/util/debug.h"
#include "linbox/ring/pir-modular-int32.h"
#include "linbox/ring/local2_32.h"
#include "linbox/ring/local-pir-modular.h"
#include "linbox/algorithms/smith-form-iliopoulos.h"
#include "linbox/algorithms/smith-form-local.h"
#include "linbox/algorithms/rational-solver-adaptive.h"
Expand All @@ -49,12 +50,6 @@
#include "linbox/algorithms/smith-form-adaptive.inl"
#include "linbox/solutions/valence.h"



#ifdef __LINBOX_HAVE_NTL
#include "linbox/ring/pir-ntl-zz_p.h"
#endif

namespace LinBox
{

Expand Down Expand Up @@ -150,7 +145,6 @@ namespace LinBox

}

#ifdef __LINBOX_HAVE_NTL
/* Compute the local smith form at prime p, when modular (p^e) doesnot fit in long
*/
template <class Matrix>
Expand All @@ -162,36 +156,23 @@ namespace LinBox
int order = (int)(A. rowdim() < A. coldim() ? A. rowdim() : A. coldim());
linbox_check ((s. size() >= (size_t) order) && (p > 0) && ( e >= 0));
integer T; T = order; T <<= 20; T = pow (T, (int) sqrt((double)order));
NTL::ZZ m; NTL::conv(m, 1); int i = 0; for (i = 0; i < e; ++ i) m *= p;
//if (m < T)
if (1) {
report << " Compute local Smith at " << p << '^' << e << " over PIR-ntl-ZZ_p\n";
PIR_ntl_ZZ_p R(m);
BlasMatrix <PIR_ntl_ZZ_p> A_local(R, A.rowdim(), A.coldim());
SmithFormLocal <PIR_ntl_ZZ_p> SF;
std::list <PIR_ntl_ZZ_p::Element> l;
MatrixHom::map (A_local, A);
SF (l, A_local, R);
std::list <PIR_ntl_ZZ_p::Element>::iterator l_p;
BlasVector<Givaro::ZRing<Integer> >::iterator s_p;
for (s_p = s. begin(), l_p = l. begin(); s_p != s. begin() +(ptrdiff_t) order; ++ s_p, ++ l_p)
R. convert(*s_p, *l_p);
report << " Done \n";
}
else {
std::cerr << "Compute local smith at " << p << '^' << e << std::endl;
std::cerr << "Not implemented yet.\n";
}
Givaro::Integer m(Integer::one);
for (int i = 0; i < e; ++ i) m *= p;
report << " Compute local Smith at " << p << '^' << e << " over ZZ_p\n";
typedef LocalPIRModular<Givaro::Integer> GZZ_p;
GZZ_p R(m);
BlasMatrix <GZZ_p> A_local(R, A.rowdim(), A.coldim());
SmithFormLocal <GZZ_p> SF;
std::list <GZZ_p::Element> l;
MatrixHom::map (A_local, A);
SF (l, A_local, R);
std::list <GZZ_p::Element>::iterator l_p;
BlasVector<Givaro::ZRing<Integer> >::iterator s_p;
for (s_p = s. begin(), l_p = l. begin(); s_p != s. begin() +(ptrdiff_t) order; ++ s_p, ++ l_p)
R. convert(*s_p, *l_p);
report << " Done \n";
return;
}
#else
template <class Matrix>
void SmithFormAdaptive::compute_local_big (BlasVector<Givaro::ZRing<Integer> >& s, const Matrix& A, int64_t p, int64_t e)
{
throw(LinBoxError("you need NTL to use SmithFormAdaptive",LB_FILE_LOC));
}

#endif

/* Compute the local smith form at prime p
*/
Expand Down Expand Up @@ -260,7 +241,6 @@ namespace LinBox
}


#ifdef __LINBOX_HAVE_NTL
/* Compute the k-rough part of the invariant factor, where k = 100.
* By EGV+ algorithm or Iliopoulos' algorithm for Smith form.
*/
Expand Down Expand Up @@ -309,8 +289,9 @@ namespace LinBox
}
else {
report << " Elimination start:\n";
PIR_ntl_ZZ_p R (m);
BlasMatrix<PIR_ntl_ZZ_p> A_ilio(R, A.rowdim(), A.coldim());
typedef LocalPIRModular<Givaro::Integer> GZZ_p;
GZZ_p R (m);
BlasMatrix<GZZ_p> A_ilio(R, A.rowdim(), A.coldim());
MatrixHom::map (A_ilio, A);
SmithFormIliopoulos::smithFormIn (A_ilio);
int i; BlasVector<Givaro::ZRing<Integer> >::iterator s_p;
Expand All @@ -320,13 +301,6 @@ namespace LinBox
}
report << "Compuation of the k-rough part of the invariant factors finishes.\n";
}
#else
template <class Matrix>
void SmithFormAdaptive::smithFormRough (BlasVector<Givaro::ZRing<Integer> >& s, const Matrix& A, integer m)
{
throw(LinBoxError("you need NTL to use SmithFormAdaptive",LB_FILE_LOC));
}
#endif

/* Compute the Smith form via valence algorithms
* Compute the local Smtih form at each possible prime
Expand Down
2 changes: 1 addition & 1 deletion linbox/algorithms/smith-form-local.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -123,7 +123,7 @@ namespace LinBox
} else {
typename Matrix::Iterator p_it;
for (p_it = A.Begin(); p_it != A.End(); ++p_it) {
R.divin(*p_it, g);
R.divin(*p_it, g);
}
typename LocalPID::Element x; R.neg(x, g);
R.gcdin(g,x);
Expand Down
127 changes: 85 additions & 42 deletions linbox/ring/local-pir-modular.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -44,29 +44,28 @@ namespace LinBox
typedef RingCategories::ModularTag categoryTag;
};

/// \ingroup ring
/// \ingroup ring
template <typename intType>
class LocalPIRModular : public Givaro::Modular<intType> {

public:

using Parent_t = Givaro::Modular<intType>;

using typename Givaro::Modular<intType>::Element;
//typedef typename Givaro::Modular<intType>::Element Element;
//typedef typename Givaro::Modular<intType>::Element Element;

using typename Givaro::Modular<intType>::RandIter;

//No default cstor
//No default cstor

LocalPIRModular (intType value, uint32_t exp = 1) :
Givaro::Modular<intType>(Givaro::power(value, exp)),
_irred(value),
Givaro::Modular<intType>(Givaro::power(value, exp)),
_exponent(exp)
{}
~LocalPIRModular() noexcept {};
{}

~LocalPIRModular() noexcept {};

using Parent_t:: zero;
using Parent_t:: one;
using Parent_t:: mOne;
Expand All @@ -76,45 +75,43 @@ namespace LinBox

using Parent_t:: characteristic;
using Parent_t:: cardinality;

using Parent_t:: maxCardinality;
using Parent_t:: minCardinality;

// ----- Checkers
// ----- Checkers
using Parent_t:: isZero;
using Parent_t:: isOne ;
using Parent_t:: isMOne;
using Parent_t:: areEqual;
using Parent_t:: length;
// ----- Ring-wise operators

// ----- Ring-wise operators
using Parent_t:: operator==;
using Parent_t:: operator!=;
using Parent_t:: operator=;

// ----- Initialisation
// ----- Initialisation
using Parent_t:: init ;

using Parent_t:: assign ;

using Parent_t:: convert;

using Parent_t:: reduce ;

using Parent_t:: mul;
using Parent_t:: div;
using Parent_t:: add;
using Parent_t:: sub;
using Parent_t:: neg;
using Parent_t:: inv;

using Parent_t:: mulin;
using Parent_t:: divin;
using Parent_t:: addin;
using Parent_t:: subin;
using Parent_t:: negin;
using Parent_t:: invin;

using Parent_t:: axpy ;
using Parent_t:: axpyin;

Expand All @@ -124,38 +121,84 @@ namespace LinBox
using Parent_t:: maxpy;
using Parent_t:: maxpyin;

// ----- Random generators
// ----- Random generators
using typename Parent_t:: NonZeroRandIter;
using Parent_t:: random;
using Parent_t:: nonzerorandom;

// --- IO methods
using Parent_t:: read ;
using Parent_t:: write;
// --- IO methods
using Parent_t:: read ;
using Parent_t:: write;

std::istream& read (std::istream& is)
{ std::string s; return Parent_t::read(is >> s)>> s >> _irred >> s >> _exponent; }
{ std::string s; return Parent_t::read(is >> s)>> s >> Parent_t::residu() >> s >> _exponent; }
std::ostream& write(std::ostream& os) const
{ return Parent_t::write(os<<"Local- ") << "irred: " << _irred << ", exponent: " << _exponent; }

Element& gcdin (Element& a, const Element& b) const
{
Givaro::gcd(a, a, b);
Givaro::gcd(a, a, Parent_t::residu());
return reduce(a);
}

{ return Parent_t::write(os<<"Local- ") << "irred: " << Parent_t::residu() << ", exponent: " << _exponent; }

Element& gcdin (Element& a, const Element& b) const {
Givaro::gcd(a, a, b);
Givaro::gcd(a, a, Parent_t::residu());
return reduce(a);
}

Element& gcd(Element& g, const Element& a, const Element& b) const {
return Givaro::gcd(g,a,b);
}

Element& xgcd(Element& g, Element& s, Element& t, const Element& a, const Element& b) const {
return Givaro::gcd(g,s,t,a,b);
}


bool isUnit(const Element& a) const {
Element g;
Givaro::gcd(g, a, Parent_t::residu());
return isOne(g);
}

bool isDivisor(const Element& a, const Element& b) const {
Integer g;
if (this->isZero(a)) return false;
else if (this->isZero(b)) return true;
else {
Givaro::gcd(g, Integer(a), Integer(Parent_t::residu()));
return this->isZero(Integer(b) % g);
}
}
Element& div(Element& r, const Element& a, const Element& b) const {
Integer g, ia(a), ib(b);
Givaro::gcd(g, ia, ib);
ia /= g;
ib /= g;
Element iv;
this->inv(iv, ib);
return this->mul(r, ia, iv);
}

Element& divin(Element& r, const Element& b) const {
Element a(r); // @enhancement: could use inplace variants
return div(r,a,b);
}

Element& normal (Element& a, const Element& b) const
{

if (this->isZero(b))
return a = this->zero;
else
return this->gcd(a, b, Parent_t::residu() );
}

Element& normalIn (Element& a) const
{
if (this->isZero(a))
return a;
else
return this->gcdin(a, Parent_t::residu() );
}

bool isUnit(const Element& a) const
{
Element g;
Givaro::gcd(g, a, _irred);
return isOne(g);
}


protected:
intType _irred;
uint32_t _exponent;
};

Expand Down
18 changes: 0 additions & 18 deletions linbox/solutions/smith-form.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -30,12 +30,8 @@
#include <iterator>
#include "linbox/util/error.h"
#include "linbox/algorithms/matrix-hom.h"
//#ifdef __LINBOX_HAVE_NTL
#include "linbox/algorithms/smith-form-adaptive.h"
//#endif
#include "givaro/zring.h"
//#include "linbox/algorithms/smith-form.h"
//#include "linbox/algorithms/smith-form-local.h"

namespace LinBox
{
Expand Down Expand Up @@ -139,18 +135,6 @@ namespace LinBox
return V;
}

#if 0
// for specialization with respect to the DomainCategory
template< class Output, class Blackbox, class SmithMethod, class DomainCategory>
Output &smithForm(Output & S,
const Blackbox &A,
const DomainCategory &tag,
const SmithMethod &M)
{
throw LinBoxError( "Smith form solution implemented only for NTL.\n Please reconfigure LinBox with NTL enabled.");
}

#endif
// The smithForm with default Method
template<class Blackbox>
SmithList<typename Blackbox::Field> &
Expand Down Expand Up @@ -214,8 +198,6 @@ namespace LinBox
}
#endif

//#ifdef __LINBOX_HAVE_NTL

/*template<>
std::list<std::pair<integer, size_t> > &
smithForm(std::list<std::pair<integer, size_t> >& S,
Expand Down
2 changes: 1 addition & 1 deletion tests/Makefile.am
View file
Original file line number Diff line number Diff line change
Expand Up @@ -79,7 +79,6 @@ BASIC_TESTS = \
test-smith-form-adaptive \
test-smith-form-iliopoulos \
test-smith-form-local \
test-one-invariant-factor \
test-last-invariant-factor \
test-qlup \
test-det \
Expand Down Expand Up @@ -136,6 +135,7 @@ FULLCHECK_TESTS = \
test-hom \
test-image-field \
test-inverse \
test-one-invariant-factor \
test-la-block-lanczos \
test-matpoly-mult \
test-matrix-domain \
Expand Down
Loading

0 comments on commit f9ef2c7

Please sign in to comment.