00001
00006 inline long zhematrix::zhesv(zgematrix& mat)
00007 {
00008 #ifdef CPPL_VERBOSE
00009 std::cerr << "# [MARK] zhematrix::zhesv(zgematrix&)"
00010 << std::endl;
00011 #endif//CPPL_VERBOSE
00012
00013 #ifdef CPPL_DEBUG
00014 if(N!=mat.N){
00015 std::cerr << "[ERROR] zhematrix::zhesv(zgematrix&) " << std::endl
00016 << "These two matrices cannot be solved." << std::endl
00017 << "Your input was (" << N << "x" << N << ") and ("
00018 << mat.N << "x" << mat.N << ")." << std::endl;
00019 exit(1);
00020 }
00021 #endif//CPPL_DEBUG
00022
00023 char UPLO('L');
00024 long NRHS(mat.N), LDA(N), *IPIV(new long[N]), LDB(mat.N), LWORK(-1), INFO(1);
00025 std::complex<double> *WORK(new std::complex<double>[1]);
00026 zhesv_(UPLO, N, NRHS, Array, LDA, IPIV, mat.Array, LDB, WORK, LWORK, INFO);
00027
00028 INFO=1;
00029 LWORK = long(std::real(WORK[0]));
00030 delete [] WORK; WORK =new std::complex<double>[LWORK];
00031 zhesv_(UPLO, N, NRHS, Array, LDA, IPIV, mat.Array, LDB, WORK, LWORK, INFO);
00032 delete [] WORK; delete [] IPIV;
00033
00034 if(INFO!=0){
00035 std::cerr << "[WARNING] zhematrix::zhesv(zhematrix&) "
00036 << "Serious trouble happend. INFO = " << INFO << "."
00037 << std::endl;
00038 }
00039 return INFO;
00040 }
00041
00042
00047 inline long zhematrix::zhesv(zcovector& vec)
00048 {
00049 #ifdef CPPL_VERBOSE
00050 std::cerr << "# [MARK] zhematrix::zhesv(zcovector&)"
00051 << std::endl;
00052 #endif//CPPL_VERBOSE
00053
00054 #ifdef CPPL_DEBUG
00055 if(N!=vec.L){
00056 std::cerr << "[ERROR] zhematrix::zhesv(zcovector&) " << std::endl
00057 << "These matrix and vector cannot be solved." << std::endl
00058 << "Your input was (" << N << "x" << N << ") and ("
00059 << vec.L << ")." << std::endl;
00060 exit(1);
00061 }
00062 #endif//CPPL_DEBUG
00063
00064 char UPLO('L');
00065 long NRHS(1), LDA(N), *IPIV(new long[N]), LDB(vec.L), LWORK(-1), INFO(1);
00066 std::complex<double> *WORK( new std::complex<double>[1] );
00067 zhesv_(UPLO, N, NRHS, Array, LDA, IPIV, vec.Array, LDB, WORK, LWORK, INFO);
00068
00069 INFO=1;
00070 LWORK = long(std::real(WORK[0]));
00071 delete WORK; WORK = new std::complex<double>[LWORK];
00072 zhesv_(UPLO, N, NRHS, Array, LDA, IPIV, vec.Array, LDB, WORK, LWORK, INFO);
00073 delete [] WORK; delete [] IPIV;
00074
00075 if(INFO!=0){
00076 std::cerr << "[WARNING] zhematrix::zhesv(zcovector&) "
00077 << "Serious trouble happend. INFO = " << INFO << "."
00078 << std::endl;
00079 }
00080 return INFO;
00081 }
00082
00086
00087
00094 inline long zhematrix::zheev(std::vector<double>& w,
00095 const bool& jobz=0)
00096 {
00097 #ifdef CPPL_VERBOSE
00098 std::cerr << "# [MARK] zhematrix::zheev(std::vector<double>&, const bool&)"
00099 << std::endl;
00100 #endif//CPPL_VERBOSE
00101
00102 w.resize(N);
00103 char JOBZ, UPLO('L');
00104 if(jobz==0){ JOBZ='N'; } else{ JOBZ='V'; }
00105 long LDA(N), INFO(1), LWORK(-1);
00106 double *RWORK(new double[max(1, 3*N-2)]);
00107 std::complex<double> *WORK(new std::complex<double>[1]);
00108 zheev_(JOBZ, UPLO, N, Array, LDA, &w[0], WORK, LWORK, RWORK, INFO);
00109
00110 INFO=1;
00111 LWORK = long(std::real(WORK[0]));
00112 delete [] WORK; WORK = new std::complex<double>[LWORK];
00113 zheev_(JOBZ, UPLO, N, Array, LDA, &w[0], WORK, LWORK, RWORK, INFO);
00114 delete [] RWORK; delete [] WORK;
00115
00116 if(INFO!=0){
00117 std::cerr << "[WARNING] zhematrix::zheev"
00118 << "(vector<std::complex<double>>&, const bool&) "
00119 << "Serious trouble happend. INFO = " << INFO << "."
00120 << std::endl;
00121 }
00122 return INFO;
00123 }
00124
00125
00132 inline long zhematrix::zheev(std::vector<double>& w,
00133 std::vector<zcovector>& v)
00134 {
00135 #ifdef CPPL_VERBOSE
00136 std::cerr << "# [MARK] zhematrix::zheev(std::vector<double>&, std::vector<zcovector>&)"
00137 << std::endl;
00138 #endif//CPPL_VERBOSE
00139
00140 w.resize(N); v.resize(N);
00141 for(long i=0; i<N; i++){ v[i].resize(N); }
00142 char JOBZ('V'), UPLO('L');
00143 long LDA(N), INFO(1), LWORK(-1);
00144 double *RWORK(new double[max(1, 3*N-2)]);
00145 std::complex<double> *WORK(new std::complex<double>[1]);
00146 zheev_(JOBZ, UPLO, N, Array, LDA, &w[0], WORK, LWORK, RWORK, INFO);
00147
00148 INFO=1;
00149 LWORK = long(std::real(WORK[0]));
00150 delete [] WORK; WORK = new std::complex<double>[LWORK];
00151 zheev_(JOBZ, UPLO, N, Array, LDA, &w[0], WORK, LWORK, RWORK, INFO);
00152 delete [] RWORK; delete [] WORK;
00153
00155 for(long i=0; i<N; i++){ for(long j=0; j<N; j++){
00156 v[j](i) = Array[i+N*j];
00157 }}
00158
00159 if(INFO!=0){
00160 std::cerr << "[WARNING] zhematrix::zheev"
00161 << "(vector<std::complex<double>>&, vector<zcovector>&) "
00162 << "Serious trouble happend. INFO = " << INFO << "."
00163 << std::endl;
00164 }
00165 return INFO;
00166 }
00167
00168
00175 inline long zhematrix::zheev(std::vector<double>& w,
00176 std::vector<zrovector>& v)
00177 {
00178 #ifdef CPPL_VERBOSE
00179 std::cerr << "# [MARK] zhematrix::zheev(std::vector<double>&, std::vector<zrovector>&)"
00180 << std::endl;
00181 #endif//CPPL_VERBOSE
00182
00183 w.resize(N); v.resize(N);
00184 for(long i=0; i<N; i++){ v[i].resize(N); }
00185 char JOBZ('V'), UPLO('L');
00186 long LDA(N), INFO(1), LWORK(-1);
00187 double *RWORK(new double[max(1, 3*N-2)]);
00188 std::complex<double> *WORK(new std::complex<double>[1]);
00189 zheev_(JOBZ, UPLO, N, Array, LDA, &w[0], WORK, LWORK, RWORK, INFO);
00190
00191 INFO=1;
00192 LWORK = long(std::real(WORK[0]));
00193 delete [] WORK; WORK = new std::complex<double>[LWORK];
00194 zheev_(JOBZ, UPLO, N, Array, LDA, &w[0], WORK, LWORK, RWORK, INFO);
00195 delete [] RWORK; delete [] WORK;
00196
00198 for(long i=0; i<N; i++){ for(long j=0; j<N; j++){
00199 v[j](i) = Array[i+N*j];
00200 }}
00201
00202 if(INFO!=0){
00203 std::cerr << "[WARNING] zhematrix::zheev"
00204 << "(vector<std::complex<double>>&, vector<zcovector>&)"
00205 << "Serious trouble happend. INFO = " << INFO << "."
00206 << std::endl;
00207 }
00208 return INFO;
00209 }