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
037 package org.jblas;
038
039 /**
040 * <p>Eigenvalue and Eigenvector related functions.</p>
041 *
042 * <p>Methods exist for working with symmetric matrices or general eigenvalues.
043 * The symmetric versions are usually much faster on symmetric matrices.</p>
044 */
045 public class Eigen {
046 private static final DoubleMatrix dummyDouble = new DoubleMatrix(1);
047
048 /** Compute the eigenvalues for a symmetric matrix. */
049 public static DoubleMatrix symmetricEigenvalues(DoubleMatrix A) {
050 A.assertSquare();
051 DoubleMatrix eigenvalues = new DoubleMatrix(A.rows);
052 int isuppz[] = new int[2*A.rows];
053 SimpleBlas.syevr('N', 'A', 'U', A.dup(), 0, 0, 0, 0, 0, eigenvalues, dummyDouble, isuppz);
054 return eigenvalues;
055 }
056
057 /**
058 * Computes the eigenvalues and eigenvectors for a symmetric matrix.
059 *
060 * @return an array of DoubleMatrix objects containing the eigenvectors
061 * stored as the columns of the first matrix, and the eigenvalues as
062 * diagonal elements of the second matrix.
063 */
064 public static DoubleMatrix[] symmetricEigenvectors(DoubleMatrix A) {
065 A.assertSquare();
066 DoubleMatrix eigenvalues = new DoubleMatrix(A.rows);
067 DoubleMatrix eigenvectors = A.dup();
068 int isuppz[] = new int[2*A.rows];
069 SimpleBlas.syevr('V', 'A', 'U', A.dup(), 0, 0, 0, 0, 0, eigenvalues, eigenvectors, isuppz);
070 return new DoubleMatrix[] { eigenvectors, DoubleMatrix.diag(eigenvalues) };
071 }
072
073 /** Computes the eigenvalues of a general matrix. */
074 public static ComplexDoubleMatrix eigenvalues(DoubleMatrix A) {
075 A.assertSquare();
076 DoubleMatrix WR = new DoubleMatrix(A.rows);
077 DoubleMatrix WI = WR.dup();
078 SimpleBlas.geev('N', 'N', A.dup(), WR, WI, dummyDouble, dummyDouble);
079
080 return new ComplexDoubleMatrix(WR, WI);
081 }
082
083 /**
084 * Computes the eigenvalues and eigenvectors of a general matrix.
085 *
086 * @return an array of ComplexDoubleMatrix objects containing the eigenvectors
087 * stored as the columns of the first matrix, and the eigenvalues as the
088 * diagonal elements of the second matrix.
089 */
090 public static ComplexDoubleMatrix[] eigenvectors(DoubleMatrix A) {
091 A.assertSquare();
092 // setting up result arrays
093 DoubleMatrix WR = new DoubleMatrix(A.rows);
094 DoubleMatrix WI = WR.dup();
095 DoubleMatrix VR = new DoubleMatrix(A.rows, A.rows);
096
097 SimpleBlas.geev('N', 'V', A.dup(), WR, WI, dummyDouble, VR);
098
099 // transferring the result
100 ComplexDoubleMatrix E = new ComplexDoubleMatrix(WR, WI);
101 ComplexDoubleMatrix V = new ComplexDoubleMatrix(A.rows, A.rows);
102 //System.err.printf("VR = %s\n", VR.toString());
103 for (int i = 0; i < A.rows; i++) {
104 if (E.get(i).isReal()) {
105 V.putColumn(i, new ComplexDoubleMatrix(VR.getColumn(i)));
106 } else {
107 ComplexDoubleMatrix v = new ComplexDoubleMatrix(VR.getColumn(i), VR.getColumn(i+1));
108 V.putColumn(i, v);
109 V.putColumn(i+1, v.conji());
110 i += 1;
111 }
112 }
113 return new ComplexDoubleMatrix[] { V, ComplexDoubleMatrix.diag(E) };
114 }
115
116 //BEGIN
117 // The code below has been automatically generated.
118 // DO NOT EDIT!
119 private static final FloatMatrix dummyFloat = new FloatMatrix(1);
120
121 /** Compute the eigenvalues for a symmetric matrix. */
122 public static FloatMatrix symmetricEigenvalues(FloatMatrix A) {
123 A.assertSquare();
124 FloatMatrix eigenvalues = new FloatMatrix(A.rows);
125 int isuppz[] = new int[2*A.rows];
126 SimpleBlas.syevr('N', 'A', 'U', A.dup(), 0, 0, 0, 0, 0, eigenvalues, dummyFloat, isuppz);
127 return eigenvalues;
128 }
129
130 /**
131 * Computes the eigenvalues and eigenvectors for a symmetric matrix.
132 *
133 * @return an array of FloatMatrix objects containing the eigenvectors
134 * stored as the columns of the first matrix, and the eigenvalues as
135 * diagonal elements of the second matrix.
136 */
137 public static FloatMatrix[] symmetricEigenvectors(FloatMatrix A) {
138 A.assertSquare();
139 FloatMatrix eigenvalues = new FloatMatrix(A.rows);
140 FloatMatrix eigenvectors = A.dup();
141 int isuppz[] = new int[2*A.rows];
142 SimpleBlas.syevr('V', 'A', 'U', A.dup(), 0, 0, 0, 0, 0, eigenvalues, eigenvectors, isuppz);
143 return new FloatMatrix[] { eigenvectors, FloatMatrix.diag(eigenvalues) };
144 }
145
146 /** Computes the eigenvalues of a general matrix. */
147 public static ComplexFloatMatrix eigenvalues(FloatMatrix A) {
148 A.assertSquare();
149 FloatMatrix WR = new FloatMatrix(A.rows);
150 FloatMatrix WI = WR.dup();
151 SimpleBlas.geev('N', 'N', A.dup(), WR, WI, dummyFloat, dummyFloat);
152
153 return new ComplexFloatMatrix(WR, WI);
154 }
155
156 /**
157 * Computes the eigenvalues and eigenvectors of a general matrix.
158 *
159 * @return an array of ComplexFloatMatrix objects containing the eigenvectors
160 * stored as the columns of the first matrix, and the eigenvalues as the
161 * diagonal elements of the second matrix.
162 */
163 public static ComplexFloatMatrix[] eigenvectors(FloatMatrix A) {
164 A.assertSquare();
165 // setting up result arrays
166 FloatMatrix WR = new FloatMatrix(A.rows);
167 FloatMatrix WI = WR.dup();
168 FloatMatrix VR = new FloatMatrix(A.rows, A.rows);
169
170 SimpleBlas.geev('N', 'V', A.dup(), WR, WI, dummyFloat, VR);
171
172 // transferring the result
173 ComplexFloatMatrix E = new ComplexFloatMatrix(WR, WI);
174 ComplexFloatMatrix V = new ComplexFloatMatrix(A.rows, A.rows);
175 //System.err.printf("VR = %s\n", VR.toString());
176 for (int i = 0; i < A.rows; i++) {
177 if (E.get(i).isReal()) {
178 V.putColumn(i, new ComplexFloatMatrix(VR.getColumn(i)));
179 } else {
180 ComplexFloatMatrix v = new ComplexFloatMatrix(VR.getColumn(i), VR.getColumn(i+1));
181 V.putColumn(i, v);
182 V.putColumn(i+1, v.conji());
183 i += 1;
184 }
185 }
186 return new ComplexFloatMatrix[] { V, ComplexFloatMatrix.diag(E) };
187 }
188
189 //END
190 }