libalgebra-demo/libalgebra/sparse_vector.h

/* *************************************************************

Copyright 2010 Terry Lyons, Stephen Buckley, Djalil Chafai, 
Greg Gyurkó and Arend Janssen. 

Distributed under the terms of the GNU General Public License, 
Version 3. (See accompanying file License.txt)

************************************************************* */




//  sparse_vector.h


// Include once wrapper
#ifndef DJC_COROPA_LIBALGEBRA_SPARSEVECTORH_SEEN
#define DJC_COROPA_LIBALGEBRA_SPARSEVECTORH_SEEN


template<class BASIS, class MAP>
class sparse_vector :
        /*private*/ MAP
        {
public:
                using MAP::begin;
                using MAP::end;
                using MAP::find;
                using MAP::erase;
                using MAP::operator [];
                using MAP::clear;
        
        void swap (typename sparse_vector & rhs) {
                        MAP::swap(rhs);
                }
                typedef typename BASIS::RATIONAL RATIONAL;
                typedef typename MAP::key_type KEY;
                typedef typename MAP::mapped_type SCALAR;
                typedef typename MAP::iterator iterator;
                typedef typename MAP::const_iterator const_iterator;
public:
                const static SCALAR zero; //0
                const static SCALAR one;  //+1
                const static SCALAR mone; //-1
        
                static BASIS basis;

public:

        sparse_vector(void) {}
        sparse_vector(const sparse_vector& v) : MAP(v) {}

        explicit sparse_vector(const KEY& k, const SCALAR& s = one)
                {
                        if (s != zero)
                                (*this)[k] = s;
                }
public:
        inline sparse_vector operator-(void) const
                {       
                        if (empty())
                                return *this;
                        const_iterator in;
                        sparse_vector result;
                        for (in = begin(); in != end(); ++in)
                                result[in->first] = -(in->second);
                        return result;
                }
        inline sparse_vector& operator*=(const SCALAR& s)
                {
                        if (s != zero)
                        {
                                iterator it;
                                if (!empty())
                                        for (it = begin(); it != end(); ++it)
                                                it->second *= s;
                        }
                        else
                                clear();
                        return *this;
                }
        inline __DECLARE_BINARY_OPERATOR(sparse_vector,*,*=,SCALAR)
        inline sparse_vector& operator/=(const RATIONAL& s)
                {
                        iterator it;
                        if (!empty())
                                for (it = begin(); it != end(); ++it)
                                {
                                        RATIONAL temp(1);
                                        it->second *= (temp / s);
                                }
                        return *this;
                }
        inline __DECLARE_BINARY_OPERATOR(sparse_vector,/,/=,RATIONAL)
        inline sparse_vector& operator+=(const sparse_vector& rhs)
                {
                        iterator it;
                        const_iterator cit;
                        if (rhs.empty())
                                return *this;
                        if (empty())
                                return *this = rhs;
                        for (cit = rhs.begin(); cit != rhs.end(); ++cit)
                        { // Instead of a bare (*this)[cit->first] += cit->second;
                                it = find(cit->first);
                                if (it == end())
                                        (*this)[cit->first] = cit->second;
                                else if ((it->second += cit->second) == zero)
                                        erase(it->first);
                        }       
                        return *this;
                }
        inline __DECLARE_BINARY_OPERATOR(sparse_vector,+,+=,sparse_vector)
        inline sparse_vector& operator-=(const sparse_vector& rhs)
                {
                        iterator it;
                        const_iterator cit;
                        if (rhs.empty())
                                return *this;
                        if (empty())
                                return *this = -rhs;
                        for (cit = rhs.begin(); cit != rhs.end(); ++cit)
                        { // Instead of a bare (*this)[cit->first] -= cit->second;
                                it = find(cit->first);
                                if (it == end())
                                        (*this)[cit->first] = -(cit->second);
                                else if ((it->second -= cit->second) == zero)
                                        erase(it->first);
                        }       
                        return *this;
                }
        inline __DECLARE_BINARY_OPERATOR(sparse_vector,-,-=,sparse_vector)

                
        inline sparse_vector& add_scal_prod(const sparse_vector& rhs,
                const SCALAR& s)
                {
                        iterator it;
                        const_iterator cit;
                        if ((s == zero) || rhs.empty())
                                return *this;
                        if (empty())
                        {
                                *this = rhs;
                                return operator *= (s);
                        }
                        for (cit = rhs.begin(); cit != rhs.end(); ++cit)
                        { // Instead of a bare (*this)[cit->first] += cit->second * s;
                                it = find(cit->first);
                                if (it == end())
                                        (*this)[cit->first] = cit->second * s;
                                else if ((it->second += cit->second * s) == zero)
                                        erase(it->first);
                        }       
                        return *this;
                }
        inline sparse_vector& sub_scal_prod(const sparse_vector& rhs,
                const SCALAR& s)
                {
                        iterator it;
                        const_iterator cit;
                        if ((s == zero) || rhs.empty())
                                return *this;
                        if (empty())
                        {
                                *this = rhs;
                                return operator *= (-s);
                        }
                        for (cit = rhs.begin(); cit != rhs.end(); ++cit)
                        { // Instead of a bare (*this)[cit->first] -= cit->second * s;
                                it = find(cit->first);
                                if (it == end())
                                        (*this)[cit->first] = cit->second * -s;
                                else if ((it->second -= cit->second * s) == zero)
                                        erase(it->first);
                        }       
                        return *this;
                }
        inline sparse_vector& add_scal_div(const sparse_vector& rhs,
                const RATIONAL& s)
                {
                        iterator it;
                        const_iterator cit;
                        if (rhs.empty())
                                return *this;
                        if (empty())
                        {
                                *this = rhs;
                                return operator /= (s);
                        }
                        for (cit = rhs.begin(); cit != rhs.end(); ++cit)
                        { // Instead of a bare (*this)[cit->first] += cit->second / s;
                                it = find(cit->first);
                                if (it == end())
                                        (*this)[cit->first] = cit->second / s;
                                else if ((it->second += (cit->second / s)) == zero)
                                        erase(it->first);
                        }       
                        return *this;
                }
        inline sparse_vector& sub_scal_div(const sparse_vector& rhs,
                const RATIONAL& s)
                {
                        iterator it;
                        const_iterator cit;
                        if (rhs.empty())
                                return *this;
                        if (empty())
                        {
                                *this = rhs;
                                return operator /= (-s);
                        }
                        for (cit = rhs.begin(); cit != rhs.end(); ++cit)
                        { // Instead of a bare (*this)[cit->first] -= cit->second / s;
                                it = find(cit->first);
                                if (it == end())
                                        (*this)[cit->first] = - cit->second / s;
                                else if ((it->second -= (cit->second / s)) == zero)
                                        erase(it->first);
                        }       
                        return *this;
                }
        bool operator==(const sparse_vector& rhs) const
                {
                        if (size() != rhs.size())
                                return false;
                        const_iterator i, j;
                        for (i = begin(); i != end(); ++i)
                        {
                                j = rhs.find(i->first);
                                if ((j == rhs.end()) || (j->second != i->second))
                                        return false;
                        }
                        return true;
                }
        bool operator < (const sparse_vector& rhs) const
                {
                        return std::lexicographical_compare(begin(), end(), rhs.begin(), rhs.end());
                }
        bool operator!=(const sparse_vector& rhs) const
                {
                        return !operator == (rhs);
                }
        inline SCALAR NormL1() const
                {
                        const_iterator i;
                        SCALAR ans(zero);
                        for (i = begin(); i != end(); ++i)
                        {
                                ans += abs(i->second);
                        }
                        return ans;
                }
        inline SCALAR NormL1(const DEG & d) const
                {
                        const_iterator i;
                        SCALAR ans(zero);
                        for (i = begin(); i != end(); ++i)
                        {
                                if (d == basis.degree(i->first))
                                        ans += abs(i->second);
                        }
                        return ans;
                }

        inline friend std::ostream& operator<<(std::ostream& os,
                const sparse_vector& rhs)
                {       
                        const_iterator cit;                     
                        std::pair<BASIS*, KEY> token;
                        token.first = &sparse_vector::basis;
                        os << '{';
                        for (cit = rhs.begin(); cit != rhs.end(); ++cit)
                        {
                                token.second = cit->first;
                                os << ' ' << cit->second << '(' << token << ')';
                        }
                        os << " }";
                        return os;
                }
        };

// Initialisation of static members of sparse_vector<>

template<class BASIS, class MAP>
BASIS sparse_vector<BASIS, MAP>::basis;

template<class BASIS, class MAP>
typename
const MAP::mapped_type sparse_vector<BASIS, MAP>::one(+1);

template<class BASIS, class MAP>
typename
const MAP::mapped_type sparse_vector<BASIS, MAP>::zero(0);

template<class BASIS, class MAP>
typename
const MAP::mapped_type sparse_vector<BASIS, MAP>::mone(-1);

// Include once wrapper
// DJC_COROPA_LIBALGEBRA_SPARSEVECTORH_SEEN
#endif
//EOF.