problem with DGELSx: variables A and B unchanged

Post here if you have a question about LAPACK or ScaLAPACK algorithm or data format

problem with DGELSx: variables A and B unchanged

hello!

I am trying to solve a linear least squares problem using LAPACK's specialized routines. Right after DGELS returns, the variables A and B are unchanged but info evaluates to 0. I am using PGI Fortran 7.2-5 with BLAS and LAPACK as released by PGI, on a Dell Precision T7400.

Has anyone else experienced this problem? Is there something I am not doing right? I even tried for a simple 2x2 matrix but the outcome was the same: A and B had the same value as before entering DGELS

Code: Select all
`call DGELS(tr, 2, 2, 1, amat, 2, bvec, 2, ym, -1, infols)`

Code: Select all
`[alexgc@predeal SWE_NEW]\$ ./sensobs  RHS before    1.000000000000000         1.000000000000000         1.000000000000000        -1.000000000000000      LHS before    1.000000000000000         0.000000000000000      info (as returned by DGELS) =             0 RHS after <=> factorization    1.000000000000000         1.000000000000000         1.000000000000000        -1.000000000000000      LHS after <=> solution    1.000000000000000         0.000000000000000     FORTRAN STOP`
alexgc

Posts: 1
Joined: Tue Dec 14, 2010 12:50 pm

Re: problem with DGELSx: variables A and B unchanged

Alex,
I am confused in the dimension of your matrices.
we need to know which matrix is which matrix

Code: Select all
`*  NRHS    (input) INTEGER*          The number of right hand sides, i.e., the number of*          columns of the matrices B and X. NRHS >=0.`

from your output your RHS has 2 colums, so the fourth arg NRHS should be two.

Also only B will be changed, A will stay unchanged.
Code: Select all
`*  B       (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS)*          On entry, the matrix B of right hand side vectors, stored*          columnwise; B is M-by-NRHS if TRANS = 'N', or N-by-NRHS*          if TRANS = 'T'.*          On exit, if INFO = 0, B is overwritten by the solution*          vectors, stored columnwise`

Below is a simple code that works, hope it helps
Julie

Code: Select all
`      program main      double precision A(2,2), B(2,1)      integer m ,n, ldb, I, J      integer info      double precision WORK(66)      lwork=-1      A(1,1)=1      A(2,1)=1      A(1,2)=1      A(2,2)=-1         WRITE(*,*) "MATRIX A"         WRITE(*,99997) ((A(I,J),J=1,2),I=1,2)      B(1,1)=1      B(2,1)=0         WRITE(*,*) "MATRIX B"         WRITE(*,99998) (B(I,1),I=1,2)      lwork=66      call dgels('T',2,2,1,A,2,B,2,work, lwork,info)      write(*,*) "INFO for dgels", INFO      IF (INFO .EQ.0) THEN         write(*,*) "res vaut :"         WRITE(*,99998) (B(I,1),I=1,2)      ELSE         write(*,*) "Problem in DGELS, INFO is not nul"      END IF99997 FORMAT ((3X,2(E16.8)))99998 FORMAT ((3X,1(E16.8)))      end `
admin
Site Admin

Posts: 600
Joined: Wed Dec 08, 2004 7:07 pm

Re: problem with DGELSx: variables A and B unchanged

Hello Alex,

call DGELS(tr, 2, 2, 1, amat, 2, bvec, 2, ym, -1, infols)

You are supplying LWORK = -1. This means that you are doing workspace query. The function just writes recommended length of WORK in the first element of YM and immediately exits. For the workspace query case the result of unchanged matrices is expected. For your case the minimal length of YM and value of LWORK is 4.

From DGELS parameters description:
* LWORK (input) INTEGER
* The dimension of the array WORK.
* LWORK >= max( 1, MN + max( MN, NRHS ) ).
* For optimal performance,
* LWORK >= max( 1, MN + max( MN, NRHS )*NB ).
* where MN = min(M,N) and NB is the optimum block size.
*
* If LWORK = -1, then a workspace query is assumed; the routine
* only calculates the optimal size of the WORK array, returns
* this value as the first entry of the WORK array, and no error
* message related to LWORK is issued by XERBLA.

W.B.R.,
Alexander
akobotov

Posts: 11
Joined: Wed Feb 03, 2010 7:38 am
Location: Intel Corp., Russia, Novosibirsk

Return to Algorithm / Data

Who is online

Users browsing this forum: Bing [Bot] and 2 guests