cmplxvec.cpp 7.2 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400
  1. //#include "stdafx.h"
  2. #include <math.h>
  3. #include "cmplxvec.h"
  4. #include <stdlib.h>
  5. /******************************* constructors ******************************/
  6. CmplxVector::CmplxVector(int n) {
  7. Complex zero(0.,0.);
  8. if (n < 0) {
  9. cerr << "complex vector: negative dimension." << endl;
  10. exit(1);
  11. }
  12. d = n;
  13. if (d > 0) {
  14. v = new Complex[d];
  15. while (n--) v[n] = zero;
  16. }
  17. else v = NULL;
  18. }
  19. CmplxVector::CmplxVector(int n, double a) {
  20. if (n <= 0) {
  21. cerr << "complex vector: undefined dimension." << endl;
  22. exit(1);
  23. }
  24. d = n;
  25. v = new Complex[d];
  26. while (n--) v[n] = a;
  27. }
  28. CmplxVector::CmplxVector(int n, Complex a) {
  29. if (n < 0) {
  30. cerr << "complex vector: negative dimension." << endl;
  31. exit(1);
  32. }
  33. d = n;
  34. v = new Complex[d];
  35. while (n--) v[n] = a;
  36. }
  37. CmplxVector::CmplxVector(const CmplxVector& p) {
  38. d = p.d;
  39. if (d > 0) {
  40. v = new Complex [d];
  41. for(int i = 0; i < d; i++) v[i] = p.v[i];
  42. }
  43. else v = NULL;
  44. }
  45. CmplxVector::CmplxVector(const Vector& v1) {
  46. d = v1.dim();
  47. v = new Complex[d];
  48. for(int i = 0; i < d; i++) v[i] = Complex(v1[i]);
  49. }
  50. /****************************** members ******************************/
  51. void CmplxVector::check_dimensions(const CmplxVector& vec) const {
  52. if (d != vec.d) {
  53. cerr << "Complex vector arguments have different dimensions." << endl;
  54. exit(1);
  55. }
  56. }
  57. /***********************************************************
  58. Complex CmplxVector::operator[](int i) const {
  59. if (i < 0 || i >= d) {
  60. cerr << "Complex vector: index out of range" << endl;
  61. exit(1);
  62. }
  63. return (v[i]);
  64. }
  65. Complex& CmplxVector::operator[](int i) {
  66. if (i < 0 || i >= d) {
  67. cerr << "Complex vector: index out of range" << endl;
  68. exit(1);
  69. }
  70. return (v[i]);
  71. }
  72. CmplxVector CmplxVector::operator+(const CmplxVector& vec) const {
  73. check_dimensions(vec);
  74. int n = d;
  75. CmplxVector result(n);
  76. while (n--) result.v[n] = v[n] + vec.v[n];
  77. return (result);
  78. }
  79. CmplxVector CmplxVector::operator+(const Vector& vec) const {
  80. check_dimensions(vec);
  81. int n = d;
  82. CmplxVector result(n);
  83. while (n--) result.v[n] = v[n] + vec[n];
  84. return (result);
  85. }
  86. CmplxVector CmplxVector::operator-(const CmplxVector& vec) const {
  87. check_dimensions(vec);
  88. int n = d;
  89. CmplxVector result(n);
  90. while (n--) result.v[n] = v[n] - vec.v[n];
  91. return (result);
  92. }
  93. CmplxVector CmplxVector::operator-(const Vector& vec) const {
  94. check_dimensions(vec);
  95. int n = d;
  96. CmplxVector result(n);
  97. while (n--) result.v[n] = v[n] - vec[n];
  98. return (result);
  99. }
  100. CmplxVector& CmplxVector::operator+=(const CmplxVector& vec) {
  101. check_dimensions(vec);
  102. int n = d;
  103. while (n--) v[n] += vec.v[n];
  104. return (*this);
  105. }
  106. CmplxVector& CmplxVector::operator+=(const Vector& vec) {
  107. check_dimensions(vec);
  108. int n = d;
  109. while (n--) v[n] += vec[n];
  110. return (*this);
  111. }
  112. CmplxVector& CmplxVector::operator-=(const CmplxVector& vec) {
  113. check_dimensions(vec);
  114. int n = d;
  115. while (n--) v[n] -= vec.v[n];
  116. return (*this);
  117. }
  118. CmplxVector& CmplxVector::operator-=(const Vector& vec) {
  119. check_dimensions(vec);
  120. int n = d;
  121. while (n--) v[n] -= vec[n];
  122. return (*this);
  123. }
  124. CmplxVector& CmplxVector::operator*=(double a) {
  125. int n = d;
  126. while (n--) v[n] *= a;
  127. return (*this);
  128. }
  129. CmplxVector& CmplxVector::operator*=(Complex a) {
  130. int n = d;
  131. while (n--) v[n] *= a;
  132. return (*this);
  133. }
  134. CmplxVector CmplxVector::operator-() {
  135. int n = d;
  136. CmplxVector result(n);
  137. while (n--) result.v[n] = -v[n];
  138. return (result);
  139. }
  140. CmplxVector CmplxVector::operator*(double x) const {
  141. int n = d;
  142. CmplxVector result(n);
  143. while (n--) result.v[n] = v[n] * x;
  144. return (result);
  145. }
  146. CmplxVector CmplxVector::operator/(double x) const {
  147. int n = d;
  148. CmplxVector result(n);
  149. while (n--) result.v[n] = v[n] / x;
  150. return (result);
  151. }
  152. CmplxVector CmplxVector::operator*(const Complex& x) const {
  153. int n = d;
  154. CmplxVector result(n);
  155. while (n--) result.v[n] = v[n] * x;
  156. return (result);
  157. }
  158. CmplxVector CmplxVector::operator/(const Complex& x) const {
  159. int n = d;
  160. CmplxVector result(n);
  161. while (n--) result.v[n] = v[n] / x;
  162. return (result);
  163. }
  164. Complex CmplxVector::operator*(const CmplxVector& vec) const {
  165. check_dimensions(vec);
  166. Complex result(0.,0.);
  167. int n = d;
  168. while (n--) result += v[n] * vec.v[n];
  169. return (result);
  170. }
  171. Complex CmplxVector::operator*(const Vector& vec) const {
  172. check_dimensions(vec);
  173. Complex result(0.,0.);
  174. int n = d;
  175. while (n--) result += v[n] * vec[n];
  176. return (result);
  177. }
  178. CmplxVector& CmplxVector::operator=(const CmplxVector& vec) {
  179. if (this != &vec) {
  180. int n = vec.d;
  181. if (n != d) {
  182. delete v;
  183. v = new Complex [n];
  184. d = n;
  185. }
  186. while (n--) v[n] = vec.v[n];
  187. }
  188. return (*this);
  189. }
  190. CmplxVector& CmplxVector::operator=(const Vector& vec) {
  191. int n = vec.dim();
  192. if (n != d) {
  193. delete v;
  194. v = new Complex [n];
  195. d = n;
  196. }
  197. while (n--) v[n] = vec[n];
  198. return (*this);
  199. }
  200. CmplxVector& CmplxVector::operator=(double a) {
  201. int n = d;
  202. while (n--) v[n] = a;
  203. return (*this);
  204. }
  205. CmplxVector& CmplxVector::operator=(Complex a) {
  206. int n = d;
  207. while (n--) v[n] = a;
  208. return (*this);
  209. }
  210. int CmplxVector::operator==(const CmplxVector& vec) const {
  211. if (vec.d != d) return (0);
  212. int i = 0;
  213. while ((i < d) && (v[i]==vec.v[i])) i++;
  214. return ((i == d) ? 1 : 0);
  215. }
  216. ***********************************/
  217. CmplxVector CmplxVector::conjg() const {
  218. int n = d;
  219. CmplxVector tmp(n);
  220. for(int i = 0; i < n; i++) tmp[i] = ::conjg(v[i]);
  221. return (tmp);
  222. }
  223. Vector CmplxVector::real() const {
  224. Vector tmp(d);
  225. int n = d;
  226. while(n--) tmp[n] = ::real(v[n]);
  227. return (tmp);
  228. }
  229. Vector CmplxVector::imag() const {
  230. Vector tmp(d);
  231. int n = d;
  232. while(n--) tmp[n] = ::imag(v[n]);
  233. return (tmp);
  234. }
  235. Vector CmplxVector::abs() const {
  236. Vector result(d);
  237. for(int i = 0; i < d; i++)
  238. result[i]=::cabs(v[i]);
  239. return (result);
  240. }
  241. CmplxVector& CmplxVector::resize(int new_dim) {
  242. if(d != new_dim) {
  243. delete v;
  244. v = new Complex [new_dim];
  245. d = new_dim;
  246. }
  247. while (new_dim--) v[new_dim] = 0.0;
  248. return (*this);
  249. }
  250. /************ friends ***********************/
  251. CmplxVector operator*(double f, const CmplxVector& v)
  252. {
  253. return v*f;
  254. }
  255. CmplxVector operator*(const Complex& f, const CmplxVector& v)
  256. {
  257. return v*f;
  258. }
  259. Complex inner(const CmplxVector& x, const CmplxVector& x0)
  260. {
  261. return(x.conjg()*x0);
  262. }
  263. ostream& operator<<(ostream& out, const CmplxVector& v)
  264. {
  265. out.precision(10);
  266. for (int i = 0; i < v.d; i++)
  267. //out << v.v[i] << endl;
  268. out << " [" << i << "]=" << v.v[i] << endl;
  269. return (out);
  270. }
  271. istream& operator>>(istream& in, CmplxVector& x)
  272. {
  273. int i = 0;
  274. while (i < x.d && in >> x.v[i++]);
  275. return (in);
  276. }
  277. double norm(const CmplxVector& x)
  278. {
  279. double sum = 0.;
  280. int n = x.dim();
  281. for(int i = 0; i < n; i++) sum += cabs(x[i]);
  282. return (sum);
  283. }
  284. double norm2(const CmplxVector& x)
  285. {
  286. double sum = 0;
  287. int n = x.dim();
  288. //while (n--) sum += real(x[n])*real(x[n])+imag(x[n])*imag(x[n]);
  289. for(int i = 0; i < n; i++) sum += cabs(x[i])*cabs(x[i]);
  290. return (sqrt(sum));
  291. }
  292. double normsa(const CmplxVector& x)
  293. {
  294. Complex sum = 0.;
  295. int n = x.dim();
  296. for(int i = 0; i < n; i++) sum += x[i];
  297. return (cabs(sum));
  298. }
  299. CmplxVector operator^(const CmplxVector& v,const CmplxVector& w) {
  300. if (v.dim() != w.dim()) {
  301. cerr << "CmplxVector=: vectors have different length" << endl;
  302. exit(1);
  303. }
  304. int n = v.dim();
  305. CmplxVector result(n, 0.0);
  306. while (n--) result[n] = v[n] * w[n];
  307. return (result);
  308. }