ThePEG/Config/PhysicalQtyComplex.h

0001 // -*- C++ -*-
0002 //
0003 // PhysicalQtyComplex.h is a part of ThePEG - Toolkit for HEP Event Generation
0004 // Copyright (C) 2006-2019 David Grellscheid, Leif Lonnblad
0005 //
0006 // ThePEG is licenced under version 3 of the GPL, see COPYING for details.
0007 // Please respect the MCnet academic guidelines, see GUIDELINES for details.
0008 //
0009 #ifndef Physical_Qty_Complex_H
0010 #define Physical_Qty_Complex_H
0011 #include "PhysicalQty.h"
0012 #include "PhysicalQtyOps.h"
0013 #include <complex>
0014
0015 /** @file PhysicalQtyComplex.h
0016  * Overloads for operations on complex physical quantities.
0017  */
0018
0019 namespace std {
0020   /**
0021    *  Template specialization for std::complex<Qty<0,0,0> >
0022    *  with conversions to complex<double>
0023    */
0024   template<>
0025   class complex<ThePEG::QtyDouble>
0026   {
0027   public:
0028     /// Default constructor
0029     constexpr complex(double r=0.0, double i=0.0)
0030       : rawValue_(r,i) {}
0031
0032     /// Constructor from complex<double>
0033     constexpr complex(complex<double> C)
0034       : rawValue_(C) {}
0035
0036     /**
0037      * The internal representation of the dimensionful quantity.
0038      * Using this will break dimension-consistency.
0039      */
0040     constexpr complex<double> rawValue() const { return rawValue_; }
0041
0042     /// Real part
0043     constexpr double real() const { return rawValue_.real(); }
0044
0045     /// Imaginary part
0046     constexpr double imag() const { return rawValue_.imag(); }
0047
0048     /// Cast to complex<double>
0049     constexpr operator complex<double>() const {
0050       return rawValue_;
0051     }
0052
0053     /// Addition-assignment
0054     complex<ThePEG::QtyDouble> &
0055     operator+=(const complex<ThePEG::QtyDouble> x) {
0056       rawValue_ += x.rawValue();
0057       return *this;
0058     }
0059
0060     /// Subtraction-assignment
0061     complex<ThePEG::QtyDouble> &
0062     operator-=(const complex<ThePEG::QtyDouble> x) {
0063       rawValue_ -= x.rawValue();
0064       return *this;
0065     }
0066
0067   private:
0068     /// Internal value of the dimensioned quantity
0069     complex<double> rawValue_;
0070   };
0071 }
0072 // =========================================
0073
0074 namespace ThePEG {
0075
0076 /// @name Overloads for mathematical operations
0077 //@{
0078 // complex qty = complex qty * complex qty
0079 template<typename L1, typename E1, typename Q1,
0080          typename L2, typename E2, typename Q2>
0081 inline constexpr auto
0082 operator*(std::complex<Qty<L1,E1,Q1>> q1,
0083           std::complex<Qty<L2,E2,Q2>> q2)
0084 -> std::complex<decltype(q1.real()*q2.real())>
0085 {
0086   return {q1.real()*q2.real() - q1.imag()*q2.imag(),
0087           q1.real()*q2.imag() + q1.imag()*q2.real()};
0088 }
0089
0090 // complex qty = complex qty * complex qty
0091 template<typename L, typename E, typename Q>
0092 inline constexpr std::complex<typename Qty<L,E,Q>::Squared>
0093 operator*(std::complex<Qty<L,E,Q>> q1,
0094           std::complex<Qty<L,E,Q>> q2)
0095 {
0096   return {q1.real()*q2.real() - q1.imag()*q2.imag(),
0097           q1.real()*q2.imag() + q1.imag()*q2.real()};
0098 }
0099
0100 // complex qty = complex double - complex qty
0101 inline constexpr std::complex<double>
0102 operator-(std::complex<double> q1, std::complex<QtyDouble> q2) {
0103   return {q1.real()-q2.real(), q1.imag()-q2.imag()};
0104 }
0105
0106 // complex qty = complex double + complex qty
0107 inline constexpr std::complex<double>
0108 operator+(std::complex<double> q1, std::complex<QtyDouble> q2) {
0109   return {q1.real()+q2.real(), q1.imag()+q2.imag()};
0110 }
0111
0112 // complex qty = complex double * complex qty
0113 template<typename L, typename E, typename Q>
0114 inline constexpr std::complex<Qty<L,E,Q>>
0115 operator*(std::complex<double> q1, std::complex<Qty<L,E,Q>> q2) {
0116   return {q1.real()*q2.real() - q1.imag()*q2.imag(),
0117           q1.real()*q2.imag() + q1.imag()*q2.real()};
0118 }
0119
0120 // complex qty = complex double / complex qty
0121 template<typename L, typename E, typename Q>
0122 inline std::complex<typename Qty<L,E,Q>::Inverse>
0123 operator/(std::complex<double> q1, std::complex<Qty<L,E,Q>> q2) {
0124   auto tmp  =  q1*conj(q2);
0125   auto norm = (q2*conj(q2)).real();
0126   return {tmp.real()/norm, tmp.imag()/norm};
0127 }
0128
0129 // complex qty = complex double / qty
0130 template<typename L, typename E, typename Q>
0131 inline constexpr std::complex<typename Qty<L,E,Q>::Inverse>
0132 operator/(std::complex<double> q1, Qty<L,E,Q> q2) {
0133   return {q1.real()/q2, q1.imag()/q2};
0134 }
0135
0136 // complex qty = complex qty / complex double
0137 template<typename L, typename E, typename Q>
0138 inline std::complex<Qty<L,E,Q>>
0139 operator/(std::complex<Qty<L,E,Q>> q1, std::complex<double> q2) {
0140   auto tmp  =  q1*conj(q2);
0141   auto norm = (q2*conj(q2)).real();
0142   return {tmp.real()/norm, tmp.imag()/norm};
0143 }
0144
0145 // complex qty = qty / complex double
0146 template<typename L, typename E, typename Q>
0147 inline std::complex<Qty<L,E,Q>>
0148 operator/(Qty<L,E,Q> q1, std::complex<double> q2) {
0149   auto tmp  =  q1*conj(q2);
0150   auto norm = (q2*conj(q2)).real();
0151   return {tmp.real()/norm, tmp.imag()/norm};
0152 }
0153
0154 // complex double = complex qty / complex qty
0155 template<typename L, typename E, typename Q>
0156 inline std::complex<double>
0157 operator/(std::complex<Qty<L,E,Q>> q1,
0158           std::complex<Qty<L,E,Q>> q2) {
0159   auto tmp  =  q1*conj(q2);
0160   auto norm = (q2*conj(q2)).real();
0161   return {tmp.real()/norm, tmp.imag()/norm};
0162 }
0163
0164 // complex double = qty / complex qty
0165 template<typename L, typename E, typename Q>
0166 inline std::complex<double>
0167 operator/(Qty<L,E,Q> q1, std::complex<Qty<L,E,Q>> q2) {
0168   auto tmp = q1*conj(q2);
0169   auto norm = (q2*conj(q2)).real();
0170   return {tmp.real()/norm, tmp.imag()/norm};
0171 }
0172
0173 // complex double = complex qty / qty
0174 template<typename L, typename E, typename Q>
0175 inline constexpr std::complex<double>
0176 operator/(std::complex<Qty<L,E,Q>> q1, Qty<L,E,Q> q2) {
0177   return {q1.real()/q2, q1.imag()/q2};
0178 }
0179
0180 // complex qty = complex qty / complex qty
0181 template<typename L1, typename E1, typename Q1,
0182          typename L2, typename E2, typename Q2>
0183 inline auto
0184 operator/(std::complex<Qty<L1,E1,Q1>> q1,
0185           std::complex<Qty<L2,E2,Q2>> q2)
0186 -> std::complex<decltype(q1.real()/q2.real())>
0187 {
0188   auto  tmp =  q1*conj(q2);
0189   auto norm = (q2*conj(q2)).real();
0190   return {tmp.real()/norm, tmp.imag()/norm};
0191 }
0192
0193 // complex qty = qty / complex qty
0194 template<typename L1, typename E1, typename Q1,
0195          typename L2, typename E2, typename Q2>
0196 inline auto
0197 operator/(Qty<L1,E1,Q1> q1,
0198           std::complex<Qty<L2,E2,Q2>> q2)
0199 -> std::complex<decltype(q1/q2.real())>
0200 {
0201   auto  tmp =  q1*conj(q2);
0202   auto norm = (q2*conj(q2)).real();
0203   return {tmp.real()/norm, tmp.imag()/norm};
0204 }
0205
0206 // complex qty = complex qty / qty
0207 template<typename L1, typename E1, typename Q1,
0208          typename L2, typename E2, typename Q2>
0209 inline constexpr auto
0210 operator/(std::complex<Qty<L1,E1,Q1>> q1, Qty<L2,E2,Q2> q2)
0211 -> std::complex<decltype(q1.real()/q2)>
0212 {
0213   return {q1.real()/q2, q1.imag()/q2};
0214 }
0215
0216
0217 // complex qty = complex qty * complex double
0218 template<typename L, typename E, typename Q>
0219 inline constexpr std::complex<Qty<L,E,Q>>
0220 operator*(std::complex<Qty<L,E,Q>> q1, std::complex<double> q2) {
0221   return q2 * q1;
0222 }
0223
0224
0225 // complex qty = qty * complex qty
0226 template<typename L1, typename E1, typename Q1,
0227          typename L2, typename E2, typename Q2>
0228 inline constexpr auto
0229 operator*(Qty<L1,E1,Q1> q1, std::complex<Qty<L2,E2,Q2>> q2)
0230 -> std::complex<decltype(q1*q2.real())>
0231 {
0232   return {q1*q2.real(), q1*q2.imag()};
0233 }
0234
0235 // complex qty = qty * complex qty
0236 template<typename L, typename E, typename Q>
0237 inline constexpr std::complex<typename Qty<L,E,Q>::Squared>
0238 operator*(Qty<L,E,Q> q1, std::complex<Qty<L,E,Q>> q2) {
0239   return {q1*q2.real(), q1*q2.imag()};
0240 }
0241
0242 // complex qty = qty * complex double
0243 template<typename L, typename E, typename Q>
0244 inline constexpr std::complex<Qty<L,E,Q>>
0245 operator*(Qty<L,E,Q> q1, std::complex<double> q2) {
0246   return {q1*q2.real(), q1*q2.imag()};
0247 }
0248
0249 // complex qty = complex double * qty
0250 template<typename L, typename E, typename Q>
0251 inline constexpr std::complex<Qty<L,E,Q>>
0252 operator*(std::complex<double> q1, Qty<L,E,Q> q2) {
0253   return q2 * q1;
0254 }
0255
0256
0257 // complex qty = complex qty * qty
0258 template<typename L1, typename E1, typename Q1,
0259          typename L2, typename E2, typename Q2>
0260 inline constexpr auto
0261 operator*(std::complex<Qty<L1,E1,Q1>> q1, Qty<L2,E2,Q2> q2)
0262 -> decltype(q2*q1)
0263 {
0264   return q2 * q1;
0265 }
0266
0267 // complex qty = complex qty * qty
0268 template<typename L, typename E, typename Q>
0269 inline constexpr std::complex<typename Qty<L,E,Q>::Squared>
0270 operator*(std::complex<Qty<L,E,Q>> q1, Qty<L,E,Q> q2) {
0271   return q2 * q1;
0272 }
0273
0274 // complex qty *= double
0275 template<typename L, typename E, typename Q>
0276 inline constexpr std::complex<Qty<L,E,Q>> &
0277 operator*=(std::complex<Qty<L,E,Q>> & q1, double q2) {
0278   return (q1 = q1 * q2);
0279 }
0280
0281 // complex qty /= double
0282 template<typename L, typename E, typename Q>
0283 inline constexpr std::complex<Qty<L,E,Q>> &
0284 operator/=(std::complex<Qty<L,E,Q>> & q1, double q2) {
0285   return (q1 = q1 / q2);
0286 }
0287 //@}
0288 }
0289
0290 #endif