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.la; 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 //BEGIN 116 // The code below has been automatically generated. 117 // DO NOT EDIT! 118 private static final FloatMatrix dummyFloat = new FloatMatrix(1); 119 120 /** Compute the eigenvalues for a symmetric matrix. */ 121 public static FloatMatrix symmetricEigenvalues(FloatMatrix A) { 122 A.assertSquare(); 123 FloatMatrix eigenvalues = new FloatMatrix(A.rows); 124 int isuppz[] = new int[2*A.rows]; 125 SimpleBlas.syevr('N', 'A', 'U', A.dup(), 0, 0, 0, 0, 0, eigenvalues, dummyFloat, isuppz); 126 return eigenvalues; 127 } 128 129 /** 130 * Computes the eigenvalues and eigenvectors for a symmetric matrix. 131 * 132 * @return an array of FloatMatrix objects containing the eigenvectors 133 * stored as the columns of the first matrix, and the eigenvalues as 134 * diagonal elements of the second matrix. 135 */ 136 public static FloatMatrix[] symmetricEigenvectors(FloatMatrix A) { 137 A.assertSquare(); 138 FloatMatrix eigenvalues = new FloatMatrix(A.rows); 139 FloatMatrix eigenvectors = A.dup(); 140 int isuppz[] = new int[2*A.rows]; 141 SimpleBlas.syevr('V', 'A', 'U', A.dup(), 0, 0, 0, 0, 0, eigenvalues, eigenvectors, isuppz); 142 return new FloatMatrix[] { eigenvectors, FloatMatrix.diag(eigenvalues) }; 143 } 144 145 /** Computes the eigenvalues of a general matrix. */ 146 public static ComplexFloatMatrix eigenvalues(FloatMatrix A) { 147 A.assertSquare(); 148 FloatMatrix WR = new FloatMatrix(A.rows); 149 FloatMatrix WI = WR.dup(); 150 SimpleBlas.geev('N', 'N', A.dup(), WR, WI, dummyFloat, dummyFloat); 151 152 return new ComplexFloatMatrix(WR, WI); 153 } 154 155 /** 156 * Computes the eigenvalues and eigenvectors of a general matrix. 157 * 158 * @return an array of ComplexFloatMatrix objects containing the eigenvectors 159 * stored as the columns of the first matrix, and the eigenvalues as the 160 * diagonal elements of the second matrix. 161 */ 162 public static ComplexFloatMatrix[] eigenvectors(FloatMatrix A) { 163 A.assertSquare(); 164 // setting up result arrays 165 FloatMatrix WR = new FloatMatrix(A.rows); 166 FloatMatrix WI = WR.dup(); 167 FloatMatrix VR = new FloatMatrix(A.rows, A.rows); 168 169 SimpleBlas.geev('N', 'V', A.dup(), WR, WI, dummyFloat, VR); 170 171 // transferring the result 172 ComplexFloatMatrix E = new ComplexFloatMatrix(WR, WI); 173 ComplexFloatMatrix V = new ComplexFloatMatrix(A.rows, A.rows); 174 System.err.printf("VR = %s\n", VR.toString()); 175 for (int i = 0; i < A.rows; i++) { 176 if (E.get(i).isReal()) { 177 V.putColumn(i, new ComplexFloatMatrix(VR.getColumn(i))); 178 } else { 179 ComplexFloatMatrix v = new ComplexFloatMatrix(VR.getColumn(i), VR.getColumn(i+1)); 180 V.putColumn(i, v); 181 V.putColumn(i+1, v.conji()); 182 i += 1; 183 } 184 } 185 return new ComplexFloatMatrix[] { V, ComplexFloatMatrix.diag(E) }; 186 } 187 //END 188 }