Eng-Tips is the largest engineering community on the Internet

Intelligent Work Forums for Engineering Professionals

doubt using dsygv

Status
Not open for further replies.

rk_19

Structural
Aug 7, 2012
71
0
0
AE
hi, I am using fortran to get the Eigen values for generalized Eigen value problem of type Ax=lamda.Bx. I am using lapack routine dsygv for this and it works perfectly for all my 5 case studies and not working as expected for the 6th one. i get first two Eigen values as -ve which was not expected. I do not see anything which could cause this error. both matrices are symmetric square matrices, all diagonals are positive. I am attaching the A and B matrices here. someone could look into them to see what is going wrong. appreciate any guidance on this.
 
 http://files.engineering.com/getfile.aspx?folder=74216901-ea42-4c9b-b8b6-895acc958e15&file=A_and_B_matrix.zip
Replies continue below

Recommended for you

instead of looking at the one 6-th case alone, I would focus on comparing this 6-th case to the other 5 and see what kind of attributes are common to those 5 and different to the 6-th; in addition to what you mentioned and other physically meaningful things regarding your problem, look to things from the computing point of view, like actual size vs allocated memory, whether they are all being passed in the same style (pointers or not), etc.

Are you solving all 6 problems in a row? If so, sometimes is good to change the order of things...maybe something is getting corrupted in memory by the time you get to 6...solve 6 first, instead or something.


 
the 6 cases are independent, and I am running them individually, same code is used for all 6 cases. they are all similar problems from numerical point of view or based on the input required for DSYGV to work. the routine requires 2 matrices which are square and symmetric and I checked these matrices before they go into the routine and it meets the criteria. is there anyway to check what is happening - what causes a -ve Eigen value ?? does that mean the B matrix is not positive definite ?? I checked B matrix (attached in earlier post), its a square matrix, symmetric about diagonal, and all diagonals are positive - in Mathcad, I tried cholesky decomposition and it successfully created lower triangular matrix for matrix B.... which confirms B is positive definite.... thanks
 
hello, to crosscheck my calcs, I took the A and B matrix coming out of fortran into Mathcad (screenshot attached) - I tried to get the eigenvalues in Mathcad and I got the e1, e2,e3 exactly as expected.so why is it not working with dsygv. i found an interesting thing when i compared mathcad output with my dsygv output(refer attachment), the eigen value e1 in mathcad matches eigenval3 from dsygv, e2 from mathcad matches eigenval 4 in fortran and so on.. this looks like the eigen value 1 and 2( which are negative) is inserted in the dsygv listing by some error ???
Link :
 
IS THERE A LAPACK/BLAS ROUTINE TO USE WHEN THE A,B MATRICES ARE NOT PRECISELY SYMMETRIC ? I FOUND THAT MY MATRICES ARE NOT EXACTLY SYMERTRIC - DUE TO PRECISION RELATED ISSUES ?
 
I had a look at this using the Python Scipy library, which includes a front end to the Lapack functions. I called that from Excel with xlwings, using the two Excel files you posted in another forum (the ones you linked to here seem to be text files, and had some problems with the numbers being corrupted when I imported them).

With the matrices taken from your files, and using the Scipy front-end to dsygv, the first two results are negative as you found. I wrote a VBA routine to copy the upper triangle to the lower triangle, to ensure that they were symmetrical, and after doing that only the first result was negative. Also the following results are approximately equal to the Matlab results, but diverge after the 3rd significant figure.

In the attached file, on sheet1, I have copied all the Scipy function results to column 2, and I have converted the function call to text in cell A1. Let me know if you want a copy of the required VBA and Python code to run the function. You would need to install Python, Scipy, and xlwings to run it.

The Scipy manual was not very helpful, it just lists the input and output for the lapack function. Scipy may have some alternative code to do what you want but it wasn't obvious to me if it does.

Doug Jenkins
Interactive Design Services
 
thanks IDS for all that effort. in your opinion, why does Mathcad show results different from the lapack results (I am taking the same matrices from fortran to excel and then using it in Mathcad).infact the Mathcad results were matching output from an independent software too. The numerical work was done to calculate natural period for a structure (for which I already have a independent FE model in sacs software). the eigvalue1 in Mathcad and sacs are matching (excuse minor differences), this is matching eig3 of fortran o/p. similary eig2 of Mathcad/sacs match eig4 of fortran o/p. eign3 is not of my interest. eig 4 of sacs/Mathcad match fortran eig 6 and so on. please see the attachment. results. txt is fortran output.once again, thanks.
the matrices A, B and the fortran program is attached
 
I don't know why Lapack gives different results, but I checked your matrices using the Scipy Eig function, and got the same results as your Mathcad results.

I had assumed the Eig function just called the Lapack function, but it seems not. Is there a reason why you don't want to use Mathcad?

If you want to call the scipy functions from Excel you can download a spreadsheet from:

Doug Jenkins
Interactive Design Services
 
I am developing this program in fortran, so cannot use Mathcad for this work. Mathcad was used just to check the discrepancy. if you see, after discarding the first two Eigen values in lapack o/p, the rest matches with Mathcad - strange !!
 
If you need to do it entirely in Fortran it should be possible to extract the Fortran code from the Scipy eig function, since it is all open source, but I doubt that would be easy.

Scipy Linalg source:
The link below may (or may not) help:


Did you get anything useful from the other forums?

If you do find a way to get Fortran code that works, I'd be interested to hear it.

Doug Jenkins
Interactive Design Services
 
Looking at the Scipy manual and Wikipedia, it looks like ARPACK is worth investigating.

From the Scipy manual:

Imagine you’d like to find the smallest and largest eigenvalues and the corresponding eigenvectors for a
large matrix. ARPACK can handle many forms of input: dense matrices such as numpy.ndarray instances,
sparse matrices such as scipy.sparse.csr_matrix, or a general linear operator derived from
scipy.sparse.linalg.LinearOperator. For this example, for simplicity, we’ll construct a symmetric,
positive-definite matrix.

Doug Jenkins
Interactive Design Services
 
Status
Not open for further replies.
Back
Top