Skip to content

Re-add the xPPTRF routines (rename xPOTRF routines from #3123) #3229

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
andrew-platt merged 3 commits into from
Mar 3, 2026
Merged

Re-add the xPPTRF routines (rename xPOTRF routines from #3123) #3229

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
75d4c40
Switch from xPPTRF to xPOTRF to improve TurbSim speed on macOS
IrisMeasure Dec 30, 2025
bed151b
nwtclib: update xPOTRF routines, re-add xPPTRF routines
andrew-platt Mar 2, 2026
1495a52
TurbSim: use xPOTRF instead of xPPTRF
andrew-platt Mar 2, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
190 changes: 190 additions & 0 deletions modules/nwtc-library/src/NetLib/lapack/NWTC_LAPACK.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -77,6 +77,12 @@ MODULE NWTC_LAPACK
MODULE PROCEDURE LAPACK_sposv
END INTERFACE

!> Compute the Cholesky factorization of a real symmetric positive definite matrix A stored in packed format (internally handled as unpacked).
INTERFACE LAPACK_potrf
MODULE PROCEDURE LAPACK_dpotrf
MODULE PROCEDURE LAPACK_spotrf
END INTERFACE

!> Compute the Cholesky factorization of a real symmetric positive definite matrix A stored in packed format.
INTERFACE LAPACK_pptrf
MODULE PROCEDURE LAPACK_dpptrf
Expand Down Expand Up @@ -1410,6 +1416,190 @@ SUBROUTINE LAPACK_SPOSV (UPLO, N, NRHS, A, B, ErrStat, ErrMsg)
END SUBROUTINE LAPACK_SPOSV
!=======================================================================
!> Compute the Cholesky factorization of a real symmetric positive definite matrix A stored in packed format.
!! use LAPACK_POTRF (nwtc_lapack::lapack_potrf) instead of this specific function.
SUBROUTINE LAPACK_DPOTRF (UPLO, N, AP, ErrStat, ErrMsg)

! DPOTRF computes the Cholesky factorization of a real symmetric
! positive definite matrix A stored in packed format.
!
! Internally, packed storage is converted to full storage and DPOTRF
! is called. The result is converted back to packed storage.
!
! The factorization has the form
! A = U**T * U, if UPLO = 'U', or
! A = L * L**T, if UPLO = 'L'.

INTEGER, INTENT(IN ) :: N
REAL(R8Ki), INTENT(INOUT) :: AP(:)
INTEGER(IntKi), INTENT( OUT) :: ErrStat
CHARACTER(*), INTENT( OUT) :: ErrMsg
CHARACTER(1), INTENT(IN ) :: UPLO

REAL(R8Ki), ALLOCATABLE :: A(:,:)
INTEGER :: I, J, IDX, INFO

ErrStat = ErrID_None
ErrMsg = ""

IF (N <= 0) RETURN

ALLOCATE(A(N,N))
A = 0.0_R8Ki

IF (UPLO == 'U') THEN
IDX = 0
DO J = 1, N
DO I = 1, J
IDX = IDX + 1
A(I,J) = AP(IDX)
END DO
END DO
ELSE IF (UPLO == 'L') THEN
IDX = 0
DO J = 1, N
DO I = J, N
IDX = IDX + 1
A(I,J) = AP(IDX)
END DO
END DO
ELSE
ErrStat = ErrID_FATAL
ErrMsg = "LAPACK_DPOTRF: invalid UPLO value."
DEALLOCATE(A)
RETURN
END IF

CALL DPOTRF (UPLO, N, A, N, INFO)

IF (INFO /= 0) THEN
ErrStat = ErrID_FATAL
WRITE(ErrMsg,*) INFO
IF (INFO < 0) THEN
ErrMsg = "LAPACK_DPOTRF: illegal value in argument "//TRIM(ErrMsg)//"."
ELSE
ErrMsg = "LAPACK_DPOTRF: Leading minor of order "//TRIM(ErrMsg)// &
" of A is not positive definite, so Cholesky factorization could not be completed."
END IF
DEALLOCATE(A)
RETURN
END IF

IF (UPLO == 'U') THEN
IDX = 0
DO J = 1, N
DO I = 1, J
IDX = IDX + 1
AP(IDX) = A(I,J)
END DO
END DO
ELSE
IDX = 0
DO J = 1, N
DO I = J, N
IDX = IDX + 1
AP(IDX) = A(I,J)
END DO
END DO
END IF

DEALLOCATE(A)
RETURN

END SUBROUTINE LAPACK_DPOTRF
!=======================================================================
!> Compute the Cholesky factorization of a real symmetric positive definite matrix A stored in packed format.
!! use LAPACK_POTRF (nwtc_lapack::lapack_potrf) instead of this specific function.
SUBROUTINE LAPACK_SPOTRF (UPLO, N, AP, ErrStat, ErrMsg)

! SPOTRF computes the Cholesky factorization of a real symmetric
! positive definite matrix A stored in packed format.
!
! Internally, packed storage is converted to full storage and SPOTRF
! is called. The result is converted back to packed storage.
!
! The factorization has the form
! A = U**T * U, if UPLO = 'U', or
! A = L * L**T, if UPLO = 'L'.

INTEGER, INTENT(IN ) :: N
REAL(SiKi), INTENT(INOUT) :: AP(:)
INTEGER(IntKi), INTENT( OUT) :: ErrStat
CHARACTER(*), INTENT( OUT) :: ErrMsg
CHARACTER(1), INTENT(IN ) :: UPLO

REAL(SiKi), ALLOCATABLE :: A(:,:)
INTEGER :: I, J, IDX, INFO

ErrStat = ErrID_None
ErrMsg = ""

IF (N <= 0) RETURN

ALLOCATE(A(N,N))
A = 0.0_SiKi

IF (UPLO == 'U') THEN
IDX = 0
DO J = 1, N
DO I = 1, J
IDX = IDX + 1
A(I,J) = AP(IDX)
END DO
END DO
ELSE IF (UPLO == 'L') THEN
IDX = 0
DO J = 1, N
DO I = J, N
IDX = IDX + 1
A(I,J) = AP(IDX)
END DO
END DO
ELSE
ErrStat = ErrID_FATAL
ErrMsg = "LAPACK_SPOTRF: invalid UPLO value."
DEALLOCATE(A)
RETURN
END IF

CALL SPOTRF (UPLO, N, A, N, INFO)

IF (INFO /= 0) THEN
ErrStat = ErrID_FATAL
WRITE(ErrMsg,*) INFO
IF (INFO < 0) THEN
ErrMsg = "LAPACK_SPOTRF: illegal value in argument "//TRIM(ErrMsg)//"."
ELSE
ErrMsg = "LAPACK_SPOTRF: Leading minor of order "//TRIM(ErrMsg)// &
" of A is not positive definite, so Cholesky factorization could not be completed."
END IF
DEALLOCATE(A)
RETURN
END IF

IF (UPLO == 'U') THEN
IDX = 0
DO J = 1, N
DO I = 1, J
IDX = IDX + 1
AP(IDX) = A(I,J)
END DO
END DO
ELSE
IDX = 0
DO J = 1, N
DO I = J, N
IDX = IDX + 1
AP(IDX) = A(I,J)
END DO
END DO
END IF

DEALLOCATE(A)
RETURN

END SUBROUTINE LAPACK_SPOTRF
!=======================================================================
!> Compute the Cholesky factorization of a real symmetric positive definite matrix A stored in packed format.
!! use LAPACK_PPTRF (nwtc_lapack::lapack_pptrf) instead of this specific function.
SUBROUTINE LAPACK_DPPTRF (UPLO, N, AP, ErrStat, ErrMsg)

Expand Down
2 changes: 1 addition & 1 deletion modules/turbsim/src/TSsubs.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -763,7 +763,7 @@ SUBROUTINE Coh2H( p, IVec, IFreq, TRH, S, ErrStat, ErrMsg )
NPts = p%grid%NPoints
END IF

CALL LAPACK_pptrf( 'L', NPts, TRH(Indx:), ErrStat, ErrMsg ) ! 'L'ower triangular 'TRH' matrix (packed form), of order 'NPoints'; returns Stat
CALL LAPACK_potrf( 'L', NPts, TRH(Indx:), ErrStat, ErrMsg ) ! 'L'ower triangular 'TRH' matrix (unpacked form), of order 'NPoints'; returns Stat

IF ( ErrStat /= ErrID_None ) THEN
IF (ErrStat < AbortErrLev) then
Expand Down