Re: do you support MTL?

On 5/10/06, Irek Szczesniak <[email protected]> wrote:
>
> Hi Peter,
>
> Thank you for your email and resolving my problems.  I implemented the
> operators for the polynomial type, and now I'm able to multiply
> vectors and matrixes of polynomials.  I'm attaching a complete program
> at the bottom of this email.
>
> Because I have the basic functionality going (a product of polynomial
> vectors and matrixes), I think I will plugh ahead harnessing MTL.
>
> I'm looking forward to the release of the new MTL.
>
>
> Best,
> Irek
>
> Peter Gottschling wrote:
>
> > Hi Irek,
> >
> > Thank you for your interest in MTL and sorry for the trouble you run
> > into.  The problem is that nobody developed MTL for some years.
> > When we restarted it (with different authors), we saw many things
> > principally different.  Therefore, we decided to restart from the
> > beginning (which unfortunately goes slower than expected).
> > Therefore I tend to keep the maintenance efforts on MTL 2.xx minimal
> > (to not delay the delivery of the new version even further).
> >
> > However, the idea of polynomials as scalar types is interesting and
> > I tried to figure out your problems.  There were severals issues and
> > not everything was MTL's fault ;-)
> > 1. print_all_matrix is implemented as std::cout << x ..., it better
> > should be 'using std::cout; cout << x; ...'
> > - there are several such functions and we will fix this one day
> > - one can work around by defining the << operator in std namespace
> > 2. print_vector used iterator -> I changed it to const_iterator in
> > the hope that all relevant vectors have this
> > 3. instead of defining * and + operator for map it works better to
> > define a new type that derives from map (which is needed by
> > 4. anyway)
> > 4. the utilization of scale(v, 0) requires that polynomial is
> > constructible from int or double
> >
> > I add some code that solved the syntactical problems.  The
> > implementation of the numerical part I leave to you although that is
> > the more interesting.  I would appreciate if you can keep me
> > up-to-date occasionally (better directly to me instead of the list).
> >
> > Best,
> > Peter
>
> ********************************************************************
>
> #include <iostream>
> #include <map>
> #include <vector>
> #include <mtl/matrix.h>
> #include <mtl/mtl.h>
> #include <mtl/utils.h>
> #include <mtl/blais.h>
> #include <mtl/sparse1D.h>
>
> using namespace mtl;
> using namespace std;
>
> struct polynomial : map<int, double>
> {
>   typedef map<int, double> base;
>
>   polynomial() : base() {}
>
>   polynomial(double c0) : base()
>   {
>     (*this)[0]= c0;
>   }
> };
>
>
> typedef matrix<polynomial,
>                 rectangle<>,
>                 array< dense<> >,
>                 row_major>::type Mat;
>
>
> polynomial
> operator *(const polynomial &p1, const polynomial &p2)
> {
>   polynomial result;
>
>   for(polynomial::const_iterator i = p1.begin(); i != p1.end(); ++i)
>     for(polynomial::const_iterator j = p2.begin(); j != p2.end(); ++j)
>       result[i->first + j->first] += i->second * j->second;
>
>   return result;
> }
>
>
> polynomial&
> operator +=(polynomial &p1, const polynomial &p2)
> {
>   for(polynomial::const_iterator i = p2.begin(); i != p2.end(); ++i)
>     p1[i->first] += i->second;
>   return p1;
> }
>
>
> polynomial
> operator +(const polynomial &p1, const polynomial &p2)
> {
>   polynomial result = p1;
>   result += p2;
>   return result;
> }
>
>
> namespace std
> {
>   ostream &operator << (ostream &out, const polynomial &p)
>   {
>     for(polynomial::const_reverse_iterator i = p.rbegin(); i !=
> p.rend();)
>       {
>        out << i->second;
>
>        if (i->first)
>          out << "*x^" << i->first;;
>
>        if (++i != p.rend())
>          cout << " + ";
>       }
>
>     return out;
>   }
> }
>
> int
> main()
> {
>   Mat A(3, 3);
>
>   A(0, 0) = polynomial();
>   A(1, 1) = polynomial();
>   A(2, 2) = polynomial();
>
>   A(0, 0)[2] = 1;
>   A(1, 1)[1] = 1;
>   A(2, 2)[0] = 1;
>
>   dense1D<polynomial> c(3), result(3);
>
>   c[0][0] = 1;
>   c[1][1] = 3;
>   c[2][2] = 2;
>
>   mult(A, c, result);
>
>   print_all_matrix(A);
>
>   print_vector(result);
>
>   return 0;
> }
> _______________________________________________
> This list is archived at http://www.osl.iu.edu/MailArchives/mtl-devel/
>



-- 
xubin

_______________________________________________
This list is archived at http://www.osl.iu.edu/MailArchives/mtl-devel/