SimplifyProduct.h

Go to the documentation of this file.

#ifndef PhysicsTools_Utilities_SimplifyProduct_h
#define PhysicsTools_Utilities_SimplifyProduct_h

#include "PhysicsTools/Utilities/interface/Product.h"
#include "PhysicsTools/Utilities/interface/Fraction.h"
#include "PhysicsTools/Utilities/interface/DecomposePower.h"
#include "PhysicsTools/Utilities/interface/ParametricTrait.h"
#include <boost/mpl/if.hpp>
#include <type_traits>

#include "PhysicsTools/Utilities/interface/Simplify_begin.h"

namespace funct {

  // a * ( b * c ) = ( a * b ) * c
  PROD_RULE(TYPT3, A, PROD_S(B, C), PROD(PROD(A, B), C), (_1 * _2._1) * _2._2);

  // 0 * a = 0
  PROD_RULE(TYPT1, NUM(0), A, NUM(0), num<0>());

  // 0 * n = 0
  PROD_RULE(TYPN1, NUM(0), NUM(n), NUM(0), num<0>());

  // 0 * ( a * b ) => 0
  PROD_RULE(TYPT2, NUM(0), PROD_S(A, B), NUM(0), num<0>());

  // 1 * a = a
  PROD_RULE(TYPT1, NUM(1), A, A, _2);

  // avoid template ambiguities
  // 1 * n = n
  PROD_RULE(TYPN1, NUM(1), NUM(n), NUM(n), _2);

  // 1 * (n/m) = (n/m)
  PROD_RULE(TYPN2, NUM(1), FRACT_S(n, m), FRACT_S(n, m), _2);

  // 1 * 1 = 1
  PROD_RULE(TYP0, NUM(1), NUM(1), NUM(1), num<1>());

  // ( - 1 ) * a = - a
  PROD_RULE(TYPT1, NUM(-1), A, MINUS_S(A), -_2);

  // ( - 1 ) * n = -n
  PROD_RULE(TYPN1, NUM(-1), NUM(n), NUM(-n), num<-n>());

  // 1 * ( a * b ) => ( a * b )
  PROD_RULE(TYPT2, NUM(1), PROD_S(A, B), PROD_S(A, B), _2);

  // a * ( -b ) =  - ( a * b )
  PROD_RULE(TYPT2, A, MINUS_S(B), MINUS(PROD(A, B)), -(_1* _2._));

  // n * ( -a ) =  - ( n * a )
  PROD_RULE(TYPN1T1, NUM(n), MINUS_S(A), MINUS(PROD(NUM(n), A)), -(_1* _2._));

  // ( a * b ) * ( -c )=  - ( ( a * b ) * c )
  PROD_RULE(TYPT3, PROD_S(A, B), MINUS_S(C), MINUS(PROD(PROD(A, B), C)), -(_1* _2._));

  // 1 * ( -a ) = -a
  PROD_RULE(TYPT1, NUM(1), MINUS_S(A), MINUS(A), _2);

  // ( - a ) * ( - b ) = a * b
  PROD_RULE(TYPT2, MINUS_S(A), MINUS_S(B), PROD(A, B), _1._* _2._);

  // ( -a ) * b = -( a * b )
  PROD_RULE(TYPT2, MINUS_S(A), B, MINUS(PROD(A, B)), -(_1._* _2));

  // a * ( b / c ) = ( a * b ) / c
  PROD_RULE(TYPT3, A, RATIO_S(B, C), RATIO(PROD(A, B), C), (_1 * _2._1) / _2._2);

  // n * ( a / b ) = ( n * a ) / b
  PROD_RULE(TYPN1T2, NUM(n), RATIO_S(A, B), RATIO(PROD(NUM(n), A), B), (_1 * _2._1) / _2._2);

  // 1 * ( a / b ) = a / b
  PROD_RULE(TYPT2, NUM(1), RATIO_S(A, B), RATIO(A, B), _2);

  // 0 * ( a / b ) = 0
  PROD_RULE(TYPT2, NUM(0), RATIO_S(A, B), NUM(0), num<0>());

  // a * n = n * a
  PROD_RULE(TYPN1T1, A, NUM(n), PROD(NUM(n), A), _2* _1);

  // ( a * b ) n = ( n * a ) * b
  PROD_RULE(TYPN1T2, PROD_S(A, B), NUM(n), PROD(PROD(NUM(n), A), B), (_2 * _1._1) * _1._2);

  // ( a * b ) * ( c * d ) => ( ( a * b ) * c ) * d
  PROD_RULE(TYPT4, PROD_S(A, B), PROD_S(C, D), PROD(PROD(PROD(A, B), C), D), (_1 * _2._1) * _2._2);

  // n/m * ( a / k ) = n/(m+k) * a
  PROD_RULE(TYPN3T1, FRACT_S(n, m), RATIO_S(A, NUM(k)), PROD(FRACT(n, m + k), A), (fract<n, m + k>() * _2._1));

  // ( a / b ) * n = ( n a ) / b
  PROD_RULE(TYPN1T2, RATIO_S(A, B), NUM(n), RATIO(PROD(NUM(n), A), B), (_2 * _1._1) / _1._2);

  // ( a / b ) * c = ( a * c ) / b
  PROD_RULE(TYPT3, RATIO_S(A, B), C, RATIO(PROD(A, C), B), (_1._1 * _2) / _1._2);

  // 0 * 1 = 0  ( avoid template ambiguity )
  PROD_RULE(TYP0, NUM(0), NUM(1), NUM(0), num<0>());

  // ( a / b ) * ( c / d )= a * c / ( b * d )
  PROD_RULE(TYPT4, RATIO_S(A, B), RATIO_S(C, D), RATIO(PROD(A, C), PROD(B, D)), (_1._1 * _2._1) / (_1._2 * _2._2));

  // a^b * a^c => a^( b + c )
  template <TYPT3, bool parametric = Parametric<A>::value>
  struct SimplifyPowerProduct {
    typedef POWER(A, B) arg1;
    typedef POWER(A, C) arg2;
    typedef PROD_S(arg1, arg2) type;
    COMBINE(arg1, arg2, type(_1, _2));
  };

  TEMPL(T3)
  struct SimplifyPowerProduct<A, B, C, false> {
    typedef POWER(A, B) arg1;
    typedef POWER(A, C) arg2;
    typedef POWER(A, SUM(B, C)) type;
    inline static type combine(const arg1& _1, const arg2& _2) {
      return pow(DecomposePower<A, B>::getBase(_1),
                 (DecomposePower<A, B>::getExp(_1) + DecomposePower<A, C>::getExp(_2)));
    }
  };

  TEMPL(T3) struct Product<POWER_S(A, B), POWER_S(A, C)> : public SimplifyPowerProduct<A, B, C> {};

  TEMPL(T2) struct Product<POWER_S(A, B), POWER_S(A, B)> : public SimplifyPowerProduct<A, B, B> {};

  TEMPL(T2) struct Product<A, POWER_S(A, B)> : public SimplifyPowerProduct<A, NUM(1), B> {};

  TEMPL(N1T1) struct Product<A, POWER_S(A, NUM(n))> : public SimplifyPowerProduct<A, NUM(1), NUM(n)> {};

  TEMPL(T2) struct Product<POWER_S(A, B), A> : public SimplifyPowerProduct<A, B, NUM(1)> {};

  TEMPL(N1T1) struct Product<POWER_S(A, NUM(n)), A> : public SimplifyPowerProduct<A, NUM(n), NUM(1)> {};

  TEMPL(T1) struct Product<A, A> : public SimplifyPowerProduct<A, NUM(1), NUM(1)> {};

  TEMPL(T2) struct Product<PROD_S(A, B), PROD_S(A, B)> : public SimplifyPowerProduct<PROD(A, B), NUM(1), NUM(1)> {};

  TEMPL(T1) struct Product<MINUS_S(A), MINUS_S(A)> : public SimplifyPowerProduct<MINUS_S(A), NUM(1), NUM(1)> {};

  // n * n = n ^ 2
  PROD_RULE(TYPN1, NUM(n), NUM(n), NUM(n* n), num<n * n>());

  // a/ b * ( c * d ) = ( a * c * d ) / b
  PROD_RULE(TYPT4, RATIO_S(A, B), PROD_S(C, D), RATIO(PROD(PROD(A, C), D), B), ((_1._1 * _2._1) * _2._2) / _1._2);

  // simplify f * g * h regardless of the order
  template <typename Prod, bool simplify = Prod::value>
  struct AuxProduct {
    typedef PROD(typename Prod::AB, typename Prod::C) type;
    COMBINE(typename Prod::AB, typename Prod::C, _1* _2);
  };

  template <typename Prod>
  struct AuxProduct<Prod, false> {
    typedef PROD_S(typename Prod::AB, typename Prod::C) type;
    COMBINE(typename Prod::AB, typename Prod::C, type(_1, _2));
  };

  template <typename F, typename G, typename H>
  struct Product<PROD_S(F, G), H> {
    struct prod0 {
      typedef F A;
      typedef G B;
      typedef H C;
      typedef PROD_S(A, B) AB;
      inline static const A& a(const F& f, const G& g, const H& h) { return f; }
      inline static const B& b(const F& f, const G& g, const H& h) { return g; }
      inline static const C& c(const F& f, const G& g, const H& h) { return h; }
      enum { value = false };
    };
    struct prod1 {
      typedef F A;
      typedef H B;
      typedef G C;
      typedef PROD_S(A, B) base;
      typedef PROD(A, B) AB;
      inline static const A& a(const F& f, const G& g, const H& h) { return f; }
      inline static const B& b(const F& f, const G& g, const H& h) { return h; }
      inline static const C& c(const F& f, const G& g, const H& h) { return g; }
      enum { value = not::std::is_same<AB, base>::value };
    };
    struct prod2 {
      typedef G A;
      typedef H B;
      typedef F C;
      typedef PROD_S(A, B) base;
      typedef PROD(A, B) AB;
      inline static const A& a(const F& f, const G& g, const H& h) { return g; }
      inline static const B& b(const F& f, const G& g, const H& h) { return h; }
      inline static const C& c(const F& f, const G& g, const H& h) { return f; }
      enum { value = not::std::is_same<AB, base>::value };
    };

    typedef typename ::boost::mpl::if_<prod1, prod1, typename ::boost::mpl::if_<prod2, prod2, prod0>::type>::type prod;
    typedef typename AuxProduct<prod>::type type;
    inline static type combine(const ProductStruct<F, G>& fg, const H& h) {
      const F& f = fg._1;
      const G& g = fg._2;
      const typename prod::A& a = prod::a(f, g, h);
      const typename prod::B& b = prod::b(f, g, h);
      const typename prod::C& c = prod::c(f, g, h);
      return AuxProduct<prod>::combine(a * b, c);
    }
  };

}  // namespace funct

#include "PhysicsTools/Utilities/interface/Simplify_end.h"

#endif