49 #include "blockmemshell/memory.h"
50 #include "scip/type_retcode.h"
54 #define BMS_CALL(x) do \
58 SCIPerrorMessage("No memory in function call\n"); \
59 return SCIP_NOMEMORY; \
65 #define DOUBLETOINT(x) ((int) (x + 0.5))
76 char* JOBZ,
char* RANGE,
char* UPLO,
77 int* N,
double* A,
int* LDA,
78 double* VL,
double* VU,
80 double* ABSTOL,
int* M,
double* W,
double* Z,
81 int* LDZ,
int* ISUPPZ,
double* WORK,
82 int* LWORK,
int* IWORK,
int* LIWORK,
87 void F77_FUNC(dgemv, DGEMV)(
char* TRANS,
int* M,
int* N,
double* ALPHA,
double* A,
int* LDA,
double* X,
int* INCX,
double* BETA,
double* Y,
int* INCY);
103 SCIP_Bool geteigenvectors,
107 SCIP_Real* eigenvalue,
108 SCIP_Real* eigenvector
134 assert( blkmem != NULL );
137 assert( 0 < i && i <= n );
138 assert( eigenvalue != NULL );
139 assert( ( ! geteigenvectors) || eigenvector != NULL );
142 JOBZ = geteigenvectors ?
'V' :
'N';
157 F77_FUNC(dsyevr, DSYEVR)( &JOBZ, &RANGE, &UPLO,
161 &ABSTOL, &M, NULL, NULL,
163 &LWORK, &WISIZE, &LIWORK,
168 SCIPerrorMessage(
"There was an error when calling DSYEVR. INFO = %d\n", INFO);
176 BMS_CALL( BMSallocBlockMemoryArray(blkmem, &WORK, LWORK) );
177 BMS_CALL( BMSallocBlockMemoryArray(blkmem, &IWORK, LIWORK) );
178 BMS_CALL( BMSallocBlockMemoryArray(blkmem, &WTMP, N) );
179 BMS_CALL( BMSallocBlockMemoryArray(blkmem, &ISUPPZ, 2) );
185 F77_FUNC(dsyevr, DSYEVR)( &JOBZ, &RANGE, &UPLO,
189 &ABSTOL, &M, WTMP, eigenvector,
191 &LWORK, IWORK, &LIWORK,
196 SCIPerrorMessage(
"There was an error when calling DSYEVR. INFO = %d\n", INFO);
201 *eigenvalue = WTMP[0];
204 BMSfreeBlockMemoryArray(blkmem, &ISUPPZ, 2);
205 BMSfreeBlockMemoryArray(blkmem, &WTMP, N);
206 BMSfreeBlockMemoryArray(blkmem, &IWORK, LIWORK);
207 BMSfreeBlockMemoryArray(blkmem, &WORK, LWORK);
245 F77_FUNC(dgemv, DGEMV)(&TRANS, &M, &N, &ALPHA, A, &LDA, X, &INCX, &BETA, Y, &INCY);
Macros for calling Fortran-functions.
SCIP_RETCODE SCIPlapackMatrixVectorMult(int nrows, int ncols, SCIP_Real *matrix, SCIP_Real *vector, SCIP_Real *result)
#define F77_FUNC(name, NAME)
interface methods for eigenvector computation and matrix multiplication using different versions of L...
EXTERN SCIP_RETCODE SCIPlapackComputeIthEigenvalue(BMS_BLKMEM *blkmem, SCIP_Bool geteigenvectors, int n, SCIP_Real *A, int i, SCIP_Real *eigenvalue, SCIP_Real *eigenvector)