00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #if defined(DO_NOT_USE_MDP_COMPLEX)
00014 #include "complex.h"
00015 typedef complex<mdp_real> mdp_complex;
00016 #else
00017
00026 class mdp_complex {
00027 public:
00028 mdp_real re;
00029 mdp_real im;
00030
00031 mdp_real& real() { return re; }
00032 mdp_real& imag() { return im; }
00033
00034 const mdp_real& real() const { return re; }
00035 const mdp_real& imag() const { return im; }
00036
00037 inline mdp_complex(const mdp_real a=0.0, const mdp_real b=0.0) {
00038 re=a;
00039 im=b;
00040 }
00041
00042 inline mdp_complex(const mdp_complex& c) {
00043 re=c.real();
00044 im=c.imag();
00045 }
00046
00047 inline bool operator== (const mdp_complex& c) {
00048 return ((re==c.real()) && (im==c.imag()));
00049 }
00050
00051 inline bool operator!= (const mdp_complex& c) {
00052 return ((re!=c.real()) || (im!=c.imag()));
00053 }
00054
00055 inline friend mdp_real real(const mdp_complex& c) {
00056 return c.real();
00057 }
00058
00059 inline friend mdp_real imag(const mdp_complex& c) {
00060 return c.imag();
00061 }
00062
00063 inline friend mdp_real abs(const mdp_complex& c) {
00064 return sqrt(c.real()*c.real()+c.imag()*c.imag());
00065 }
00066
00067 inline friend mdp_real arg(const mdp_complex& c) {
00068 return phase(c);
00069 }
00070
00071 inline friend mdp_complex pow(const mdp_complex& c, mdp_real z) {
00072 mdp_real r=pow((double) abs(c),(double) z),p=(mdp_real) z*phase(c);
00073 return mdp_complex(r*cos(p),r*sin(p));
00074 }
00075
00076 inline friend mdp_complex sqrt(const mdp_complex& c) {
00077 mdp_real r=sqrt(abs(c)),p=0.5*arg(c);
00078 return mdp_complex(r*cos(p),r*sin(p));
00079 }
00080
00081 inline friend mdp_complex exp(const mdp_complex& c) {
00082 mdp_real exp_c_re=exp(c.real());
00083 return mdp_complex(exp_c_re*cos(c.imag()), exp_c_re*sin(c.imag()));
00084 }
00085
00086 inline friend mdp_complex sin(const mdp_complex& c) {
00087 return mdp_complex(cosh(c.imag())*sin(c.real()), sinh(c.imag())*cos(c.real()));
00088 }
00089
00090 inline friend mdp_complex cos(const mdp_complex& c) {
00091 return mdp_complex(cosh(c.imag())*cos(c.real()), -sinh(c.imag())*sin(c.real()));
00092 }
00093
00094 inline friend mdp_complex times_i(const mdp_complex& c) {
00095 return mdp_complex(-c.imag(), c.real());
00096 }
00097
00098 inline friend mdp_complex times_minus_i(const mdp_complex& c) {
00099 return mdp_complex(c.imag(), -c.real());
00100 }
00101
00102 inline friend mdp_complex operator- (const mdp_complex& c) {
00103 return mdp_complex(-c.real(), -c.imag());
00104 }
00105
00106 inline friend mdp_complex operator+ (const mdp_complex& c) {
00107 return c;
00108 }
00109
00110 inline void operator+= (const mdp_complex& c) {
00111 re+=c.real();
00112 im+=c.imag();
00113 }
00114
00115 inline void operator-= (const mdp_complex& c) {
00116 re-=c.real();
00117 im-=c.imag();
00118 }
00119
00120 inline void operator*= (const mdp_complex& c);
00121
00122 inline void operator/= (const mdp_complex& c);
00123
00124 inline friend mdp_real phase(const mdp_complex& c) {
00125 return atan2(c.imag(), c.real());
00126 }
00127
00128 inline friend mdp_complex conj(const mdp_complex& a) {
00129 return mdp_complex(a.real(), -a.imag());
00130 }
00131
00132 inline void operator+= (const mdp_real c) {
00133 re+=c;
00134 }
00135
00136 inline void operator-= (const mdp_real c) {
00137 re-=c;
00138 }
00139
00140 inline void operator*= (const mdp_real c) {
00141 re*=c;
00142 im*=c;
00143 }
00144
00145 inline void operator/= (const mdp_real c) {
00146 re/=c;
00147 im/=c;
00148 }
00149
00150 };
00151
00152
00153 inline mdp_complex operator+(const mdp_complex& a, const mdp_complex& b) {
00154 return mdp_complex(a.real()+b.real(), a.imag()+b.imag());
00155 }
00156
00157 inline mdp_complex operator-(const mdp_complex& a, const mdp_complex& b) {
00158 return mdp_complex(a.real()-b.real(), a.imag()-b.imag());
00159 }
00160
00161 inline mdp_complex operator*(const mdp_complex& a, const mdp_complex& b) {
00162 return mdp_complex(a.real()*b.real()-a.imag()*b.imag(),a.real()*b.imag()+a.imag()*b.real());
00163 }
00164
00165 inline mdp_complex operator/(const mdp_complex& a, const mdp_complex& b) {
00166 mdp_real den=b.real()*b.real()+b.imag()*b.imag();
00167 return mdp_complex((a.real()*b.real()+a.imag()*b.imag())/den,(a.imag()*b.real()-a.real()*b.imag())/den);
00168 }
00169
00170 inline void mdp_complex::operator*= (const mdp_complex& c) {
00171 (*this)=(*this)*c;
00172 }
00173
00174 inline void mdp_complex::operator/= (const mdp_complex& c) {
00175 (*this)=(*this)/c;
00176 }
00177
00178
00179
00180 inline mdp_complex operator+ (const mdp_complex& a, const int c) {
00181 return mdp_complex(a.real()+c,a.imag());
00182 }
00183
00184 inline mdp_complex operator- (const mdp_complex& a, const int c) {
00185 return mdp_complex(a.real()-c,a.imag());
00186 }
00187
00188 inline mdp_complex operator* (const mdp_complex& a, const int c) {
00189 return mdp_complex(a.real()*c,a.imag()*c);
00190 }
00191
00192 inline mdp_complex operator/ (const mdp_complex& a, const int c) {
00193 return mdp_complex(a.real()/c,a.imag()/c);
00194 }
00195
00196 inline mdp_complex operator+ (const int a, const mdp_complex& c) {
00197 return mdp_complex(c.real()+a,c.imag());
00198 }
00199
00200 inline mdp_complex operator- (const int a, const mdp_complex& c) {
00201 return mdp_complex(c.real()-a,c.imag());
00202 }
00203
00204 inline mdp_complex operator* (const int a, const mdp_complex& c) {
00205 return mdp_complex(c.real()*a,c.imag()*a);
00206 }
00207
00208 inline mdp_complex operator/ (const int a, const mdp_complex& c) {
00209 mdp_real den=(c.real()*c.real()+c.imag()*c.imag())/a;
00210 return mdp_complex((c.real())/den,(-c.imag())/den);
00211 }
00212
00213
00214
00215
00216 inline mdp_complex operator+ (const mdp_complex& a, const long c) {
00217 return mdp_complex(a.real()+c,a.imag());
00218 }
00219
00220 inline mdp_complex operator- (const mdp_complex& a, const long c) {
00221 return mdp_complex(a.real()-c,a.imag());
00222 }
00223
00224 inline mdp_complex operator* (const mdp_complex& a, const long c) {
00225 return mdp_complex(a.real()*c,a.imag()*c);
00226 }
00227
00228 inline mdp_complex operator/ (const mdp_complex& a, const long c) {
00229 return mdp_complex(a.real()/c,a.imag()/c);
00230 }
00231
00232 inline mdp_complex operator+ (const long a, const mdp_complex& c) {
00233 return mdp_complex(c.real()+a,c.imag());
00234 }
00235
00236 inline mdp_complex operator- (const long a, const mdp_complex& c) {
00237 return mdp_complex(c.real()-a,c.imag());
00238 }
00239
00240 inline mdp_complex operator* (const long a, const mdp_complex& c) {
00241 return mdp_complex(c.real()*a,c.imag()*a);
00242 }
00243
00244 inline mdp_complex operator/ (const long a, const mdp_complex& c) {
00245 mdp_real den=(c.real()*c.real()+c.imag()*c.imag())/a;
00246 return mdp_complex((c.real())/den,(-c.imag())/den);
00247 }
00248
00249
00250
00251 inline mdp_complex operator+ (const mdp_complex& a, const float c) {
00252 return mdp_complex(a.real()+c,a.imag());
00253 }
00254
00255 inline mdp_complex operator- (const mdp_complex& a, const float c) {
00256 return mdp_complex(a.real()-c,a.imag());
00257 }
00258
00259 inline mdp_complex operator* (const mdp_complex& a, const float c) {
00260 return mdp_complex(a.real()*c,a.imag()*c);
00261 }
00262
00263 inline mdp_complex operator/ (const mdp_complex& a, const float c) {
00264 return mdp_complex(a.real()/c,a.imag()/c);
00265 }
00266
00267 inline mdp_complex operator+ (const float a, const mdp_complex& c) {
00268 return mdp_complex(c.real()+a,c.imag());
00269 }
00270
00271 inline mdp_complex operator- (const float a, const mdp_complex& c) {
00272 return mdp_complex(c.real()-a,c.imag());
00273 }
00274
00275 inline mdp_complex operator* (const float a, const mdp_complex& c) {
00276 return mdp_complex(c.real()*a,c.imag()*a);
00277 }
00278
00279 inline mdp_complex operator/ (const float a, const mdp_complex& c) {
00280 mdp_real den=(c.real()*c.real()+c.imag()*c.imag())/a;
00281 return mdp_complex((c.real())/den,(-c.imag())/den);
00282 }
00283
00284
00285 inline mdp_complex operator+ (const mdp_complex& a, const double c) {
00286 return mdp_complex(a.real()+c,a.imag());
00287 }
00288
00289 inline mdp_complex operator- (const mdp_complex& a, const double c) {
00290 return mdp_complex(a.real()-c,a.imag());
00291 }
00292
00293 inline mdp_complex operator* (const mdp_complex& a, const double c) {
00294 return mdp_complex(a.real()*c,a.imag()*c);
00295 }
00296
00297 inline mdp_complex operator/ (const mdp_complex& a, const double c) {
00298 return mdp_complex(a.real()/c,a.imag()/c);
00299 }
00300
00301 inline mdp_complex operator+ (const double a, const mdp_complex& c) {
00302 return mdp_complex(c.real()+a,c.imag());
00303 }
00304
00305 inline mdp_complex operator- (const double a, const mdp_complex& c) {
00306 return mdp_complex(c.real()-a,c.imag());
00307 }
00308
00309 inline mdp_complex operator* (const double a, const mdp_complex& c) {
00310 return mdp_complex(c.real()*a,c.imag()*a);
00311 }
00312
00313 inline mdp_complex operator/ (const double a, const mdp_complex& c) {
00314 mdp_real den=(c.real()*c.real()+c.imag()*c.imag())/a;
00315 return mdp_complex((c.real())/den,(-c.imag())/den);
00316 }
00317
00318 #endif
00319
00320
00321 const mdp_complex I=mdp_complex(0,1);
00322
00323 inline mdp_real abs2(const mdp_complex &a) {
00324 return real(a)*real(a)+imag(a)*imag(a);
00325 }
00326
00327 ostream& operator<< (ostream& os, const mdp_complex& a) {
00328 os << a.real() << "+" << a.imag() << "I";
00329 return os;
00330 }
00331
00332
