lapack icon indicating copy to clipboard operation
lapack copied to clipboard

This a draft for dgecx.f.

Open scr2016 opened this issue 3 months ago • 2 comments

This a draft for dgecx.f. to REVIEW and COMMENT

scr2016 avatar Oct 15 '25 07:10 scr2016

by mistake 2 branches (one with dgecx.f and one with DGECX.f ) were joined in a single pull request, the files are identical. DGECX.f can be removed.

On Wed, Oct 15, 2025 at 11:15 AM Mark Gates @.***> wrote:

@.**** requested changes on this pull request.

Some initial comments.

On SRC/DGECX.f https://github.com/Reference-LAPACK/lapack/pull/1161#discussion_r2433413718 :

There is both DGECX.f (incorrect) and dgecx.f (correct) file? This file should be removed.

In SRC/dgecx.f https://github.com/Reference-LAPACK/lapack/pull/1161#discussion_r2433415475 :

+> array as the first K elements. +> b) If the flag COMP_X is set, then the routine also explicitly +> computes and returns the factor X = pseudoinv(C) * A. +> +> \endverbatim + +* Arguments: +* ========== +* +> \param[in] FACT +> \verbatim +> FACT is CHARACTER1 +> Specifies how the factors of CX factorization +> are returned. +> +> = 'P' or 'p' : return only the column permutaion matrix P

Pedantic: LAPACK doesn't normally specify both Uppercase and lowercase. Cf.:

lapack/SRC> grep "^*> += '\w'" lapack/SRC/zpotrf.f *> = 'U': Upper triangle of A is stored; *> = 'L': Lower triangle of A is stored.

lapack/SRC> grep "^*> += '\w'" lapack/SRC/*.f`


In SRC/dgecx.f https://github.com/Reference-LAPACK/lapack/pull/1161#discussion_r2433417924 :

+> +> 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. +> \endverbatim +> +> \param[out] IWORK +> \verbatim +> IWORK is INTEGER array, dimension (N). +> Is a work array. ( IWORK is used by DGEQP3RK to store indices +> of "bad" columns for norm downdating in the residual +> matrix in the blocked step auxiliary subroutine DLAQP3RK ). +> \endverbatim +> +> \param[out] LIWORK

Input, right? \param[in]

In SRC/dgecx.f https://github.com/Reference-LAPACK/lapack/pull/1161#discussion_r2433420002 :

+> +> \param[out] IWORK +> \verbatim +> IWORK is INTEGER array, dimension (N). +> Is a work array. ( IWORK is used by DGEQP3RK to store indices +> of "bad" columns for norm downdating in the residual +> matrix in the blocked step auxiliary subroutine DLAQP3RK ). +> \endverbatim +> +> \param[out] LIWORK +> \verbatim +> LIWORK is INTEGER +> The dimension of the array LIWORK. LIWORK >= N +> +> If LIWORK = -1, then a workspace query is assumed; the routine +> only calculates the optimal size of the WORK array, returns

Typo: IWORK array.

In SRC/dgecx.f https://github.com/Reference-LAPACK/lapack/pull/1161#discussion_r2433428427 :

+* Copy matrix C into work array WORK. +*

  •     CALL DLACPY( 'F', M, K, C, LDC, WORK, M )
    

+*

  •     CALL DGELS( 'N', M, K, N, WORK, M, X, LDX,
    
  • $               WORK( M*K+1 ), LWORK,
    
  • $               IINFO )
    
  •     INFO = IINFO
    

+*

  •  END IF
    

+*

  •  WORK( 1 ) = DBLE( LWKOPT )
    

+* +* DGECX +*

  •  END
    

Files should end with a newline. Notice ⊖ on Github.

In SRC/dgecx.f https://github.com/Reference-LAPACK/lapack/pull/1161#discussion_r2433448149 :

+* code. +*

  •  IF( INFO.EQ.0 ) THEN
    
  •     MINMN = MIN( M, N )
    
  •     IF( MINMN.EQ.0 ) THEN
    
  •        LWKMIN = 1
    
  •        LWKOPT = 1
    
  •     ELSE
    
  •        LWKMIN = 3*N + 1
    

+* +* Assign to NBOPT optimal block size. +*

  •        NBOPT = ILAENV( INB, 'DGEQRF', ' ', M, N, -1, -1 )
    
  •        LWKOPT =  1000
    
  •     END IF
    
  •     WORK( 1 ) = DBLE( LWKOPT )
    

Add: LWORK( 1 ) = N or whatever the formula is. Also elsewhere.

In SRC/dgecx.f https://github.com/Reference-LAPACK/lapack/pull/1161#discussion_r2433452890 :

+* that uses QR factorization. For that purpose, we store C into +* WORK array WORK(1:MK), and the matrix A was copied into + the array X at the begining of the routine. +* +* Copy matrix C into work array WORK. +*

  •     CALL DLACPY( 'F', M, K, C, LDC, WORK, M )
    

+*

  •     CALL DGELS( 'N', M, K, N, WORK, M, X, LDX,
    
  • $               WORK( M*K+1 ), LWORK,
    
  • $               IINFO )
    
  •     INFO = IINFO
    

+*

  •  END IF
    

+*

  •  WORK( 1 ) = DBLE( LWKOPT )
    

Add: LWORK( 1 ) = N as above.

In SRC/dgecx.f https://github.com/Reference-LAPACK/lapack/pull/1161#discussion_r2433471710 :

+*

  •  MNSUB = MIN(MSUB, NSUB)
    
  •  MRESID = MSUB-NSEL
    
  •  NRESID = NSUB-NSEL
    
  •  IF( NSEL.GT.0 ) THEN
    
  •     IF( MSUB.LT.NSEL ) THEN
    

+* TODO: Move this part to the top of the routine. +* a) Case MSUB < NSEL. +* When the number of preselected columns +* is larger than MSUB, hence the factorization of all NSEL +* columns cannot be completed. Return from the routine with the +* error of COL_SEL_DESEL parameter. NSEL cannot be larger than +* MSUB. +*

  •        INFO = -6
    
  •        WORK( 1 ) = DBLE( LWKOPT )
    

Add: LWORK( 1 ) = N or whatever the formula is. Also elsewhere.

— Reply to this email directly, view it on GitHub https://github.com/Reference-LAPACK/lapack/pull/1161#pullrequestreview-3341523325, or unsubscribe https://github.com/notifications/unsubscribe-auth/AAHYAZG2IG6CWKRVGDULGWT3X2FMPAVCNFSM6AAAAACJHEAYW6VHI2DSMVQWIX3LMV43YUDVNRWFEZLROVSXG5CSMV3GSZLXHMZTGNBRGUZDGMZSGU . You are receiving this because you authored the thread.Message ID: @.***>

scr2016 avatar Oct 16 '25 05:10 scr2016

We've discussed that this is the expert interface, and we would have a simplified interface (e.g., without all the column selection). Elsewhere in LAPACK, expert interfaces have an extra x, such as gesvx, gesvxx, posvx, posvxx. Should we rename this one gecxx and create the simplified interface as gecx?

mgates3 avatar Nov 05 '25 22:11 mgates3