001// --- BEGIN LICENSE BLOCK ---
002/* 
003 * Copyright (c) 2009, Mikio L. Braun
004 * All rights reserved.
005 * 
006 * Redistribution and use in source and binary forms, with or without
007 * modification, are permitted provided that the following conditions are
008 * met:
009 * 
010 *     * Redistributions of source code must retain the above copyright
011 *       notice, this list of conditions and the following disclaimer.
012 * 
013 *     * Redistributions in binary form must reproduce the above
014 *       copyright notice, this list of conditions and the following
015 *       disclaimer in the documentation and/or other materials provided
016 *       with the distribution.
017 * 
018 *     * Neither the name of the Technische Universit?t Berlin nor the
019 *       names of its contributors may be used to endorse or promote
020 *       products derived from this software without specific prior
021 *       written permission.
022 * 
023 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
024 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
025 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
026 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
027 * HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
028 * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
029 * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
030 * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
031 * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
032 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
033 * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
034 */
035// --- END LICENSE BLOCK ---
036
037package org.jblas;
038
039import org.jblas.exceptions.SizeException;
040
041import java.io.DataInputStream;
042import java.io.DataOutputStream;
043import java.io.FileInputStream;
044import java.io.FileOutputStream;
045import java.io.IOException;
046
047public class ComplexDoubleMatrix {
048        
049        public int rows;
050        public int columns;
051        public int length;
052        public double[] data = null; // rows are contiguous
053
054        /**************************************************************************
055         * 
056         * Constructors and factory functions
057         * 
058         **************************************************************************/
059
060        /**
061   * Create a new matrix with <i>newRows</i> rows, <i>newColumns</i> columns
062         * using <i>newData></i> as the data.
063         */
064        public ComplexDoubleMatrix(int newRows, int newColumns, double... newData) {
065                rows = newRows;
066                columns = newColumns;
067                length = rows * columns;
068
069    if (newData.length != 2 * newRows * newColumns)
070                        throw new IllegalArgumentException(
071                                        "Passed data must match matrix dimensions.");
072                data = newData;
073        }
074        
075        /**
076         * Creates a new <i>n</i> times <i>m</i> <tt>ComplexDoubleMatrix</tt>.
077         * @param newRows the number of rows (<i>n</i>) of the new matrix.
078         * @param newColumns the number of columns (<i>m</i>) of the new matrix.
079         */
080        public ComplexDoubleMatrix(int newRows, int newColumns) {
081                this(newRows, newColumns, new double[2 * newRows * newColumns]);
082        }
083        
084        /**
085         * Creates a new <tt>ComplexDoubleMatrix</tt> of size 0 times 0.
086         */
087        public ComplexDoubleMatrix() {
088                this(0, 0, null);
089        }
090
091        /**
092         * Create a Matrix of length <tt>len</tt>. By default, this creates a row vector.
093         * @param len
094         */
095        public ComplexDoubleMatrix(int len) {
096                this(len, 1, new double[2 * len]);
097        }
098        
099        public ComplexDoubleMatrix(double[] newData) {
100                this(newData.length/2);
101                                
102                data = newData;
103        }
104
105        public ComplexDoubleMatrix(ComplexDouble[] newData) {
106                this(newData.length);
107                                
108                for (int i = 0; i < newData.length; i++)
109                        put(i, newData[i]);
110        }
111                
112        
113  /** Construct a complex matrix from a real matrix. */
114  public ComplexDoubleMatrix(DoubleMatrix m) {
115    this(m.rows, m.columns);
116
117    NativeBlas.dcopy(m.length, m.data, 0, 1, data, 0, 2);
118  }
119
120  /** Construct a complex matrix from separate real and imaginary parts. Either
121   * part can be set to null in which case it will be ignored.
122   */
123  public ComplexDoubleMatrix(DoubleMatrix real, DoubleMatrix imag) {
124      this(real.rows, real.columns);
125      real.assertSameSize(imag);
126
127      if (real != null)
128          NativeBlas.dcopy(length, real.data, 0, 1, data, 0, 2);
129      if (imag != null)
130          NativeBlas.dcopy(length, imag.data, 0, 1, data, 1, 2);
131  }
132        
133        /**
134         * Creates a new matrix by reading it from a file.
135         * @param filename the path and name of the file to read the matrix from
136         * @throws IOException 
137         */
138        public ComplexDoubleMatrix(String filename) throws IOException {
139                load(filename);
140        }
141        
142        /**
143         * Creates a new <i>n</i> times <i>m</i> <tt>ComplexDoubleMatrix</tt> from
144         * the given <i>n</i> times <i>m</i> 2D data array. The first dimension of the array makes the
145         * rows (<i>n</i>) and the second dimension the columns (<i>m</i>). For example, the
146         * given code <br/><br/>
147         * <code>new ComplexDoubleMatrix(new double[][]{{1d, 2d, 3d}, {4d, 5d, 6d}, {7d, 8d, 9d}}).print();</code><br/><br/>
148         * will constructs the following matrix:
149         * <pre>
150         * 1.0  2.0     3.0
151         * 4.0  5.0     6.0
152         * 7.0  8.0     9.0
153         * </pre>.
154         * @param data <i>n</i> times <i>m</i> data array
155         */ 
156        public ComplexDoubleMatrix(double[][] data) {
157                this(data.length, data[0].length);
158                                                
159                for (int r = 0; r < rows; r++)
160                        assert(data[r].length == columns);
161                
162                for (int r = 0; r < rows; r++)
163                        for (int c = 0; c < columns; c++)
164                                put(r, c, data[r][c]);
165        }
166        
167        /**
168         * Creates a new matrix in which all values are equal 0.
169         * @param rows number of rows
170         * @param columns number of columns
171         * @return new matrix
172         */
173        public static ComplexDoubleMatrix zeros(int rows, int columns) {
174                return new ComplexDoubleMatrix(rows, columns);
175        }
176        
177        public static ComplexDoubleMatrix zeros(int length) {
178                return zeros(length, 1);
179        }
180
181        /**
182         * Creates a new matrix in which all values are equal 1.
183         * @param rows number of rows
184         * @param columns number of columns
185         * @return new matrix
186         */
187        public static ComplexDoubleMatrix ones(int rows, int columns) {
188                ComplexDoubleMatrix m = new ComplexDoubleMatrix(rows, columns);
189                
190                for (int i = 0; i < rows * columns; i++)
191                        m.put(i, 1.0);
192                
193                return m;
194        }
195        
196        public static ComplexDoubleMatrix ones(int length) {
197                return ones(length, 1);
198        }
199        
200        /**
201         * Creates a new matrix where the values of the given vector are the diagonal values of
202         * the matrix.
203         * @param x the diagonal values
204         * @return new matrix
205         */
206        public static ComplexDoubleMatrix diag(ComplexDoubleMatrix x) {
207                ComplexDoubleMatrix m = new ComplexDoubleMatrix(x.length, x.length);
208                
209                for (int i = 0; i < x.length; i++)
210                        m.put(i, i, x.get(i));
211                
212                return m;
213        }
214
215  /**
216   * Construct a matrix of arbitrary shape and set the diagonal according
217   * to a passed vector.
218   *
219   * length of needs to be smaller than rows or columns.
220   *
221   * @param x vector to fill the diagonal with
222   * @param rows number of rows of the resulting matrix
223   * @param columns number of columns of the resulting matrix
224   * @return a matrix with dimensions rows * columns whose diagonal elements are filled by x
225   */
226  public static ComplexDoubleMatrix diag(ComplexDoubleMatrix x, int rows, int columns) {
227    if (x.length > rows || x.length > columns) {
228      throw new SizeException("Length of diagonal matrix must be larger than both rows and columns.");
229    }
230    
231    ComplexDoubleMatrix m = new ComplexDoubleMatrix(rows, columns);
232
233    for (int i = 0; i < x.length; i++)
234      m.put(i, i, x.get(i));
235
236    return m;
237  }
238        
239        /**
240         * Create a 1 * 1 - matrix. For many operations, this matrix functions like a
241         * normal double
242         * @param s value of the matrix
243         * @return the constructed ComplexDoubleMatrix 
244         */
245        public static ComplexDoubleMatrix scalar(double s) {
246                ComplexDoubleMatrix m = new ComplexDoubleMatrix(1, 1);
247                m.put(0, 0, s);
248                return m;
249        }
250        
251        /** Test whether a matrix is scalar */
252        public boolean isScalar() {
253                return length == 1;
254        }
255        
256        /** Return the first element of the matrix */
257        public ComplexDouble scalar() {
258                return get(0);
259        }
260        
261        public static ComplexDoubleMatrix concatHorizontally(ComplexDoubleMatrix A, ComplexDoubleMatrix B) {
262                if (A.rows != B.rows)
263                        throw new SizeException("Matrices don't have same number of rows.");
264                
265                ComplexDoubleMatrix result = new ComplexDoubleMatrix(A.rows, A.columns + B.columns);
266                SimpleBlas.copy(A, result);
267                NativeBlas.zcopy(B.length, B.data, 0, 1, result.data, A.length, 1);
268                return result;
269        }
270
271        public static ComplexDoubleMatrix concatVertically(ComplexDoubleMatrix A, ComplexDoubleMatrix B) {
272                if (A.columns != B.columns)
273                        throw new SizeException("Matrices don't have same number of columns.");
274                
275                ComplexDoubleMatrix result = new ComplexDoubleMatrix(A.rows + B.rows, A.columns);
276
277                for (int i = 0; i < A.columns; i++) {
278                        NativeBlas.zcopy(A.rows, A.data, A.index(0, i), 1, result.data, result.index(0, i), 1);
279                        NativeBlas.zcopy(B.rows, B.data, B.index(0, i), 1, result.data, result.index(A.rows, i), 1);
280                }
281                
282                return result;
283        }
284        
285        /**************************************************************************
286         * Working with slices (Man! 30+ methods just to make this a bit flexible...) 
287         */
288
289        public ComplexDoubleMatrix get(int[] indices) {
290                ComplexDoubleMatrix result = new ComplexDoubleMatrix(indices.length);
291                
292                for (int i = 0; i < indices.length; i++)
293                        result.put(i, get(indices[i]));
294                
295                return result;
296        }
297        
298        public ComplexDoubleMatrix get(int r, int[] indices) {
299                ComplexDoubleMatrix result = new ComplexDoubleMatrix(1, indices.length);
300                
301                for (int i = 0; i < indices.length; i++)
302                        result.put(i, get(r, indices[i]));
303                
304                return result;
305        }
306        
307        public ComplexDoubleMatrix get(int[] indices, int c) {
308                ComplexDoubleMatrix result = new ComplexDoubleMatrix(indices.length, 1);
309                
310                for (int i = 0; i < indices.length; i++)
311                        result.put(i, get(indices[i], c));
312                
313                return result;
314        }
315        
316        public ComplexDoubleMatrix get(int[] rindices, int[] cindices) {
317                ComplexDoubleMatrix result = new ComplexDoubleMatrix(rindices.length, cindices.length);
318                
319                for (int i = 0; i < rindices.length; i++)
320                        for (int j = 0; j < cindices.length; j++)
321                                result.put(i, j, get(rindices[i], cindices[j]));
322                
323                return result;
324        }
325        
326        public ComplexDoubleMatrix get(ComplexDoubleMatrix indices) {
327                return get(indices.findIndices());
328        }
329
330        public ComplexDoubleMatrix get(int r, ComplexDoubleMatrix indices) {
331                return get(r, indices.findIndices());
332        }
333        
334        public ComplexDoubleMatrix get(ComplexDoubleMatrix indices, int c) {
335                return get(indices.findIndices(), c);
336        }
337
338        public ComplexDoubleMatrix get(ComplexDoubleMatrix rindices, ComplexDoubleMatrix cindices) {
339                return get(rindices.findIndices(), cindices.findIndices());
340        }
341        
342        private void checkLength(int l) {
343                if (length != l)
344                        throw new SizeException("Matrix does not have the necessary length (" + length + " != " + l + ").");
345        }
346
347        private void checkRows(int r) {
348                if (rows != r)
349                        throw new SizeException("Matrix does not have the necessary length (" + length + " != " + r + ").");
350        }
351        
352        private void checkColumns(int c) {
353                if (columns != c)
354                        throw new SizeException("Matrix does not have the necessary length (" + length + " != " + c + ").");
355        }
356
357        public ComplexDoubleMatrix put(int[] indices, ComplexDoubleMatrix x) {
358                if (x.isScalar())
359                        return put(indices, x.scalar());
360                x.checkLength(indices.length);
361                
362                for (int i = 0; i < indices.length; i++)
363                        put(indices[i], x.get(i));
364                
365                return this;
366        }
367        
368        public ComplexDoubleMatrix put(int r, int[] indices, ComplexDoubleMatrix x) {
369                if (x.isScalar())
370                        return put(r, indices, x.scalar());
371                x.checkColumns(indices.length);
372                
373                for (int i = 0; i < indices.length; i++)
374                        put(r, indices[i], x.get(i));
375                
376                return this;
377        }
378        
379        public ComplexDoubleMatrix put(int[] indices, int c, ComplexDoubleMatrix x) {
380                if (x.isScalar())
381                        return put(indices, c, x.scalar());             
382                x.checkRows(indices.length);
383                
384                for (int i = 0; i < indices.length; i++)
385                        put(indices[i], c, x.get(i));
386                
387                return this;
388        }
389        
390        public ComplexDoubleMatrix put(int[] rindices, int[] cindices, ComplexDoubleMatrix x) {
391                if (x.isScalar())
392                        return put(rindices, cindices, x.scalar());             
393                x.checkRows(rindices.length);
394                x.checkColumns(cindices.length);
395                
396                for (int i = 0; i < rindices.length; i++)
397                        for (int j = 0; j < cindices.length; j++)
398                                put(rindices[i], cindices[j], x.get(i,j));
399                
400                return this;
401        }
402        
403        public ComplexDoubleMatrix put(int[] indices, double v) {
404                for (int i = 0; i < indices.length; i++)
405                        put(indices[i], v);
406                
407                return this;
408        }
409
410        public ComplexDoubleMatrix putReal(int[] indices, double v) {
411                return put(indices, v);
412        }
413
414        public ComplexDoubleMatrix putImag(int[] indices, double v) {
415                for (int i = 0; i < indices.length; i++)
416                        putImag(indices[i], v);
417                
418                return this;
419        }
420
421        public ComplexDoubleMatrix put(int[] indices, ComplexDouble v) {
422                for (int i = 0; i < indices.length; i++)
423                        put(indices[i], v);
424                
425                return this;
426        }
427
428        public ComplexDoubleMatrix put(int r, int[] indices, double v) {
429                for (int i = 0; i < indices.length; i++)
430                        put(r, indices[i], v);
431                
432                return this;
433        }
434
435        public ComplexDoubleMatrix putReal(int r, int[] indices, double v) {
436                return put(r, indices, v);
437        }
438
439        public ComplexDoubleMatrix putImag(int r, int[] indices, double v) {
440                for (int i = 0; i < indices.length; i++)
441                        putImag(r, indices[i], v);
442                
443                return this;
444        }
445
446        public ComplexDoubleMatrix put(int r, int[] indices, ComplexDouble v) {
447                for (int i = 0; i < indices.length; i++)
448                        put(r, indices[i], v);
449                
450                return this;
451        }
452
453        public ComplexDoubleMatrix put(int[] indices, int c, double v) {
454                for (int i = 0; i < indices.length; i++)
455                        put(indices[i], c, v);
456                
457                return this;
458        }
459        
460        public ComplexDoubleMatrix putReal(int[] indices, int c, double v) {
461                return put(indices, c, v);
462        }
463        
464        public ComplexDoubleMatrix putImag(int[] indices, int c, double v) {
465                for (int i = 0; i < indices.length; i++)
466                        putImag(indices[i], c, v);
467                
468                return this;
469        }
470        
471        public ComplexDoubleMatrix put(int[] indices, int c, ComplexDouble v) {
472                for (int i = 0; i < indices.length; i++)
473                        put(indices[i], c, v);
474                
475                return this;
476        }
477        
478        public ComplexDoubleMatrix put(int[] rindices, int[] cindices, double v) {
479                for (int i = 0; i < rindices.length; i++)
480                        for (int j = 0; j < cindices.length; j++)
481                                put(rindices[i], cindices[j], v);
482                
483                return this;
484        }
485        
486        public ComplexDoubleMatrix putReal(int[] rindices, int[] cindices, double v) {
487                return put(rindices, cindices, v);
488        }
489        
490        public ComplexDoubleMatrix putImag(int[] rindices, int[] cindices, double v) {
491                for (int i = 0; i < rindices.length; i++)
492                        for (int j = 0; j < cindices.length; j++)
493                                put(rindices[i], cindices[j], v);
494                
495                return this;
496        }
497
498        public ComplexDoubleMatrix put(int[] rindices, int[] cindices, ComplexDouble v) {
499                for (int i = 0; i < rindices.length; i++)
500                        for (int j = 0; j < cindices.length; j++)
501                                put(rindices[i], cindices[j], v);
502                
503                return this;
504        }
505
506        public ComplexDoubleMatrix put(ComplexDoubleMatrix indices, ComplexDoubleMatrix v) {
507                return put(indices.findIndices(), v);
508        }
509
510        public ComplexDoubleMatrix put(int r, ComplexDoubleMatrix indices, ComplexDoubleMatrix v) {
511                return put(r, indices.findIndices(), v);
512        }
513        
514        public ComplexDoubleMatrix put(ComplexDoubleMatrix indices, int c, ComplexDoubleMatrix v) {
515                return put(indices.findIndices(), c, v);
516        }
517
518        public ComplexDoubleMatrix put(ComplexDoubleMatrix rindices, ComplexDoubleMatrix cindices, ComplexDoubleMatrix v) {
519                return put(rindices.findIndices(), cindices.findIndices(), v);
520        }
521
522        public ComplexDoubleMatrix put(ComplexDoubleMatrix indices, double v) {
523                return put(indices.findIndices(), v);
524        }
525
526        public ComplexDoubleMatrix putReal(ComplexDoubleMatrix indices, double v) {
527                return put(indices, v);
528        }
529
530        public ComplexDoubleMatrix putImag(ComplexDoubleMatrix indices, double v) {
531                return putImag(indices.findIndices(), v);
532        }
533
534        public ComplexDoubleMatrix put(ComplexDoubleMatrix indices, ComplexDouble v) {
535                return put(indices.findIndices(), v);
536        }
537        
538        public ComplexDoubleMatrix put(int r, ComplexDoubleMatrix indices, double v) {
539                return put(r, indices.findIndices(), v);
540        }
541        
542        public ComplexDoubleMatrix putReal(int r, ComplexDoubleMatrix indices, double v) {
543                return put(r, indices, v);
544        }
545
546        public ComplexDoubleMatrix putImag(int r, ComplexDoubleMatrix indices, double v) {
547                return putImag(r, indices.findIndices(), v);
548        }
549
550        public ComplexDoubleMatrix put(int r, ComplexDoubleMatrix indices, ComplexDouble v) {
551                return put(r, indices.findIndices(), v);
552        }
553
554        public ComplexDoubleMatrix put(ComplexDoubleMatrix indices, int c, double v) {
555                return put(indices.findIndices(), c, v);
556        }
557
558        public ComplexDoubleMatrix putReal(ComplexDoubleMatrix indices, int c, double v) {
559                return put(indices, c, v);
560        }
561
562        public ComplexDoubleMatrix putImag(ComplexDoubleMatrix indices, int c, double v) {
563                return putImag(indices.findIndices(), c, v);
564        }
565
566        public ComplexDoubleMatrix put(ComplexDoubleMatrix indices, int c, ComplexDouble v) {
567                return put(indices.findIndices(), c, v);
568        }
569
570        public ComplexDoubleMatrix put(ComplexDoubleMatrix rindices, ComplexDoubleMatrix cindices, double v) {
571                return put(rindices.findIndices(), cindices.findIndices(), v);
572        }
573
574        public ComplexDoubleMatrix putReal(ComplexDoubleMatrix rindices, ComplexDoubleMatrix cindices, double v) {
575                return putReal(rindices, cindices, v);
576        }
577
578        public ComplexDoubleMatrix putImag(ComplexDoubleMatrix rindices, ComplexDoubleMatrix cindices, double v) {
579                return putImag(rindices.findIndices(), cindices.findIndices(), v);
580        }
581
582        public ComplexDoubleMatrix put(ComplexDoubleMatrix rindices, ComplexDoubleMatrix cindices, ComplexDouble v) {
583                return put(rindices.findIndices(), cindices.findIndices(), v);
584        }
585
586        
587        public int[] findIndices() {
588                int len = 0;
589                for (int i = 0; i < length; i++)
590                        if (!get(i).isZero())
591                                len++;
592                
593                int[] indices = new int[len];
594                int c = 0;
595                
596                for (int i = 0; i < length; i++)
597                        if (!get(i).isZero())
598                                indices[c++] = i;
599                
600                return indices;
601        }
602        
603        /**************************************************************************
604         * Basic operations (copying, resizing, element access)
605         */
606        
607        /** Return transposed copy of this matrix */
608        public ComplexDoubleMatrix transpose() {
609                ComplexDoubleMatrix result = new ComplexDoubleMatrix(columns, rows);
610
611                ComplexDouble c = new ComplexDouble(0);
612
613                for (int i = 0; i < rows; i++)
614                        for (int j = 0; j < columns; j++)
615                                result.put(j, i, get(i, j, c));
616                
617                return result;
618        }
619
620        public ComplexDoubleMatrix hermitian() {
621            ComplexDoubleMatrix result = new ComplexDoubleMatrix(columns, rows);
622
623            ComplexDouble c = new ComplexDouble(0);
624
625            for (int i = 0; i < rows; i++)
626                for (int j = 0; j < columns; j++)
627                    result.put(j, i, get(i, j, c).conji());
628            return result;
629        }
630
631        /**
632         * Compute complex conjugate (in-place).
633         */
634        public ComplexDoubleMatrix conji() {
635            ComplexDouble c = new ComplexDouble(0.0);
636            for (int i = 0; i < length; i++)
637                put(i, get(i, c).conji());
638            return this;
639        }
640
641        /**
642         * Compute complex conjugate.
643         */
644        public ComplexDoubleMatrix conj() {
645            return dup().conji();
646        }
647
648                
649        /** Compare two matrices.
650         * @param o Object to compare to
651         * @return true if and only if other is also a ComplexDoubleMatrix which has the same size and the
652         * maximal absolute difference in matrix elements is smaller thatn 1e-6.  */
653        public boolean equals(Object o) {
654                if (!(o instanceof ComplexDoubleMatrix))
655                        return false;
656
657                ComplexDoubleMatrix other = (ComplexDoubleMatrix) o;
658
659                if (!sameSize(other))
660                        return false;
661                
662                DoubleMatrix diff = MatrixFunctions.absi(sub(other)).getReal();
663                
664                return diff.max() / (rows * columns) < 1e-6;
665        }
666
667        
668        /** Resize the matrix. All elements will be set to zero. */
669        public void resize(int newRows, int newColumns) {
670                rows = newRows;
671                columns = newColumns;
672                length = newRows * newColumns;
673                data = new double[2 * rows * columns];
674        }
675
676        
677        /** Reshape the matrix. Number of elements must not change. */
678        public ComplexDoubleMatrix reshape(int newRows, int newColumns) {
679                if (length != newRows * newColumns)
680                        throw new IllegalArgumentException(
681                                        "Number of elements must not change.");
682
683                rows = newRows;
684                columns = newColumns;
685                
686                return this;
687        }
688
689        /** Checks whether two matrices have the same size. */
690        public boolean sameSize(ComplexDoubleMatrix a) {
691                return rows == a.rows && columns == a.columns;
692        }
693
694        /** 
695         * Assert that two matrices have the same size.
696         * 
697         * @param a the other matrix
698         * @throws SizeException if matrix sizes don't match. 
699         * */
700        public void assertSameSize(ComplexDoubleMatrix a) {
701                if (!sameSize(a))
702                        throw new SizeException("Matrices must have the same size.");
703        }
704        
705        /** 
706         * Check whether this can be multiplied with a. 
707         * 
708         * @param a right-hand-side of the multiplication.
709         * @return true iff <tt>this.columns == a.rows</tt>
710         */
711        public boolean multipliesWith(ComplexDoubleMatrix a) {
712                return columns == a.rows;
713        }
714        
715        public void assertMultipliesWith(ComplexDoubleMatrix a) {
716                if (!multipliesWith(a))
717                        throw new SizeException("Number of columns of left matrix must be equal to number of rows of right matrix.");
718        }
719        
720        public boolean sameLength(ComplexDoubleMatrix a) {
721                return length == a.length;
722        }
723        
724        public void assertSameLength(ComplexDoubleMatrix a) {
725                if (!sameLength(a))
726                        throw new SizeException("Matrices must have same length (is: " + length + " and " + a.length + ")");
727        }
728        
729        /** Copy ComplexDoubleMatrix a to this. this a is resized if necessary. */
730        public ComplexDoubleMatrix copy(ComplexDoubleMatrix a) {
731                if (!sameSize(a))
732                        resize(a.rows, a.columns);
733                
734                SimpleBlas.copy(a, this);
735                return a;
736        }
737        
738        /** Returns a duplicate of this matrix. Geometry is the same (including offsets, transpose, etc.),
739         * but the buffer is not shared.
740         */
741        public ComplexDoubleMatrix dup() {
742                ComplexDoubleMatrix out = new ComplexDoubleMatrix(rows, columns);
743
744                JavaBlas.rcopy(2*length, data, 0, 1, out.data, 0, 1);
745                
746                return out;
747        }
748        
749        public ComplexDoubleMatrix swapColumns(int i, int j) {
750                NativeBlas.zswap(rows, data, index(0, i), 1, data, index(0, j), 1);
751                return this;
752        }
753        
754        public ComplexDoubleMatrix swapRows(int i, int j) {
755                NativeBlas.zswap(columns, data, index(i, 0), rows, data, index(j, 0), rows);
756                return this;
757        }
758                
759        /** Set matrix element */
760        public ComplexDoubleMatrix put(int rowIndex, int columnIndex, double value) {
761                data[2*index(rowIndex, columnIndex)] =  value;
762                return this;
763        }
764
765        public ComplexDoubleMatrix put(int rowIndex, int columnIndex, double realValue, double complexValue) {
766                data[2*index(rowIndex, columnIndex)] =  realValue;
767                data[2*index(rowIndex, columnIndex)+1] =  complexValue;
768                return this;
769        }
770
771        public ComplexDoubleMatrix put(int rowIndex, int columnIndex, ComplexDouble value) {
772                int i = 2*index(rowIndex, columnIndex);
773                data[i] = value.real(); data[i+1] = value.imag();
774                return this;
775        }
776
777        public ComplexDoubleMatrix putReal(int rowIndex, int columnIndex, double value) {
778                data[2*index(rowIndex, columnIndex)] = value;
779                return this;
780        }
781
782        public ComplexDoubleMatrix putImag(int rowIndex, int columnIndex, double value) {
783                data[2*index(rowIndex, columnIndex)+1] = value;
784                return this;
785        }
786        
787        /** Retrieve matrix element */
788        public ComplexDouble get(int rowIndex, int columnIndex) {
789            int i = 2*index(rowIndex, columnIndex);
790            return new ComplexDouble(data[i], data[i+1]);
791        }
792
793        /** Get matrix element, passing the variable to store the result. */
794        public ComplexDouble get(int rowIndex, int columnIndex, ComplexDouble result) {
795            return get(index(rowIndex, columnIndex), result);
796        }
797        
798        public DoubleMatrix getReal() {
799                DoubleMatrix result = new DoubleMatrix(rows, columns);
800                
801                NativeBlas.dcopy(length, data, 0, 2, result.data, 0, 1);
802                
803                return result;
804        }
805
806        /** Get index of an element */
807        public int index(int rowIndex, int columnIndex) {
808                //System.out.printf("Index for (%d, %d) -> %d\n", rowIndex, columnIndex, (rows * columnIndex + rowIndex) * 2);
809                return rows * columnIndex + rowIndex;
810        }
811
812        public ComplexDouble get(int i) {
813                return new ComplexDouble(data[i * 2], data[i * 2 + 1]);
814        }
815        
816        public ComplexDouble get(int i, ComplexDouble result) {
817            return result.set(data[i * 2], data[i*2+1]);
818        }
819        
820        public double getReal(int i) {
821                return data[2*i];
822        }
823        
824        public double getImag(int i) {
825                return data[2*i + 1]; 
826        }
827
828        public ComplexDoubleMatrix put(int i, double v) {
829                data[2*i] = v;
830                return this;
831        }
832
833        public ComplexDoubleMatrix put(int i, double r, double c) {
834            data[2*i] = r;
835            data[2*i+1] = c;
836            return this;
837        }
838        
839        public ComplexDoubleMatrix put(int i, ComplexDouble v) {
840                data[2*i] = v.real();
841                data[2*i+1] = v.imag();
842                return this;
843        }
844        
845        public ComplexDoubleMatrix putReal(int i, double v) {
846                return put(i, v);
847        }
848        
849        public ComplexDoubleMatrix putImag(int i, double v) {
850                data[2*i+1] = v;
851                return this;
852        }
853
854        public int getRows() {
855                return rows;
856        }
857        
858        public int getColumns() {
859                return columns;
860        }
861        
862        public int getLength() {
863                return length;
864        }
865        
866        /** Checks whether the matrix is empty. */
867        public boolean isEmpty() {
868                return columns == 0 || rows == 0;
869        }
870        
871        /** Checks whether the matrix is square. */
872        public boolean isSquare() {
873                return columns == rows;
874        }
875        
876        public void assertSquare() {
877                if (!isSquare())
878                        throw new SizeException("Matrix must be square!");
879        }
880        
881        /** Checks whether the matrix is a vector. */
882        public boolean isVector() {
883                return columns == 1 || rows == 1;
884        }
885        
886        public boolean isRowVector() {
887                return columns == 1;
888        }
889        
890        public boolean isColumnVector() {
891                return rows == 1;
892        }
893                
894        /** Get diagonal of the matrix. */
895        public ComplexDoubleMatrix diag() {
896                ComplexDoubleMatrix d = new ComplexDoubleMatrix(rows);
897                NativeBlas.zcopy(rows, data, 0, rows + 1, d.data, 0, 1);
898                return d;
899        }
900        
901        /** Get real part of the matrix. */
902        public DoubleMatrix real() {
903            DoubleMatrix result = new DoubleMatrix(rows, columns);
904            NativeBlas.dcopy(length, data, 0, 2, result.data, 0, 1);
905            return result;
906        }
907        
908        /** Get imaginary part of the matrix. */
909        public DoubleMatrix imag() {
910            DoubleMatrix result = new DoubleMatrix(rows, columns);
911            NativeBlas.dcopy(length, data, 1, 2, result.data, 0, 1);
912            return result;            
913        }
914
915        
916        /** 
917         * Pretty-print this matrix to <tt>System.out</tt>. 
918         * */
919        public void print() {
920                System.out.println(toString());
921        }
922
923        /** 
924         * Generate string representation of this matrix 
925         * (multi-line).
926         * */
927        public String toString() {
928                StringBuilder s = new StringBuilder();
929
930                s.append("[");
931                
932                for (int i = 0; i < rows; i++) {
933                        for (int j = 0; j < columns; j++) {
934                                s.append(get(i, j));
935                                if (j < columns - 1)
936                                        s.append(", ");
937                        }
938                        if (i < rows - 1)
939                                s.append("; ");
940                }
941
942                s.append("]");
943                
944                return s.toString();
945        }
946
947        public double[] toDoubleArray() {
948                double[] array = new double[2*length];
949                
950                for (int i = 0; i < 2*length; i++)
951                        array[i] = data[i];
952                
953                return array;
954        }
955        
956        public ComplexDouble[] toArray() {
957                ComplexDouble[] array = new ComplexDouble[length];
958                
959                for (int i = 0; i < length; i++)
960                        array[i] = get(i);
961                
962                return array;           
963        }
964        
965        public ComplexDouble[][] toArray2() {
966                ComplexDouble[][] array = new ComplexDouble[rows][columns];
967                
968                for (int r = 0; r < rows; r++)
969                        for (int c = 0; c < columns; c++)
970                                array[r][c] = get(r, c);
971                                
972                return array;
973        }
974        
975        public boolean[] toBooleanArray() {
976                boolean[] array = new boolean[length];
977                
978                for (int i = 0; i < length; i++)
979                        array[i] = get(i).isZero() ? false : true;
980                
981                return array;
982        }
983        
984        public boolean[][] toBooleanArray2() {
985                boolean[][] array = new boolean[rows][columns];
986                
987                for (int r = 0; r < rows; r++)
988                        for (int c = 0; c < columns; c++)
989                                array[r][c] = get(r, c).isZero() ? false : true;
990                                
991                return array;
992        }
993
994        /**************************************************************************
995         * Arithmetic Operations
996         */
997
998        /** 
999         * Ensures that the result vector has the same length as this. If not,
1000         * resizing result is tried, which fails if result == this or result == other.
1001         */
1002        private void ensureResultLength(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
1003                if (!sameLength(result)) {
1004                        if (result == this || result == other)
1005                                throw new SizeException("Cannot resize result matrix because it is used in-place.");
1006                        result.resize(rows, columns);
1007                }
1008        }
1009
1010        /** Add two matrices. */
1011        public ComplexDoubleMatrix addi(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
1012                if (other.isScalar())
1013                        return addi(other.scalar(), result);
1014                
1015                assertSameLength(other);
1016                ensureResultLength(other, result);
1017                
1018                if (result == this)
1019                        SimpleBlas.axpy(ComplexDouble.UNIT, other, result);
1020                else if (result == other)
1021                        SimpleBlas.axpy(ComplexDouble.UNIT, this, result);
1022                else {
1023                        SimpleBlas.copy(this, result);
1024                        SimpleBlas.axpy(ComplexDouble.UNIT, other, result);
1025                }
1026
1027                return result;
1028        }
1029        
1030        /** Add a scalar to a matrix. */
1031        public ComplexDoubleMatrix addi(ComplexDouble v, ComplexDoubleMatrix result) {
1032                ensureResultLength(null, result);
1033                
1034                for (int i = 0; i < length; i++)
1035                        result.put(i, get(i).add(v));
1036                return result;
1037        }
1038        
1039        public ComplexDoubleMatrix addi(double v, ComplexDoubleMatrix result) {
1040                return addi(new ComplexDouble(v), result);
1041        }
1042
1043        /** Subtract two matrices. */
1044        public ComplexDoubleMatrix subi(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
1045                if (other.isScalar())
1046                        return subi(other.scalar(), result);
1047                
1048                assertSameLength(other);
1049                ensureResultLength(other, result);
1050                
1051                if (result == this)
1052                        SimpleBlas.axpy(ComplexDouble.NEG_UNIT, other, result);
1053                else if (result == other) {
1054                        SimpleBlas.scal(ComplexDouble.NEG_UNIT, result);
1055                        SimpleBlas.axpy(ComplexDouble.UNIT, this, result);
1056                }
1057                else {
1058                        SimpleBlas.copy(this, result);
1059                        SimpleBlas.axpy(ComplexDouble.NEG_UNIT, other, result);
1060                }
1061                return result;
1062        }
1063        
1064        /** Subtract a scalar from a matrix */
1065        public ComplexDoubleMatrix subi(ComplexDouble v, ComplexDoubleMatrix result) {
1066                ensureResultLength(null, result);
1067                
1068                for (int i = 0; i < length; i++)
1069                        result.put(i, get(i).sub(v));
1070                return result;
1071        }
1072        
1073        public ComplexDoubleMatrix subi(double v, ComplexDoubleMatrix result) {
1074                return subi(new ComplexDouble(v), result);
1075        }
1076
1077        /** 
1078         * Subtract two matrices, but subtract first from second matrix, that is, 
1079         * compute <em>result = other - this</em>. 
1080         * */
1081        public ComplexDoubleMatrix rsubi(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
1082                return other.subi(this, result);
1083        }
1084        
1085        /** Subtract a matrix from a scalar */
1086        public ComplexDoubleMatrix rsubi(ComplexDouble a, ComplexDoubleMatrix result) {
1087                ensureResultLength(null, result);
1088                
1089                for (int i = 0; i < length; i++)
1090                        result.put(i, a.sub(get(i)));
1091                return result;
1092        }
1093
1094        public ComplexDoubleMatrix rsubi(double a, ComplexDoubleMatrix result) {
1095                return rsubi(new ComplexDouble(a), result);
1096        }
1097
1098        /** (Elementwise) Multiplication */ 
1099        public ComplexDoubleMatrix muli(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
1100                if (other.isScalar())
1101                        return muli(other.scalar(), result);
1102                
1103                assertSameLength(other);
1104                ensureResultLength(other, result);
1105                
1106                ComplexDouble c = new ComplexDouble(0.0);
1107                ComplexDouble d = new ComplexDouble(0.0);
1108                
1109                for (int i = 0; i < length; i++)
1110                        result.put(i, get(i, c).muli(other.get(i, d)));
1111                return result;
1112        }
1113        
1114        /** (Elementwise) Multiplication with a scalar */
1115        public ComplexDoubleMatrix muli(ComplexDouble v, ComplexDoubleMatrix result) {
1116                ensureResultLength(null, result);
1117                
1118                ComplexDouble c = new ComplexDouble(0.0);
1119                
1120                for (int i = 0; i < length; i++)
1121                        result.put(i, get(i, c).muli(v));
1122                return result;
1123        }
1124
1125        public ComplexDoubleMatrix muli(double v, ComplexDoubleMatrix result) {
1126                return muli(new ComplexDouble(v), result);
1127        }
1128
1129        /** Matrix-Matrix Multiplication */
1130        public ComplexDoubleMatrix mmuli(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
1131                if (other.isScalar())
1132                        return muli(other.scalar(), result);
1133
1134                /* check sizes and resize if necessary */
1135                assertMultipliesWith(other);
1136                if (result.rows != rows || result.columns != other.columns) {
1137                        if (result != this && result != other)
1138                                result.resize(rows, other.columns);
1139                        else
1140                                throw new SizeException("Cannot resize result matrix because it is used in-place.");
1141                }
1142                
1143                if (result == this || result == other) {
1144                        /* actually, blas cannot do multiplications in-place. Therefore, we will fake by
1145                         * allocating a temporary object on the side and copy the result later.
1146                         */
1147                        ComplexDoubleMatrix temp = new ComplexDoubleMatrix(result.rows, result.columns);
1148                        SimpleBlas.gemm(ComplexDouble.UNIT, this, other, ComplexDouble.ZERO, temp);
1149                        SimpleBlas.copy(temp, result);
1150                }
1151                else {
1152                        SimpleBlas.gemm(ComplexDouble.UNIT, this, other, ComplexDouble.ZERO, result);
1153                }               
1154                return result;
1155        }
1156        
1157        /** Matrix-Matrix Multiplication with a scalar (for symmetry, does the
1158         * same as muli(scalar)
1159         */
1160        public ComplexDoubleMatrix mmuli(ComplexDouble v, ComplexDoubleMatrix result) {
1161                return muli(v, result);
1162        }
1163
1164        public ComplexDoubleMatrix mmuli(double v, ComplexDoubleMatrix result) {
1165                return muli(v, result);
1166        }
1167        
1168        /** (Elementwise) division */
1169        public ComplexDoubleMatrix divi(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
1170                if (other.isScalar())
1171                        return divi(other.scalar(), result);
1172                
1173                assertSameLength(other);
1174                ensureResultLength(other, result);
1175                
1176                ComplexDouble c1 = new ComplexDouble(0.0);
1177                ComplexDouble c2 = new ComplexDouble(0.0);
1178                
1179                for (int i = 0; i < length; i++)
1180                        result.put(i, get(i, c1).divi(other.get(i, c2)));
1181                return result;
1182        }
1183                
1184        /** (Elementwise) division with a scalar */
1185        public ComplexDoubleMatrix divi(ComplexDouble a, ComplexDoubleMatrix result) {
1186                ensureResultLength(null, result);
1187                
1188                ComplexDouble c = new ComplexDouble(0.0);
1189                
1190                for (int i = 0; i < length; i++)
1191                        result.put(i, get(i, c).divi(a));
1192                return result;
1193        }       
1194
1195        public ComplexDoubleMatrix divi(double a, ComplexDoubleMatrix result) {
1196                return divi(new ComplexDouble(a), result);
1197        }
1198
1199        /** 
1200         * (Elementwise) division, with operands switched. Computes
1201         * <em>result = other / this</em>. */
1202        public ComplexDoubleMatrix rdivi(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
1203                if (other.isScalar())
1204                        return divi(other.scalar(), result);
1205                
1206                assertSameLength(other);
1207                ensureResultLength(other, result);
1208
1209                ComplexDouble c1 = new ComplexDouble(0.0);
1210                ComplexDouble c2 = new ComplexDouble(0.0);
1211
1212                for (int i = 0; i < length; i++)
1213                        result.put(i, other.get(i, c1).divi(get(i, c2)));
1214                return result;
1215        }
1216                
1217        /** (Elementwise) division with a scalar, with operands switched. Computes
1218         * <em>result = a / this</em>.*/
1219        public ComplexDoubleMatrix rdivi(ComplexDouble a, ComplexDoubleMatrix result) {
1220                ensureResultLength(null, result);
1221
1222                ComplexDouble c1 = new ComplexDouble(0.0);
1223                ComplexDouble c2 = new ComplexDouble(0.0);
1224
1225                for (int i = 0; i < length; i++) {
1226                    c1.copy(a);
1227                    result.put(i, c1.divi(get(i, c2)));                    
1228                }
1229                return result;
1230        }
1231
1232        public ComplexDoubleMatrix rdivi(double a, ComplexDoubleMatrix result) {
1233                return rdivi(new ComplexDouble(a), result);
1234        }
1235        
1236        public ComplexDoubleMatrix negi() {
1237                ComplexDouble c = new ComplexDouble(0.0);
1238                for (int i = 0; i < length; i++)
1239                        put(i, get(i, c).negi());
1240                return this;
1241        }
1242        
1243        public ComplexDoubleMatrix neg() {
1244                return dup().negi();
1245        }
1246
1247        public ComplexDoubleMatrix noti() {
1248                ComplexDouble c = new ComplexDouble(0.0);
1249                for (int i = 0; i < length; i++)
1250                        put(i, get(i, c).isZero() ? 1.0 : 0.0);
1251                return this;
1252        }
1253        
1254        public ComplexDoubleMatrix not() {
1255                return dup().noti();
1256        }
1257        
1258        public ComplexDoubleMatrix truthi() {
1259                ComplexDouble c = new ComplexDouble(0.0);
1260                for (int i = 0; i < length; i++)
1261                        put(i, get(i, c).isZero() ? 0.0 : 1.0);
1262                return this;
1263        }
1264        
1265        public ComplexDoubleMatrix truth() {
1266                return dup().truthi();
1267        }
1268
1269        /****************************************************************
1270         * Rank one-updates
1271         */
1272        
1273        /** Computes a rank-1-update A = A + alpha * x * y'. */ 
1274        public ComplexDoubleMatrix rankOneUpdate(ComplexDouble alpha, ComplexDoubleMatrix x, ComplexDoubleMatrix y) {
1275                if (rows != x.length)
1276                        throw new SizeException("Vector x has wrong length (" + x.length + " != " + rows + ").");
1277                if (columns != y.length)
1278                        throw new SizeException("Vector y has wrong length (" + x.length + " != " + columns + ").");                    
1279                
1280                SimpleBlas.gerc(alpha, x, y, this);
1281                return this;
1282        }
1283
1284        public ComplexDoubleMatrix rankOneUpdate(double alpha, ComplexDoubleMatrix x, ComplexDoubleMatrix y) {
1285                return rankOneUpdate(new ComplexDouble(alpha), x, y);
1286        }
1287
1288        /** Computes a rank-1-update A = A + alpha * x * x'. */ 
1289        public ComplexDoubleMatrix rankOneUpdate(double alpha, ComplexDoubleMatrix x) {
1290                return rankOneUpdate(new ComplexDouble(alpha), x, x);
1291        }
1292
1293        /** Computes a rank-1-update A = A + alpha * x * x'. */ 
1294        public ComplexDoubleMatrix rankOneUpdate(ComplexDouble alpha, ComplexDoubleMatrix x) {
1295                return rankOneUpdate(alpha, x, x);
1296        }
1297
1298        /** Computes a rank-1-update A = A + x * x'. */ 
1299        public ComplexDoubleMatrix rankOneUpdate(ComplexDoubleMatrix x) {
1300                return rankOneUpdate(1.0, x, x);
1301        }
1302
1303        /** Computes a rank-1-update A = A + x * y'. */ 
1304        public ComplexDoubleMatrix rankOneUpdate(ComplexDoubleMatrix x, ComplexDoubleMatrix y) {
1305                return rankOneUpdate(1.0, x, y);
1306        }
1307
1308        /****************************************************************
1309         * Logical operations
1310         */
1311        
1312        public ComplexDouble sum() {
1313                ComplexDouble s = new ComplexDouble(0.0);
1314                ComplexDouble c = new ComplexDouble(0.0);
1315                for (int i = 0; i < length; i++)
1316                        s.addi(get(i, c));
1317                return s;
1318        }
1319        
1320        public ComplexDouble mean() {
1321                return sum().div((double)length);
1322        }
1323        
1324        /** Computes this^T * other */
1325        public ComplexDouble dotc(ComplexDoubleMatrix other) {
1326                return SimpleBlas.dotc(this, other);
1327        }
1328        
1329        /** Computes this^H * other */
1330        public ComplexDouble dotu(ComplexDoubleMatrix other) {
1331                return SimpleBlas.dotu(this, other);
1332        }
1333
1334        public double norm2() {
1335                return SimpleBlas.nrm2(this);
1336        }
1337        
1338        public double normmax() {
1339                int i = SimpleBlas.iamax(this);
1340                return get(i).abs();
1341        }
1342
1343        public double norm1() {
1344                return SimpleBlas.asum(this);
1345        }
1346                
1347        /** Return a vector containing the sums of the columns (having number of columns many entries) */
1348        public ComplexDoubleMatrix columnSums() {
1349                ComplexDoubleMatrix v =
1350                        new ComplexDoubleMatrix(1, columns);
1351
1352                for (int c = 0; c < columns; c++)
1353                        v.put(c, getColumn(c).sum());
1354
1355                return v;
1356        }
1357
1358        public ComplexDoubleMatrix columnMeans() {
1359                return columnSums().divi(rows);
1360        }
1361        
1362        public ComplexDoubleMatrix rowSums() {
1363                ComplexDoubleMatrix v = new ComplexDoubleMatrix(rows);
1364
1365                for (int r = 0; r < rows; r++)
1366                        v.put(r, getRow(r).sum());
1367
1368                return v;
1369        }
1370
1371        public ComplexDoubleMatrix rowMeans() {
1372                return rowSums().divi(columns);
1373        }
1374
1375        public ComplexDoubleMatrix getColumn(int c) {
1376                ComplexDoubleMatrix result = new ComplexDoubleMatrix(rows, 1);
1377                NativeBlas.zcopy(rows, data, index(0, c), 1, result.data, 0, 1);
1378                return result;
1379        }
1380        
1381        public void putColumn(int c, ComplexDoubleMatrix v) {
1382                NativeBlas.zcopy(rows, v.data, 0, 1, data, index(0, c), 1);
1383        }
1384
1385        public ComplexDoubleMatrix getRow(int r) {
1386                ComplexDoubleMatrix result = new ComplexDoubleMatrix(1, columns);
1387                NativeBlas.zcopy(columns, data, index(r, 0), rows, result.data, 0, 1);
1388                return result;
1389        }
1390        
1391        public void putRow(int r, ComplexDoubleMatrix v) {
1392                NativeBlas.zcopy(columns, v.data, 0, 1, data, index(r, 0), rows);
1393        }
1394
1395        /**************************************************************************
1396         * Elementwise Functions
1397         */
1398
1399        /** Add a row vector to all rows of the matrix */
1400        public void addRowVector(ComplexDoubleMatrix x) {
1401                for (int r = 0; r < rows; r++) {
1402                        NativeBlas.zaxpy(columns, ComplexDouble.UNIT, x.data, 0, 1, data, index(r, 0), rows);
1403                }
1404        }
1405
1406        /** Add a vector to all columns of the matrix */
1407        public void addColumnVector(ComplexDoubleMatrix x) {
1408                for (int c = 0; c < columns; c++) {
1409                        NativeBlas.zaxpy(rows, ComplexDouble.UNIT, x.data, 0, 1, data, index(0, c), 1);
1410                }
1411        }
1412
1413        /** Add a row vector to all rows of the matrix */
1414        public void subRowVector(ComplexDoubleMatrix x) {
1415                for (int r = 0; r < rows; r++) {
1416                        NativeBlas.zaxpy(columns, ComplexDouble.NEG_UNIT, x.data, 0, 1, data, index(r, 0), rows);
1417                }
1418        }
1419
1420        /** Add a vector to all columns of the matrix */
1421        public void subColumnVector(ComplexDoubleMatrix x) {
1422                for (int c = 0; c < columns; c++) {
1423                        NativeBlas.zaxpy(rows, ComplexDouble.NEG_UNIT, x.data, 0, 1, data, index(0, c), 1);
1424                }
1425        }
1426
1427        /**
1428         * Writes out this matrix to the given data stream.
1429         * @param dos the data output stream to write to.
1430         * @throws IOException 
1431         */
1432        public void out(DataOutputStream dos) throws IOException {
1433                dos.writeUTF("double");
1434                dos.writeInt(columns);
1435                dos.writeInt(rows);
1436                
1437                dos.writeInt(data.length);
1438                for(int i=0; i < data.length;i++)
1439                        dos.writeDouble(data[i]);
1440        }
1441        
1442        /**
1443         * Reads in a matrix from the given data stream. Note
1444         * that the old data of this matrix will be discarded.
1445         * @param dis the data input stream to read from.
1446         * @throws IOException 
1447         */
1448        public void in(DataInputStream dis) throws IOException {
1449                if(!dis.readUTF().equals("double")) 
1450                        throw new IllegalStateException("The matrix in the specified file is not of the correct type!");
1451                
1452                this.columns    = dis.readInt();
1453                this.rows               = dis.readInt();
1454
1455                final int MAX = dis.readInt();
1456                data = new double[MAX];
1457                for(int i=0; i < MAX;i++)
1458                        data[i] = dis.readDouble();
1459        }       
1460        
1461        /**
1462         * Saves this matrix to the specified file.
1463         * @param filename the file to write the matrix in.
1464         * @throws IOException thrown on errors while writing the matrix to the file
1465         */
1466        public void save(String filename) throws IOException {
1467                DataOutputStream dos = new DataOutputStream(new FileOutputStream(filename, false));
1468                this.out(dos);
1469        }
1470        
1471        /**
1472         * Loads a matrix from a file into this matrix. Note that the old data
1473         * of this matrix will be discarded.
1474         * @param filename the file to read the matrix from
1475         * @throws IOException thrown on errors while reading the matrix
1476         */
1477        public void load(String filename) throws IOException {
1478                DataInputStream dis = new DataInputStream(new FileInputStream(filename));
1479                this.in(dis);
1480        }
1481
1482        /****************************************************************
1483         * Autogenerated code
1484         */
1485        
1486        /***** Code for operators ***************************************/ 
1487
1488        /* Overloads for the usual arithmetic operations */
1489        /*#
1490         def gen_overloads(base, result_rows, result_cols); <<-EOS
1491        public ComplexDoubleMatrix #{base}i(ComplexDoubleMatrix other) {
1492                return #{base}i(other, this);
1493        }
1494                
1495        public ComplexDoubleMatrix #{base}(ComplexDoubleMatrix other) {
1496                return #{base}i(other, new ComplexDoubleMatrix(#{result_rows}, #{result_cols}));
1497        }
1498
1499        public ComplexDoubleMatrix #{base}i(ComplexDouble v) {
1500                return #{base}i(v, this);
1501        }
1502        
1503        public ComplexDoubleMatrix #{base}i(double v) {
1504                return #{base}i(new ComplexDouble(v), this);
1505        }
1506
1507        public ComplexDoubleMatrix #{base}(ComplexDouble v) {
1508                return #{base}i(v, new ComplexDoubleMatrix(rows, columns));
1509        }       
1510
1511        public ComplexDoubleMatrix #{base}(double v) {
1512                return #{base}i(new ComplexDouble(v), new ComplexDoubleMatrix(rows, columns));
1513        }       
1514                
1515                EOS
1516          end
1517        #*/
1518
1519        /* Generating code for logical operators. This not only generates the stubs 
1520         * but really all of the code.
1521         */
1522        
1523        /*#
1524         def gen_compare(name, op); <<-EOS
1525         public ComplexDoubleMatrix #{name}i(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
1526            if (other.isScalar())
1527               return #{name}i(other.scalar(), result);
1528               
1529                assertSameLength(other);
1530                ensureResultLength(other, result);
1531                
1532                ComplexDouble c1 = new ComplexDouble(0.0);
1533                ComplexDouble c2 = new ComplexDouble(0.0);
1534          
1535                for (int i = 0; i < length; i++)
1536                    result.put(i, get(i, c1).#{op}(other.get(i, c2)) ? 1.0 : 0.0);
1537           return result;
1538         }
1539         
1540         public ComplexDoubleMatrix #{name}i(ComplexDoubleMatrix other) {
1541           return #{name}i(other, this);
1542         }
1543         
1544         public ComplexDoubleMatrix #{name}(ComplexDoubleMatrix other) {
1545           return #{name}i(other, new ComplexDoubleMatrix(rows, columns));
1546         }
1547         
1548         public ComplexDoubleMatrix #{name}i(ComplexDouble value, ComplexDoubleMatrix result) {
1549           ensureResultLength(null, result);
1550           ComplexDouble c = new ComplexDouble(0.0);
1551           for (int i = 0; i < length; i++)
1552             result.put(i, get(i, c).#{op}(value) ? 1.0 : 0.0);
1553           return result;
1554         }
1555
1556         public ComplexDoubleMatrix #{name}i(double value, ComplexDoubleMatrix result) {
1557           return #{name}i(new ComplexDouble(value), result);
1558         }
1559
1560         public ComplexDoubleMatrix #{name}i(ComplexDouble value) {
1561           return #{name}i(value, this);
1562         }
1563         
1564         public ComplexDoubleMatrix #{name}i(double value) {
1565           return #{name}i(new ComplexDouble(value));
1566         }
1567         
1568         public ComplexDoubleMatrix #{name}(ComplexDouble value) {
1569           return #{name}i(value, new ComplexDoubleMatrix(rows, columns));
1570         }
1571         
1572         public ComplexDoubleMatrix #{name}(double value) {
1573           return #{name}i(new ComplexDouble(value));
1574         }
1575
1576         EOS
1577         end
1578         #*/
1579        
1580        /*#
1581         def gen_logical(name, op); <<-EOS
1582         public ComplexDoubleMatrix #{name}i(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
1583                assertSameLength(other);
1584                ensureResultLength(other, result);
1585                
1586                ComplexDouble t1 = new ComplexDouble(0.0);
1587                ComplexDouble t2 = new ComplexDouble(0.0);
1588         
1589               for (int i = 0; i < length; i++)
1590                  result.put(i, (!get(i, t1).isZero()) #{op} (!other.get(i, t2).isZero()) ? 1.0 : 0.0);
1591           return result;
1592         }
1593         
1594         public ComplexDoubleMatrix #{name}i(ComplexDoubleMatrix other) {
1595           return #{name}i(other, this);
1596         }
1597         
1598         public ComplexDoubleMatrix #{name}(ComplexDoubleMatrix other) {
1599           return #{name}i(other, new ComplexDoubleMatrix(rows, columns));
1600         }
1601         
1602         public ComplexDoubleMatrix #{name}i(ComplexDouble value, ComplexDoubleMatrix result) {
1603                ensureResultLength(null, result);
1604                boolean val = !value.isZero();
1605                ComplexDouble t = new ComplexDouble(0.0);
1606                for (int i = 0; i < length; i++)
1607                     result.put(i, !get(i, t).isZero() #{op} val ? 1.0 : 0.0);
1608           return result;
1609         }
1610
1611         public ComplexDoubleMatrix #{name}i(double value, ComplexDoubleMatrix result) {
1612           return #{name}i(new ComplexDouble(value), result);
1613         }
1614
1615         public ComplexDoubleMatrix #{name}i(ComplexDouble value) {
1616           return #{name}i(value, this);
1617         }
1618
1619         public ComplexDoubleMatrix #{name}i(double value) {
1620           return #{name}i(new ComplexDouble(value), this);
1621         }
1622
1623         public ComplexDoubleMatrix #{name}(ComplexDouble value) {
1624           return #{name}i(value, new ComplexDoubleMatrix(rows, columns));
1625         }
1626         
1627         public ComplexDoubleMatrix #{name}(double value) {
1628           return #{name}i(new ComplexDouble(value));
1629         }
1630         EOS
1631         end
1632         #*/
1633
1634        /*# collect(gen_overloads('add', 'rows', 'columns'),
1635          gen_overloads('sub', 'rows', 'columns'),
1636          gen_overloads('rsub', 'rows', 'columns'),
1637          gen_overloads('div', 'rows', 'columns'),
1638          gen_overloads('rdiv', 'rows', 'columns'),
1639          gen_overloads('mul', 'rows', 'columns'),
1640          gen_overloads('mmul', 'rows', 'other.columns'),
1641          gen_compare('eq', 'eq'),
1642          gen_compare('ne', 'eq'),
1643          gen_logical('and', '&'),
1644          gen_logical('or', '|'),
1645          gen_logical('xor', '^'))
1646         #*/
1647//RJPP-BEGIN------------------------------------------------------------
1648        public ComplexDoubleMatrix addi(ComplexDoubleMatrix other) {
1649                return addi(other, this);
1650        }
1651                
1652        public ComplexDoubleMatrix add(ComplexDoubleMatrix other) {
1653                return addi(other, new ComplexDoubleMatrix(rows, columns));
1654        }
1655
1656        public ComplexDoubleMatrix addi(ComplexDouble v) {
1657                return addi(v, this);
1658        }
1659        
1660        public ComplexDoubleMatrix addi(double v) {
1661                return addi(new ComplexDouble(v), this);
1662        }
1663
1664        public ComplexDoubleMatrix add(ComplexDouble v) {
1665                return addi(v, new ComplexDoubleMatrix(rows, columns));
1666        }       
1667
1668        public ComplexDoubleMatrix add(double v) {
1669                return addi(new ComplexDouble(v), new ComplexDoubleMatrix(rows, columns));
1670        }       
1671                
1672
1673        public ComplexDoubleMatrix subi(ComplexDoubleMatrix other) {
1674                return subi(other, this);
1675        }
1676                
1677        public ComplexDoubleMatrix sub(ComplexDoubleMatrix other) {
1678                return subi(other, new ComplexDoubleMatrix(rows, columns));
1679        }
1680
1681        public ComplexDoubleMatrix subi(ComplexDouble v) {
1682                return subi(v, this);
1683        }
1684        
1685        public ComplexDoubleMatrix subi(double v) {
1686                return subi(new ComplexDouble(v), this);
1687        }
1688
1689        public ComplexDoubleMatrix sub(ComplexDouble v) {
1690                return subi(v, new ComplexDoubleMatrix(rows, columns));
1691        }       
1692
1693        public ComplexDoubleMatrix sub(double v) {
1694                return subi(new ComplexDouble(v), new ComplexDoubleMatrix(rows, columns));
1695        }       
1696                
1697
1698        public ComplexDoubleMatrix rsubi(ComplexDoubleMatrix other) {
1699                return rsubi(other, this);
1700        }
1701                
1702        public ComplexDoubleMatrix rsub(ComplexDoubleMatrix other) {
1703                return rsubi(other, new ComplexDoubleMatrix(rows, columns));
1704        }
1705
1706        public ComplexDoubleMatrix rsubi(ComplexDouble v) {
1707                return rsubi(v, this);
1708        }
1709        
1710        public ComplexDoubleMatrix rsubi(double v) {
1711                return rsubi(new ComplexDouble(v), this);
1712        }
1713
1714        public ComplexDoubleMatrix rsub(ComplexDouble v) {
1715                return rsubi(v, new ComplexDoubleMatrix(rows, columns));
1716        }       
1717
1718        public ComplexDoubleMatrix rsub(double v) {
1719                return rsubi(new ComplexDouble(v), new ComplexDoubleMatrix(rows, columns));
1720        }       
1721                
1722
1723        public ComplexDoubleMatrix divi(ComplexDoubleMatrix other) {
1724                return divi(other, this);
1725        }
1726                
1727        public ComplexDoubleMatrix div(ComplexDoubleMatrix other) {
1728                return divi(other, new ComplexDoubleMatrix(rows, columns));
1729        }
1730
1731        public ComplexDoubleMatrix divi(ComplexDouble v) {
1732                return divi(v, this);
1733        }
1734        
1735        public ComplexDoubleMatrix divi(double v) {
1736                return divi(new ComplexDouble(v), this);
1737        }
1738
1739        public ComplexDoubleMatrix div(ComplexDouble v) {
1740                return divi(v, new ComplexDoubleMatrix(rows, columns));
1741        }       
1742
1743        public ComplexDoubleMatrix div(double v) {
1744                return divi(new ComplexDouble(v), new ComplexDoubleMatrix(rows, columns));
1745        }       
1746                
1747
1748        public ComplexDoubleMatrix rdivi(ComplexDoubleMatrix other) {
1749                return rdivi(other, this);
1750        }
1751                
1752        public ComplexDoubleMatrix rdiv(ComplexDoubleMatrix other) {
1753                return rdivi(other, new ComplexDoubleMatrix(rows, columns));
1754        }
1755
1756        public ComplexDoubleMatrix rdivi(ComplexDouble v) {
1757                return rdivi(v, this);
1758        }
1759        
1760        public ComplexDoubleMatrix rdivi(double v) {
1761                return rdivi(new ComplexDouble(v), this);
1762        }
1763
1764        public ComplexDoubleMatrix rdiv(ComplexDouble v) {
1765                return rdivi(v, new ComplexDoubleMatrix(rows, columns));
1766        }       
1767
1768        public ComplexDoubleMatrix rdiv(double v) {
1769                return rdivi(new ComplexDouble(v), new ComplexDoubleMatrix(rows, columns));
1770        }       
1771                
1772
1773        public ComplexDoubleMatrix muli(ComplexDoubleMatrix other) {
1774                return muli(other, this);
1775        }
1776                
1777        public ComplexDoubleMatrix mul(ComplexDoubleMatrix other) {
1778                return muli(other, new ComplexDoubleMatrix(rows, columns));
1779        }
1780
1781        public ComplexDoubleMatrix muli(ComplexDouble v) {
1782                return muli(v, this);
1783        }
1784        
1785        public ComplexDoubleMatrix muli(double v) {
1786                return muli(new ComplexDouble(v), this);
1787        }
1788
1789        public ComplexDoubleMatrix mul(ComplexDouble v) {
1790                return muli(v, new ComplexDoubleMatrix(rows, columns));
1791        }       
1792
1793        public ComplexDoubleMatrix mul(double v) {
1794                return muli(new ComplexDouble(v), new ComplexDoubleMatrix(rows, columns));
1795        }       
1796                
1797
1798        public ComplexDoubleMatrix mmuli(ComplexDoubleMatrix other) {
1799                return mmuli(other, this);
1800        }
1801                
1802        public ComplexDoubleMatrix mmul(ComplexDoubleMatrix other) {
1803                return mmuli(other, new ComplexDoubleMatrix(rows, other.columns));
1804        }
1805
1806        public ComplexDoubleMatrix mmuli(ComplexDouble v) {
1807                return mmuli(v, this);
1808        }
1809        
1810        public ComplexDoubleMatrix mmuli(double v) {
1811                return mmuli(new ComplexDouble(v), this);
1812        }
1813
1814        public ComplexDoubleMatrix mmul(ComplexDouble v) {
1815                return mmuli(v, new ComplexDoubleMatrix(rows, columns));
1816        }       
1817
1818        public ComplexDoubleMatrix mmul(double v) {
1819                return mmuli(new ComplexDouble(v), new ComplexDoubleMatrix(rows, columns));
1820        }       
1821                
1822
1823         public ComplexDoubleMatrix eqi(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
1824            if (other.isScalar())
1825               return eqi(other.scalar(), result);
1826               
1827                assertSameLength(other);
1828                ensureResultLength(other, result);
1829                
1830                ComplexDouble c1 = new ComplexDouble(0.0);
1831                ComplexDouble c2 = new ComplexDouble(0.0);
1832          
1833                for (int i = 0; i < length; i++)
1834                    result.put(i, get(i, c1).eq(other.get(i, c2)) ? 1.0 : 0.0);
1835           return result;
1836         }
1837         
1838         public ComplexDoubleMatrix eqi(ComplexDoubleMatrix other) {
1839           return eqi(other, this);
1840         }
1841         
1842         public ComplexDoubleMatrix eq(ComplexDoubleMatrix other) {
1843           return eqi(other, new ComplexDoubleMatrix(rows, columns));
1844         }
1845         
1846         public ComplexDoubleMatrix eqi(ComplexDouble value, ComplexDoubleMatrix result) {
1847           ensureResultLength(null, result);
1848           ComplexDouble c = new ComplexDouble(0.0);
1849           for (int i = 0; i < length; i++)
1850             result.put(i, get(i, c).eq(value) ? 1.0 : 0.0);
1851           return result;
1852         }
1853
1854         public ComplexDoubleMatrix eqi(double value, ComplexDoubleMatrix result) {
1855           return eqi(new ComplexDouble(value), result);
1856         }
1857
1858         public ComplexDoubleMatrix eqi(ComplexDouble value) {
1859           return eqi(value, this);
1860         }
1861         
1862         public ComplexDoubleMatrix eqi(double value) {
1863           return eqi(new ComplexDouble(value));
1864         }
1865         
1866         public ComplexDoubleMatrix eq(ComplexDouble value) {
1867           return eqi(value, new ComplexDoubleMatrix(rows, columns));
1868         }
1869         
1870         public ComplexDoubleMatrix eq(double value) {
1871           return eqi(new ComplexDouble(value));
1872         }
1873
1874
1875         public ComplexDoubleMatrix nei(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
1876            if (other.isScalar())
1877               return nei(other.scalar(), result);
1878               
1879                assertSameLength(other);
1880                ensureResultLength(other, result);
1881                
1882                ComplexDouble c1 = new ComplexDouble(0.0);
1883                ComplexDouble c2 = new ComplexDouble(0.0);
1884          
1885                for (int i = 0; i < length; i++)
1886                    result.put(i, get(i, c1).eq(other.get(i, c2)) ? 1.0 : 0.0);
1887           return result;
1888         }
1889         
1890         public ComplexDoubleMatrix nei(ComplexDoubleMatrix other) {
1891           return nei(other, this);
1892         }
1893         
1894         public ComplexDoubleMatrix ne(ComplexDoubleMatrix other) {
1895           return nei(other, new ComplexDoubleMatrix(rows, columns));
1896         }
1897         
1898         public ComplexDoubleMatrix nei(ComplexDouble value, ComplexDoubleMatrix result) {
1899           ensureResultLength(null, result);
1900           ComplexDouble c = new ComplexDouble(0.0);
1901           for (int i = 0; i < length; i++)
1902             result.put(i, get(i, c).eq(value) ? 1.0 : 0.0);
1903           return result;
1904         }
1905
1906         public ComplexDoubleMatrix nei(double value, ComplexDoubleMatrix result) {
1907           return nei(new ComplexDouble(value), result);
1908         }
1909
1910         public ComplexDoubleMatrix nei(ComplexDouble value) {
1911           return nei(value, this);
1912         }
1913         
1914         public ComplexDoubleMatrix nei(double value) {
1915           return nei(new ComplexDouble(value));
1916         }
1917         
1918         public ComplexDoubleMatrix ne(ComplexDouble value) {
1919           return nei(value, new ComplexDoubleMatrix(rows, columns));
1920         }
1921         
1922         public ComplexDoubleMatrix ne(double value) {
1923           return nei(new ComplexDouble(value));
1924         }
1925
1926
1927         public ComplexDoubleMatrix andi(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
1928                assertSameLength(other);
1929                ensureResultLength(other, result);
1930                
1931                ComplexDouble t1 = new ComplexDouble(0.0);
1932                ComplexDouble t2 = new ComplexDouble(0.0);
1933         
1934               for (int i = 0; i < length; i++)
1935                  result.put(i, (!get(i, t1).isZero()) & (!other.get(i, t2).isZero()) ? 1.0 : 0.0);
1936           return result;
1937         }
1938         
1939         public ComplexDoubleMatrix andi(ComplexDoubleMatrix other) {
1940           return andi(other, this);
1941         }
1942         
1943         public ComplexDoubleMatrix and(ComplexDoubleMatrix other) {
1944           return andi(other, new ComplexDoubleMatrix(rows, columns));
1945         }
1946         
1947         public ComplexDoubleMatrix andi(ComplexDouble value, ComplexDoubleMatrix result) {
1948                ensureResultLength(null, result);
1949                boolean val = !value.isZero();
1950                ComplexDouble t = new ComplexDouble(0.0);
1951                for (int i = 0; i < length; i++)
1952                     result.put(i, !get(i, t).isZero() & val ? 1.0 : 0.0);
1953           return result;
1954         }
1955
1956         public ComplexDoubleMatrix andi(double value, ComplexDoubleMatrix result) {
1957           return andi(new ComplexDouble(value), result);
1958         }
1959
1960         public ComplexDoubleMatrix andi(ComplexDouble value) {
1961           return andi(value, this);
1962         }
1963
1964         public ComplexDoubleMatrix andi(double value) {
1965           return andi(new ComplexDouble(value), this);
1966         }
1967
1968         public ComplexDoubleMatrix and(ComplexDouble value) {
1969           return andi(value, new ComplexDoubleMatrix(rows, columns));
1970         }
1971         
1972         public ComplexDoubleMatrix and(double value) {
1973           return andi(new ComplexDouble(value));
1974         }
1975
1976         public ComplexDoubleMatrix ori(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
1977                assertSameLength(other);
1978                ensureResultLength(other, result);
1979                
1980                ComplexDouble t1 = new ComplexDouble(0.0);
1981                ComplexDouble t2 = new ComplexDouble(0.0);
1982         
1983               for (int i = 0; i < length; i++)
1984                  result.put(i, (!get(i, t1).isZero()) | (!other.get(i, t2).isZero()) ? 1.0 : 0.0);
1985           return result;
1986         }
1987         
1988         public ComplexDoubleMatrix ori(ComplexDoubleMatrix other) {
1989           return ori(other, this);
1990         }
1991         
1992         public ComplexDoubleMatrix or(ComplexDoubleMatrix other) {
1993           return ori(other, new ComplexDoubleMatrix(rows, columns));
1994         }
1995         
1996         public ComplexDoubleMatrix ori(ComplexDouble value, ComplexDoubleMatrix result) {
1997                ensureResultLength(null, result);
1998                boolean val = !value.isZero();
1999                ComplexDouble t = new ComplexDouble(0.0);
2000                for (int i = 0; i < length; i++)
2001                     result.put(i, !get(i, t).isZero() | val ? 1.0 : 0.0);
2002           return result;
2003         }
2004
2005         public ComplexDoubleMatrix ori(double value, ComplexDoubleMatrix result) {
2006           return ori(new ComplexDouble(value), result);
2007         }
2008
2009         public ComplexDoubleMatrix ori(ComplexDouble value) {
2010           return ori(value, this);
2011         }
2012
2013         public ComplexDoubleMatrix ori(double value) {
2014           return ori(new ComplexDouble(value), this);
2015         }
2016
2017         public ComplexDoubleMatrix or(ComplexDouble value) {
2018           return ori(value, new ComplexDoubleMatrix(rows, columns));
2019         }
2020         
2021         public ComplexDoubleMatrix or(double value) {
2022           return ori(new ComplexDouble(value));
2023         }
2024
2025         public ComplexDoubleMatrix xori(ComplexDoubleMatrix other, ComplexDoubleMatrix result) {
2026                assertSameLength(other);
2027                ensureResultLength(other, result);
2028                
2029                ComplexDouble t1 = new ComplexDouble(0.0);
2030                ComplexDouble t2 = new ComplexDouble(0.0);
2031         
2032               for (int i = 0; i < length; i++)
2033                  result.put(i, (!get(i, t1).isZero()) ^ (!other.get(i, t2).isZero()) ? 1.0 : 0.0);
2034           return result;
2035         }
2036         
2037         public ComplexDoubleMatrix xori(ComplexDoubleMatrix other) {
2038           return xori(other, this);
2039         }
2040         
2041         public ComplexDoubleMatrix xor(ComplexDoubleMatrix other) {
2042           return xori(other, new ComplexDoubleMatrix(rows, columns));
2043         }
2044         
2045         public ComplexDoubleMatrix xori(ComplexDouble value, ComplexDoubleMatrix result) {
2046                ensureResultLength(null, result);
2047                boolean val = !value.isZero();
2048                ComplexDouble t = new ComplexDouble(0.0);
2049                for (int i = 0; i < length; i++)
2050                     result.put(i, !get(i, t).isZero() ^ val ? 1.0 : 0.0);
2051           return result;
2052         }
2053
2054         public ComplexDoubleMatrix xori(double value, ComplexDoubleMatrix result) {
2055           return xori(new ComplexDouble(value), result);
2056         }
2057
2058         public ComplexDoubleMatrix xori(ComplexDouble value) {
2059           return xori(value, this);
2060         }
2061
2062         public ComplexDoubleMatrix xori(double value) {
2063           return xori(new ComplexDouble(value), this);
2064         }
2065
2066         public ComplexDoubleMatrix xor(ComplexDouble value) {
2067           return xori(value, new ComplexDoubleMatrix(rows, columns));
2068         }
2069         
2070         public ComplexDoubleMatrix xor(double value) {
2071           return xori(new ComplexDouble(value));
2072         }
2073//RJPP-END--------------------------------------------------------------
2074}