2012-09-23

SCIP in C++11 ― 2.5.3節その3


代数操作その6:有理関数(不完全版)と問題2.932.95

rationalに対する整数と多項式の場合分けで少し手間がある以外は素直。

----

class RationalArithmetic{
public:
    (中略)

    const List makeRational
    (const List& numerator,const List& denominator)const
    {
        const List g(Generic::greatestCommonDivisor(numerator,denominator));
        if(typeTag(numerator)!="polynomial"
           || typeTag(denominator)!="polynomial"){
            if(typeTag(numerator)=="number"
               && typeTag(denominator)=="number"){
                if(contents(denominator)<makeLeaf(0)){
                    return(makeList
                           (Generic::div(Generic::negate(numerator),g),
                            Generic::div(Generic::negate(denominator),g)));
                }else{
                    return(makeList(Generic::div(numerator,g),
                                    Generic::div(denominator,g)));
                }
            }else{
                return(makeList(numerator,denominator));
            }
        }
        return(makeList(car(Generic::div(numerator,g)),
                        car(Generic::div(denominator,g))));
    }
    (略)
};

//---------abstraction barrier---------

class SparseTermList{
public:
    SparseTermList(const TagType tagIn="sparse"):tagString(tagIn)
    {
        (中略)
        put(makeList("remainder",makeList(this->getTag(),this->getTag())),
            makeLeaf(function<List(List,List)>
                     ([this](const List& x,const List& y)
                      {return(this->tag(this->remainderTerms(x,y)));})));
        put(makeList("gcd",makeList(this->getTag(),this->getTag())),
            makeLeaf(function<List(List,List)>
                     ([this](const List& x,const List& y)
                      {return(this->tag(this->gcdTerms(x,y)));})));
       
    const bool isConstant(const List& termList)const
    {
        return(this->termListOrder(termList)==makeLeaf(0));
    }

    (中略)
    const List remainderTerms(const List& L1, const List& L2)
    {
        return(cadr(this->divTerms(L1,L2)));
    }
   
    const List gcdTerms(const List& L1, const List& L2)
    {
        if(this->isZero(L2)){
            return(L1);
        }else if(this->isConstant(L1)||this->isConstant(L2)){
            return(makeList(this->makeTerm(makeLeaf(0),Generic::makeNumber(1))));
        }
        return(this->gcdTerms(L2,this->remainderTerms(L1,L2)));
    }

   (略)
}:

class DenseTermList{
public:
    DenseTermList(const TagType tagIn="dense"):tagString(tagIn)
    {
        (中略)
        put(makeList("remainder",makeList(this->getTag(),this->getTag())),
            makeLeaf(function<List(List,List)>
                     ([this](const List& x,const List& y)
                      {return(this->tag(this->remainderTerms(x,y)));})));
        put(makeList("gcd",makeList(this->getTag(),this->getTag())),
            makeLeaf(function<List(List,List)>
                     ([this](const List& x,const List& y)
                      {return(this->tag(this->gcdTerms(x,y)));})));
       
(中略)
    const bool isConstant(const List& termList)const
    {
        return(this->termListOrder(termList)==makeLeaf(0));
    }

  
    const List remainderTerms(const List& L1, const List& L2)
    {
        return(cadr(this->divTerms(L1,L2)));
    }
   
    const List gcdTerms(const List& L1, const List& L2)
    {
        if(this->isZero(L2)){
            return(L1);
        }else if(this->isConstant(L1)||this->isConstant(L2)){
            return(makeList(Generic::makeNumber(1)));
        }
        return(this->gcdTerms(L2,this->remainderTerms(L1,L2)));
    }
(略)
};

//---------abstraction barrier---------
class PolynomialArithmetic{
public:
    PolynomialArithmetic(const TagType tagIn="polynomial")
        :tagString(tagIn),
         _sparseTermList(new SparseTermList()),
         _denseTermList(new DenseTermList())
    {
        (略)           put(makeList("remainder",makeList(this->getTag(),this->getTag())),
            makeLeaf(function<List(List,List)>
                     ([this](const List& x,const List& y)
                      {return(this->tag(this->remainderPolynomial(x,y)));})));
        put(makeList("gcd",makeList(this->getTag(),this->getTag())),
            makeLeaf(function<List(List,List)>
                     ([this](const List& x,const List& y)
                      {return(this->tag(this->gcdPolynomial(x,y)));})));
        (略)
    }
   (略)

    const List remainderPolynomial(const List& p1, const List& p2)
    {
        if(this->isSameVariable
           (this->variable(p1),this->variable(p2))){
            return(cadr(this->divPolynomial(p1,p2)));
        }
        cerr<<"Polynomials not in same variable -- REMAINDER-POLY "
            <<listString(p1)<<listString(p2)<<endl;
        exit(1);
        return(makeList());
    }
   
    const List gcdPolynomial(const List& p1, const List& p2)
    {
        if(this->isSameVariable
           (this->variable(p1),this->variable(p2))){
            return(this->makePolynomial
                   (this->variable(p1),
                    Generic::greatestCommonDivisor
                    (this->termList(p1),this->termList(p2))));
        }
        cerr<<"Polynomials not in same variable -- GCD-POLY "
            <<listString(p1)<<listString(p2)<<endl;
        exit(1);
        return(makeList());
    }
(略)
};
//---------abstraction barrier---------


int main(int argc, char** argv)
{
    installNumberPackage();
    installRationalPackage();
    installRealPackage();
    installComplexPackage();
    installPolynomialPackage();
    installCoercion();

    using namespace Generic;

    cout<<"Excersize 2.93:"<<endl;
    const auto p1(makePolynomial
                  ("x",
                   makeSparseTermList
                   (makeList
                    (
                     makeList(2,makeNumber(1)),
                     makeList(0,makeNumber(1))))));
    const auto p2(makePolynomial
                  ("x",
                   makeSparseTermList
                   (makeList
                    (makeList(3,makeNumber(1)),
                     makeList(0,makeNumber(1))))));
    const auto r1(makeRational(p1,p2));
    cout<<"p1 = "<<expressionString(p1)<<endl;
    cout<<"p2 = "<<expressionString(p2)<<endl;
    cout<<"r1 = "<<expressionString(r1)<<endl;
    cout<<"(add r1 r1) = "<<expressionString(add(r1,r1))<<endl;
   
   
    cout<<endl<<"Excersize 2.94:"<<endl;
    const auto p3(makePolynomial
                  ("x",
                   makeSparseTermList
                   (makeList
                    (makeList(4,makeNumber(1)),
                     makeList(3,makeNumber(-1)),
                     makeList(2,makeNumber(-2)),
                     makeList(1,makeNumber(2))))));
    const auto p4(makePolynomial
                  ("x",
                   makeSparseTermList
                   (makeList
                    (makeList(3,makeNumber(1)),
                     makeList(1,makeNumber(-1))))));
    cout<<"p3 = "<<expressionString(p3)<<endl;
    cout<<"p4 = "<<expressionString(p4)<<endl;
    cout<<"(greatest-common-divisor p3 p4) = "
        <<expressionString(Generic::greatestCommonDivisor(p3,p4))<<endl;
   
   
    cout<<endl<<"Excersize 2.95:"<<endl;
    const auto p5(makePolynomial
                  ("x",
                   makeSparseTermList
                   (makeList
                    (makeList(2,makeNumber(1)),
                     makeList(1,makeNumber(-2)),
                     makeList(0,makeNumber(1))))));
    const auto p6(makePolynomial
                  ("x",
                   makeSparseTermList
                   (makeList
                    (makeList(2,makeNumber(11)),
                     makeList(0,makeNumber(7))))));
    const auto p7(makePolynomial
                  ("x",
                   makeSparseTermList
                   (makeList
                    (makeList(1,makeNumber(13)),
                     makeList(0,makeNumber(5))))));
    const auto q1(Generic::mul(p5,p6));
    const auto q2(Generic::mul(p5,p7));
    cout<<"p5 = "<<expressionString(p5)<<endl;
    cout<<"p6 = "<<expressionString(p6)<<endl;
    cout<<"p7 = "<<expressionString(p7)<<endl;
    cout<<"q1 = "<<expressionString(q1)<<endl;
    cout<<"q2 = "<<expressionString(q2)<<endl;
    cout<<"(greatest-common-divisor q1 q2) = "
        <<expressionString(Generic::greatestCommonDivisor(q1,q2))<<endl;

    uninstallNumberPackage();
    uninstallRationalPackage();
    uninstallRealPackage();
    uninstallComplexPackage();
    uninstallPolynomialPackage();
    return(0);
}

----
出力
----
Excersize 2.93:
p1 = x^2+1
p2 = x^3+1
r1 = (x^2+1)/(x^3+1)
(add r1 r1) = ((1/2)x^2+(1/2))/((1/4)x^3+(1/4))

Excersize 2.94:
p3 = x^4-x^3-2x^2+2x
p4 = x^3-x
(greatest-common-divisor p3 p4) = -x^2+x

Excersize 2.95:
p5 = x^2-2x+1
p6 = 11x^2+7
p7 = 13x+5
q1 = 11x^4-22x^3+18x^2-14x+7
q2 = 13x^3-21x^2+3x+5
(greatest-common-divisor q1 q2) = (1458/169)x^2-(2916/169)x-(-1458/169)
 

0 件のコメント :

コメントを投稿