From 75d4c408aaf2b9e9485266297d599329931e8737 Mon Sep 17 00:00:00 2001 From: Typhaceae Date: Tue, 30 Dec 2025 20:28:58 +0800 Subject: [PATCH 01/27] Switch from xPPTRF to xPOTRF to improve TurbSim speed on macOS --- .../src/NetLib/lapack/NWTC_LAPACK.f90 | 280 +++++++++++------- 1 file changed, 169 insertions(+), 111 deletions(-) diff --git a/modules/nwtc-library/src/NetLib/lapack/NWTC_LAPACK.f90 b/modules/nwtc-library/src/NetLib/lapack/NWTC_LAPACK.f90 index bf48d36268..a9e898f8be 100644 --- a/modules/nwtc-library/src/NetLib/lapack/NWTC_LAPACK.f90 +++ b/modules/nwtc-library/src/NetLib/lapack/NWTC_LAPACK.f90 @@ -1412,127 +1412,185 @@ END SUBROUTINE LAPACK_SPOSV !> 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) - - ! DPPTRF computes the Cholesky factorization of a real symmetric - ! positive definite matrix A stored in packed format. - ! - ! The factorization has the form - ! A = U**T * U, if UPLO = 'U', or - ! A = L * L**T, if UPLO = 'L', - ! where U is an upper triangular matrix and L is lower triangular. - - ! passed parameters - - INTEGER, intent(in ) :: N !< The order of the matrix A. N >= 0. - - ! .. Array Arguments .. - REAL(R8Ki) ,intent(inout) :: AP( : ) !< AP is REAL array, dimension (N*(N+1)/2) - !! On entry, the upper or lower triangle of the symmetric matrix A, packed columnwise in a linear array. The j-th column of A - !! is stored in the array AP as follows: - !! if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j; - !! if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n. - !! See below for further details. - !! On exit, if INFO = 0, the triangular factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T, in the same storage format as A. - !! - !! Further details: - !! The packed storage scheme is illustrated by the following example - !! when N = 4, UPLO = 'U': - !! - !! Two-dimensional storage of the symmetric matrix A: - !! - !! a11 a12 a13 a14 - !! a22 a23 a24 - !! a33 a34 (aij = aji) - !! a44 - !! - !! Packed storage of the upper triangle of A: - !! - !! AP = [ a11, a12, a22, a13, a23, a33, a14, a24, a34, a44 ] - - - INTEGER(IntKi), intent( out) :: ErrStat !< Error level - CHARACTER(*), intent( out) :: ErrMsg !< Message describing error - CHARACTER(1), intent(in ) :: UPLO !< 'U': Upper triangle of A is stored; 'L': Lower triangle of A is stored. - - - - ! local variables - INTEGER :: INFO ! = 0: successful exit; < 0: if INFO = -i, the i-th argument had an illegal value; - ! > 0: if INFO = i, the leading minor of order i is not positive definite, and the factorization could not be completed. - - - ErrStat = ErrID_None - ErrMsg = "" - - CALL DPPTRF (UPLO, N, AP, INFO) - - IF (INFO /= 0) THEN - ErrStat = ErrID_FATAL - WRITE( ErrMsg, * ) INFO - IF (INFO < 0) THEN - ErrMsg = "LAPACK_DPPTRF: illegal value in argument "//TRIM(ErrMsg)//"." + ! DPPTRF 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 - ErrMsg = 'LAPACK_DPPTRF: Leading minor order '//TRIM(ErrMsg)//' of A is not positive definite, so Cholesky factorization could not be completed.' + ErrStat = ErrID_FATAL + ErrMsg = "LAPACK_DPPTRF: invalid UPLO value." + DEALLOCATE(A) + RETURN END IF - END IF - - - RETURN - END SUBROUTINE LAPACK_DPPTRF + + CALL DPOTRF (UPLO, N, A, N, INFO) + + IF (INFO /= 0) THEN + ErrStat = ErrID_FATAL + WRITE(ErrMsg,*) INFO + IF (INFO < 0) THEN + ErrMsg = "LAPACK_DPPTRF: illegal value in argument "//TRIM(ErrMsg)//"." + ELSE + ErrMsg = "LAPACK_DPPTRF: 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_DPPTRF !======================================================================= !> 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_SPPTRF (UPLO, N, AP, ErrStat, ErrMsg) - - ! SPPTRF computes the Cholesky factorization of a real symmetric - ! positive definite matrix A stored in packed format. - ! - ! The factorization has the form - ! A = U**T * U, if UPLO = 'U', or - ! A = L * L**T, if UPLO = 'L', - ! where U is an upper triangular matrix and L is lower triangular. - - ! passed parameters - - INTEGER, intent(in ) :: N !< The order of the matrix A. N >= 0. - - ! .. Array Arguments .. - REAL(SiKi) ,intent(inout) :: AP( : ) !< AP is REAL array, dimension (N*(N+1)/2) - !! On entry, the upper or lower triangle of the symmetric matrix A, packed columnwise in a linear array. The j-th column of A - !! is stored in the array AP as follows: - !! if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j; - !! if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n. - !! See LAPACK_DPPTRF for further details. - !! On exit, if INFO = 0, the triangular factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T, in the same storage format as A. - - INTEGER(IntKi), intent( out) :: ErrStat !< Error level - CHARACTER(*), intent( out) :: ErrMsg !< Message describing error - CHARACTER(1), intent(in ) :: UPLO !< 'U': Upper triangle of A is stored; 'L': Lower triangle of A is stored. - - ! local variables - INTEGER :: INFO ! = 0: successful exit; < 0: if INFO = -i, the i-th argument had an illegal value; - ! > 0: if INFO = i, the leading minor of order i is not positive definite, and the factorization could not be completed - - - ErrStat = ErrID_None - ErrMsg = "" - - CALL SPPTRF (UPLO, N, AP, INFO) - - IF (INFO /= 0) THEN - ErrStat = ErrID_FATAL - WRITE( ErrMsg, * ) INFO - IF (INFO < 0) THEN - ErrMsg = "LAPACK_SPPTRF: illegal value in argument "//TRIM(ErrMsg)//"." + ! SPPTRF 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 - ErrMsg = 'LAPACK_SPPTRF: Leading minor order '//TRIM(ErrMsg)//' of A is not positive definite, so Cholesky factorization could not be completed.' + ErrStat = ErrID_FATAL + ErrMsg = "LAPACK_SPPTRF: invalid UPLO value." + DEALLOCATE(A) + RETURN END IF - END IF - - - RETURN + + CALL SPOTRF (UPLO, N, A, N, INFO) + + IF (INFO /= 0) THEN + ErrStat = ErrID_FATAL + WRITE(ErrMsg,*) INFO + IF (INFO < 0) THEN + ErrMsg = "LAPACK_SPPTRF: illegal value in argument "//TRIM(ErrMsg)//"." + ELSE + ErrMsg = "LAPACK_SPPTRF: 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_SPPTRF !======================================================================= !> Compute singular value decomposition (SVD) for a general matrix, A. From ab6c179cbf423c7ae2b65b018bbc76493cc983df Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Thu, 19 Feb 2026 15:13:41 -0700 Subject: [PATCH 02/27] ADI c-bind: comment out code that causes segfaults - seg fault if `debug_level > 0` and the file names are passed in --- .../aerodyn/src/AeroDyn_Inflow_C_Binding.f90 | 18 ++++++++++-------- 1 file changed, 10 insertions(+), 8 deletions(-) diff --git a/modules/aerodyn/src/AeroDyn_Inflow_C_Binding.f90 b/modules/aerodyn/src/AeroDyn_Inflow_C_Binding.f90 index 22ecc7d1b6..df76ffa0e5 100644 --- a/modules/aerodyn/src/AeroDyn_Inflow_C_Binding.f90 +++ b/modules/aerodyn/src/AeroDyn_Inflow_C_Binding.f90 @@ -935,18 +935,20 @@ subroutine ShowPassedData() call WrScr(" ADinputFilePassed_C "//TmpFlag ) call WrScr(" ADinputFileString_C (ptr addr) "//trim(Num2LStr(LOC(ADinputFileString_C))) ) call WrScr(" ADinputFileStringLength_C "//trim(Num2LStr( ADinputFileStringLength_C )) ) - if (ADinputFilePassed==0_c_int) then - i = index(ADinputFileString, char(0)) ! skip anything after c_null_char - call WrScr(" ADinputFileString_C "//ADinputFileString(1:i)) - endif +!FIXME: This causes a seg fault +! if (ADinputFilePassed==0_c_int) then +! i = index(ADinputFileString, char(0)) ! skip anything after c_null_char +! call WrScr(" ADinputFileString_C "//ADinputFileString(1:i)) +! endif TmpFlag="F"; if (IfWinputFilePassed==1_c_int) TmpFlag="T" call WrScr(" IfWinputFilePassed_C "//TmpFlag ) call WrScr(" IfWinputFileString_C (ptr addr)"//trim(Num2LStr(LOC(IfWinputFileString_C))) ) call WrScr(" IfWinputFileStringLength_C "//trim(Num2LStr( IfWinputFileStringLength_C )) ) - if (IfWinputFilePassed==0_c_int) then - i = index(IfWinputFileString, char(0)) ! skip anything after c_null_char - call WrScr(" IfWinputFileString_C "//trim(IfWinputFileString(1:i))) - endif +!FIXME: This causes a seg fault +! if (IfWinputFilePassed==0_c_int) then +! i = index(IfWinputFileString, char(0)) ! skip anything after c_null_char +! call WrScr(" IfWinputFileString_C "//trim(IfWinputFileString(1:i))) +! endif call WrScr(" OutRootName "//trim(OutRootName) ) call WrScr(" Interpolation") call WrScr(" InterpOrder_C "//trim(Num2LStr( InterpOrder_C )) ) From 3434712225163aff8fb80f567aca82339805ae94 Mon Sep 17 00:00:00 2001 From: Derek Slaughter Date: Fri, 27 Feb 2026 16:06:42 +0000 Subject: [PATCH 03/27] Move !$OMP THREADPRIVATE(TRH) in TSsubs.f90 to after the variable declaration for IFX compiler --- modules/turbsim/src/TSsubs.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/modules/turbsim/src/TSsubs.f90 b/modules/turbsim/src/TSsubs.f90 index fd8753c5ab..d0bc39fa27 100644 --- a/modules/turbsim/src/TSsubs.f90 +++ b/modules/turbsim/src/TSsubs.f90 @@ -55,8 +55,8 @@ SUBROUTINE CalcFourierCoeffs_IEC( p, U, PhaseAngles, S, V, TRH_in, ErrStat, ErrM INTEGER(IntKi), INTENT(OUT) :: ErrStat CHARACTER(*), INTENT(OUT) :: ErrMsg -!$OMP THREADPRIVATE(TRH) REAL(ReKi), allocatable, save :: TRH(:) ! Each OMP thread gets its own copy of this array +!$OMP THREADPRIVATE(TRH) REAL(ReKi), ALLOCATABLE :: Dist(:) ! The distance between points REAL(ReKi), ALLOCATABLE :: DistU(:) From 458a11172d1ee0faeca925e11034308206bd6071 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Mon, 2 Mar 2026 10:00:15 -0700 Subject: [PATCH 04/27] [BugFix] AD: DTAero from input file was ignored --- modules/aerodyn/src/AeroDyn.f90 | 2 +- modules/aerodyn/src/AeroDyn_IO.f90 | 5 +---- 2 files changed, 2 insertions(+), 5 deletions(-) diff --git a/modules/aerodyn/src/AeroDyn.f90 b/modules/aerodyn/src/AeroDyn.f90 index f70b8d71b7..641c17c6cc 100644 --- a/modules/aerodyn/src/AeroDyn.f90 +++ b/modules/aerodyn/src/AeroDyn.f90 @@ -370,7 +370,7 @@ subroutine AD_Init( InitInp, u, p, x, xd, z, OtherState, y, m, Interval, InitOut if (Failed()) return; ! ----------------------------------------------------------------- ! Read the AeroDyn blade files, or copy from passed input - call ReadInputFiles( InitInp%InputFile, InputFileData, interval, p%RootName, NumBlades, AeroProjMod, UnEcho, calcCrvAngle, ErrStat2, ErrMsg2 ) + call ReadInputFiles( InitInp%InputFile, InputFileData, p%RootName, NumBlades, AeroProjMod, UnEcho, calcCrvAngle, ErrStat2, ErrMsg2 ) if (Failed()) return; ! override some parameters to simplify for aero maps diff --git a/modules/aerodyn/src/AeroDyn_IO.f90 b/modules/aerodyn/src/AeroDyn_IO.f90 index 3e14564239..fa0717bf7e 100644 --- a/modules/aerodyn/src/AeroDyn_IO.f90 +++ b/modules/aerodyn/src/AeroDyn_IO.f90 @@ -588,14 +588,12 @@ end subroutine Calc_WriteOutput_FVW END SUBROUTINE Calc_WriteOutput !---------------------------------------------------------------------------------------------------------------------------------- -SUBROUTINE ReadInputFiles( InputFileName, InputFileData, Default_DT, OutFileRoot, NumBlades, AeroProjMod, UnEcho, calcCrvAngle, ErrStat, ErrMsg ) +SUBROUTINE ReadInputFiles( InputFileName, InputFileData, OutFileRoot, NumBlades, AeroProjMod, UnEcho, calcCrvAngle, ErrStat, ErrMsg ) ! This subroutine reads the input file and stores all the data in the AD_InputFile structure. ! It does not perform data validation. !.................................................................................................................................. ! Passed variables - REAL(DbKi), INTENT(IN) :: Default_DT ! The default DT (from glue code) - CHARACTER(*), INTENT(IN) :: InputFileName ! Name of the input file CHARACTER(*), INTENT(IN) :: OutFileRoot ! The rootname of all the output files written by this routine. @@ -624,7 +622,6 @@ SUBROUTINE ReadInputFiles( InputFileName, InputFileData, Default_DT, OutFileRoot ErrStat = ErrID_None ErrMsg = '' - InputFileData%DTAero = Default_DT ! the glue code's suggested DT for the module (may be overwritten in ReadPrimaryFile()) calcCrvAngle = .false. ! initialize in case of early return ! get the blade input-file data From 872c1dbb809116e03e577c5326996bd22eef1626 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Mon, 2 Mar 2026 15:36:48 -0700 Subject: [PATCH 05/27] regtest: check if destination dirs exist before copyTree call --- reg_tests/executeAerodiskRegressionCase.py | 3 ++- reg_tests/executeAerodynPyRegressionCase.py | 3 ++- reg_tests/executeAerodynRegressionCase.py | 7 ++++--- reg_tests/executeHydrodynPyRegressionCase.py | 3 ++- reg_tests/executeMoordynPyRegressionCase.py | 3 ++- reg_tests/executeOpenfastLinearRegressionCase.py | 4 ++-- reg_tests/executeSimpleElastodynRegressionCase.py | 3 ++- reg_tests/executeSubdynRegressionCase.py | 3 ++- reg_tests/executeUnsteadyAeroRegressionCase.py | 3 ++- 9 files changed, 20 insertions(+), 12 deletions(-) diff --git a/reg_tests/executeAerodiskRegressionCase.py b/reg_tests/executeAerodiskRegressionCase.py index ae1dcd5eec..06b8624c88 100644 --- a/reg_tests/executeAerodiskRegressionCase.py +++ b/reg_tests/executeAerodiskRegressionCase.py @@ -88,7 +88,8 @@ if not os.path.isdir(inputsDirectory): rtl.exitWithError("The test data inputs directory, {}, does not exist. Verify your local repository is up to date.".format(inputsDirectory)) -rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict={'adsk_driver.outb':'adsk_driver_ref.outb'}) +if not os.path.isdir(testBuildDirectory): + rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict={'adsk_driver.outb':'adsk_driver_ref.outb'}) ### Run aerodisk on the test case if not noExec: diff --git a/reg_tests/executeAerodynPyRegressionCase.py b/reg_tests/executeAerodynPyRegressionCase.py index e9de4f3bb0..93a654fd90 100644 --- a/reg_tests/executeAerodynPyRegressionCase.py +++ b/reg_tests/executeAerodynPyRegressionCase.py @@ -90,7 +90,8 @@ # create the local output directory and initialize it with input files -rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict={'py_ad_driver.out':'py_ad_driver_ref.out'}) +if not os.path.isdir(dst): + rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict={'py_ad_driver.out':'py_ad_driver_ref.out'}) # , excludeExt=['.out','.outb']) ### Run aerodyn on the test case diff --git a/reg_tests/executeAerodynRegressionCase.py b/reg_tests/executeAerodynRegressionCase.py index 3c6cbeac8f..efb71d54ca 100644 --- a/reg_tests/executeAerodynRegressionCase.py +++ b/reg_tests/executeAerodynRegressionCase.py @@ -92,16 +92,17 @@ # Special case, copy the BAR Baseline files dst = os.path.join(buildDirectory, "BAR_Baseline") src = os.path.join(moduleDirectory, "BAR_Baseline") -try: +if not os.path.isdir(dst): rtl.copyTree(src, dst) -except: +else: # This can fail if two processes are copying the file at the same time print('>>> Copy failed') import time time.sleep(1) # create the local output directory and initialize it with input files -rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict={'ad_driver.outb':'ad_driver_ref.outb'}) +if not os.path.isdir(testBuildDirectory): + rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict={'ad_driver.outb':'ad_driver_ref.outb'}) # , excludeExt=['.out','.outb']) ### Run aerodyn on the test case diff --git a/reg_tests/executeHydrodynPyRegressionCase.py b/reg_tests/executeHydrodynPyRegressionCase.py index 56cf5a8bf3..d121932b6f 100644 --- a/reg_tests/executeHydrodynPyRegressionCase.py +++ b/reg_tests/executeHydrodynPyRegressionCase.py @@ -90,7 +90,8 @@ # create the local output directory and initialize it with input files -rtl.copyTree(inputsDirectory, testBuildDirectory) +if not os.path.isdir(testBuildDirectory): + rtl.copyTree(inputsDirectory, testBuildDirectory) # , excludeExt=['.out','.outb']) ### Run HydroDyn on the test case diff --git a/reg_tests/executeMoordynPyRegressionCase.py b/reg_tests/executeMoordynPyRegressionCase.py index 0f43de1e66..912d90d2fd 100644 --- a/reg_tests/executeMoordynPyRegressionCase.py +++ b/reg_tests/executeMoordynPyRegressionCase.py @@ -90,7 +90,8 @@ rtl.exitWithError("The test data inputs directory, {}, does not exist. Verify your local repository is up to date.".format(inputsDirectory)) # create the local output directory and initialize it with input files -rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict={'MD.out':'MD_ref.out','MDroot.MD.out':'MDroot.MD_ref.out'}) +if not os.path.isdir(testBuildDirectory): + rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict={'MD.out':'MD_ref.out','MDroot.MD.out':'MDroot.MD_ref.out'}) # , excludeExt=['.out','.outb']) ### Run MoorDyn on the test case diff --git a/reg_tests/executeOpenfastLinearRegressionCase.py b/reg_tests/executeOpenfastLinearRegressionCase.py index 91f43062b3..c5bb641041 100644 --- a/reg_tests/executeOpenfastLinearRegressionCase.py +++ b/reg_tests/executeOpenfastLinearRegressionCase.py @@ -163,8 +163,8 @@ def indent(msg, sindent='\t'): shutil.copy2(srcname, dstname) # # Copying the actual test directory -# if not os.path.isdir(testBuildDirectory): -rtl.copyTree(inputsDirectory, testBuildDirectory, excludeExt=excludeExt, renameExtDict={'.lin':'.ref_lin'}) +if not os.path.isdir(testBuildDirectory): + rtl.copyTree(inputsDirectory, testBuildDirectory, excludeExt=excludeExt, renameExtDict={'.lin':'.ref_lin'}) ### Run openfast on the test case if not noExec: diff --git a/reg_tests/executeSimpleElastodynRegressionCase.py b/reg_tests/executeSimpleElastodynRegressionCase.py index e271872574..96af0cea29 100644 --- a/reg_tests/executeSimpleElastodynRegressionCase.py +++ b/reg_tests/executeSimpleElastodynRegressionCase.py @@ -88,7 +88,8 @@ if not os.path.isdir(inputsDirectory): rtl.exitWithError("The test data inputs directory, {}, does not exist. Verify your local repository is up to date.".format(inputsDirectory)) -rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict={'sed_driver.outb':'sed_driver_ref.outb'}) +f not os.path.isdir(testBuildDirectory): + rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict={'sed_driver.outb':'sed_driver_ref.outb'}) ### Run SED on the test case if not noExec: diff --git a/reg_tests/executeSubdynRegressionCase.py b/reg_tests/executeSubdynRegressionCase.py index 04363472d4..b4a10fad06 100644 --- a/reg_tests/executeSubdynRegressionCase.py +++ b/reg_tests/executeSubdynRegressionCase.py @@ -89,7 +89,8 @@ rtl.exitWithError("The test data inputs directory, {}, does not exist. Verify your local repository is up to date.".format(inputsDirectory)) # create the local output directory and initialize it with input files (overwrite if exists) -rtl.copyTree(inputsDirectory, testBuildDirectory, renameExtDict={'.out':'_ref.out'}, includeExt=['.dvr','.dat','.out','.csv']) +f not os.path.isdir(testBuildDirectory): + rtl.copyTree(inputsDirectory, testBuildDirectory, renameExtDict={'.out':'_ref.out'}, includeExt=['.dvr','.dat','.out','.csv']) ### Run SubDyn on the test case diff --git a/reg_tests/executeUnsteadyAeroRegressionCase.py b/reg_tests/executeUnsteadyAeroRegressionCase.py index 53c8bf567e..f57be33779 100644 --- a/reg_tests/executeUnsteadyAeroRegressionCase.py +++ b/reg_tests/executeUnsteadyAeroRegressionCase.py @@ -92,7 +92,8 @@ # create the local output directory and initialize it with input files renameDict={'UA'+str(i)+'.outb':'UA'+str(i)+'_ref.outb' for i in [2,3,4,5,6,7]} -rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict=renameDict +f not os.path.isdir(testBuildDirectory): + rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict=renameDict , excludeExt=['.sum']) # , excludeExt=['.out','.outb','.sum']) From e99a8aa2bd56d7dcf6c8bf81b9852a3033205831 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Mon, 2 Mar 2026 16:01:55 -0700 Subject: [PATCH 06/27] fix typos in previous commit --- reg_tests/executeSimpleElastodynRegressionCase.py | 2 +- reg_tests/executeSubdynRegressionCase.py | 2 +- reg_tests/executeUnsteadyAeroRegressionCase.py | 2 +- 3 files changed, 3 insertions(+), 3 deletions(-) diff --git a/reg_tests/executeSimpleElastodynRegressionCase.py b/reg_tests/executeSimpleElastodynRegressionCase.py index 96af0cea29..a5232038b4 100644 --- a/reg_tests/executeSimpleElastodynRegressionCase.py +++ b/reg_tests/executeSimpleElastodynRegressionCase.py @@ -88,7 +88,7 @@ if not os.path.isdir(inputsDirectory): rtl.exitWithError("The test data inputs directory, {}, does not exist. Verify your local repository is up to date.".format(inputsDirectory)) -f not os.path.isdir(testBuildDirectory): +if not os.path.isdir(testBuildDirectory): rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict={'sed_driver.outb':'sed_driver_ref.outb'}) ### Run SED on the test case diff --git a/reg_tests/executeSubdynRegressionCase.py b/reg_tests/executeSubdynRegressionCase.py index b4a10fad06..ef96940713 100644 --- a/reg_tests/executeSubdynRegressionCase.py +++ b/reg_tests/executeSubdynRegressionCase.py @@ -89,7 +89,7 @@ rtl.exitWithError("The test data inputs directory, {}, does not exist. Verify your local repository is up to date.".format(inputsDirectory)) # create the local output directory and initialize it with input files (overwrite if exists) -f not os.path.isdir(testBuildDirectory): +if not os.path.isdir(testBuildDirectory): rtl.copyTree(inputsDirectory, testBuildDirectory, renameExtDict={'.out':'_ref.out'}, includeExt=['.dvr','.dat','.out','.csv']) diff --git a/reg_tests/executeUnsteadyAeroRegressionCase.py b/reg_tests/executeUnsteadyAeroRegressionCase.py index f57be33779..94a11e613f 100644 --- a/reg_tests/executeUnsteadyAeroRegressionCase.py +++ b/reg_tests/executeUnsteadyAeroRegressionCase.py @@ -92,7 +92,7 @@ # create the local output directory and initialize it with input files renameDict={'UA'+str(i)+'.outb':'UA'+str(i)+'_ref.outb' for i in [2,3,4,5,6,7]} -f not os.path.isdir(testBuildDirectory): +if not os.path.isdir(testBuildDirectory): rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict=renameDict , excludeExt=['.sum']) # , excludeExt=['.out','.outb','.sum']) From dc6520e31b7b4e2d559e643e821f638bde1d3d07 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Mon, 2 Mar 2026 16:25:10 -0700 Subject: [PATCH 07/27] regtest: ad -- remove else from the copyTree wrap --- reg_tests/executeAerodynRegressionCase.py | 5 ----- 1 file changed, 5 deletions(-) diff --git a/reg_tests/executeAerodynRegressionCase.py b/reg_tests/executeAerodynRegressionCase.py index efb71d54ca..30fe3c2b90 100644 --- a/reg_tests/executeAerodynRegressionCase.py +++ b/reg_tests/executeAerodynRegressionCase.py @@ -94,11 +94,6 @@ src = os.path.join(moduleDirectory, "BAR_Baseline") if not os.path.isdir(dst): rtl.copyTree(src, dst) -else: - # This can fail if two processes are copying the file at the same time - print('>>> Copy failed') - import time - time.sleep(1) # create the local output directory and initialize it with input files if not os.path.isdir(testBuildDirectory): From bed151b6fd8843086d546a8ee4a57b80e84af526 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Mon, 2 Mar 2026 16:50:58 -0700 Subject: [PATCH 08/27] nwtclib: update xPOTRF routines, re-add xPPTRF routines commit 75d4c408a changed the xPPTRF routines to call the xPOTRF lapack routines, but we want to keep both. --- .../src/NetLib/lapack/NWTC_LAPACK.f90 | 158 ++++++++++++++++-- 1 file changed, 145 insertions(+), 13 deletions(-) diff --git a/modules/nwtc-library/src/NetLib/lapack/NWTC_LAPACK.f90 b/modules/nwtc-library/src/NetLib/lapack/NWTC_LAPACK.f90 index a9e898f8be..3561c98162 100644 --- a/modules/nwtc-library/src/NetLib/lapack/NWTC_LAPACK.f90 +++ b/modules/nwtc-library/src/NetLib/lapack/NWTC_LAPACK.f90 @@ -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 @@ -1410,10 +1416,10 @@ 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_PPTRF (nwtc_lapack::lapack_pptrf) instead of this specific function. - SUBROUTINE LAPACK_DPPTRF (UPLO, N, AP, ErrStat, ErrMsg) +!! use LAPACK_POTRF (nwtc_lapack::lapack_potrf) instead of this specific function. + SUBROUTINE LAPACK_DPOTRF (UPLO, N, AP, ErrStat, ErrMsg) - ! DPPTRF computes the Cholesky factorization of a real symmetric + ! 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 @@ -1458,7 +1464,7 @@ SUBROUTINE LAPACK_DPPTRF (UPLO, N, AP, ErrStat, ErrMsg) END DO ELSE ErrStat = ErrID_FATAL - ErrMsg = "LAPACK_DPPTRF: invalid UPLO value." + ErrMsg = "LAPACK_DPOTRF: invalid UPLO value." DEALLOCATE(A) RETURN END IF @@ -1469,9 +1475,9 @@ SUBROUTINE LAPACK_DPPTRF (UPLO, N, AP, ErrStat, ErrMsg) ErrStat = ErrID_FATAL WRITE(ErrMsg,*) INFO IF (INFO < 0) THEN - ErrMsg = "LAPACK_DPPTRF: illegal value in argument "//TRIM(ErrMsg)//"." + ErrMsg = "LAPACK_DPOTRF: illegal value in argument "//TRIM(ErrMsg)//"." ELSE - ErrMsg = "LAPACK_DPPTRF: Leading minor of order "//TRIM(ErrMsg)// & + 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) @@ -1499,13 +1505,13 @@ SUBROUTINE LAPACK_DPPTRF (UPLO, N, AP, ErrStat, ErrMsg) DEALLOCATE(A) RETURN - END SUBROUTINE LAPACK_DPPTRF + END SUBROUTINE LAPACK_DPOTRF !======================================================================= !> 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_SPPTRF (UPLO, N, AP, ErrStat, ErrMsg) +!! use LAPACK_POTRF (nwtc_lapack::lapack_potrf) instead of this specific function. + SUBROUTINE LAPACK_SPOTRF (UPLO, N, AP, ErrStat, ErrMsg) - ! SPPTRF computes the Cholesky factorization of a real symmetric + ! 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 @@ -1550,7 +1556,7 @@ SUBROUTINE LAPACK_SPPTRF (UPLO, N, AP, ErrStat, ErrMsg) END DO ELSE ErrStat = ErrID_FATAL - ErrMsg = "LAPACK_SPPTRF: invalid UPLO value." + ErrMsg = "LAPACK_SPOTRF: invalid UPLO value." DEALLOCATE(A) RETURN END IF @@ -1561,9 +1567,9 @@ SUBROUTINE LAPACK_SPPTRF (UPLO, N, AP, ErrStat, ErrMsg) ErrStat = ErrID_FATAL WRITE(ErrMsg,*) INFO IF (INFO < 0) THEN - ErrMsg = "LAPACK_SPPTRF: illegal value in argument "//TRIM(ErrMsg)//"." + ErrMsg = "LAPACK_SPOTRF: illegal value in argument "//TRIM(ErrMsg)//"." ELSE - ErrMsg = "LAPACK_SPPTRF: Leading minor of order "//TRIM(ErrMsg)// & + 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) @@ -1591,6 +1597,132 @@ SUBROUTINE LAPACK_SPPTRF (UPLO, N, AP, ErrStat, ErrMsg) 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) + + ! DPPTRF computes the Cholesky factorization of a real symmetric + ! positive definite matrix A stored in packed format. + ! + ! The factorization has the form + ! A = U**T * U, if UPLO = 'U', or + ! A = L * L**T, if UPLO = 'L', + ! where U is an upper triangular matrix and L is lower triangular. + + + ! passed parameters + + INTEGER, intent(in ) :: N !< The order of the matrix A. N >= 0. + + ! .. Array Arguments .. + REAL(R8Ki) ,intent(inout) :: AP( : ) !< AP is REAL array, dimension (N*(N+1)/2) + !! On entry, the upper or lower triangle of the symmetric matrix A, packed columnwise in a linear array. The j-th column of A + !! is stored in the array AP as follows: + !! if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j; + !! if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n. + !! See below for further details. + !! On exit, if INFO = 0, the triangular factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T, in the same storage format as A. + !! + !! Further details: + !! The packed storage scheme is illustrated by the following example + !! when N = 4, UPLO = 'U': + !! + !! Two-dimensional storage of the symmetric matrix A: + !! + !! a11 a12 a13 a14 + !! a22 a23 a24 + !! a33 a34 (aij = aji) + !! a44 + !! + !! Packed storage of the upper triangle of A: + !! + !! AP = [ a11, a12, a22, a13, a23, a33, a14, a24, a34, a44 ] + + + INTEGER(IntKi), intent( out) :: ErrStat !< Error level + CHARACTER(*), intent( out) :: ErrMsg !< Message describing error + CHARACTER(1), intent(in ) :: UPLO !< 'U': Upper triangle of A is stored; 'L': Lower triangle of A is stored. + + + + ! local variables + INTEGER :: INFO ! = 0: successful exit; < 0: if INFO = -i, the i-th argument had an illegal value; + ! > 0: if INFO = i, the leading minor of order i is not positive definite, and the factorization could not be completed. + + + ErrStat = ErrID_None + ErrMsg = "" + + CALL DPPTRF (UPLO, N, AP, INFO) + + IF (INFO /= 0) THEN + ErrStat = ErrID_FATAL + WRITE( ErrMsg, * ) INFO + IF (INFO < 0) THEN + ErrMsg = "LAPACK_DPPTRF: illegal value in argument "//TRIM(ErrMsg)//"." + ELSE + ErrMsg = 'LAPACK_DPPTRF: Leading minor order '//TRIM(ErrMsg)//' of A is not positive definite, so Cholesky factorization could not be completed.' + END IF + END IF + + + RETURN + END SUBROUTINE LAPACK_DPPTRF +!======================================================================= +!> 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_SPPTRF (UPLO, N, AP, ErrStat, ErrMsg) + + ! SPPTRF computes the Cholesky factorization of a real symmetric + ! positive definite matrix A stored in packed format. + ! + ! The factorization has the form + ! A = U**T * U, if UPLO = 'U', or + ! A = L * L**T, if UPLO = 'L', + ! where U is an upper triangular matrix and L is lower triangular. + + + ! passed parameters + + INTEGER, intent(in ) :: N !< The order of the matrix A. N >= 0. + + ! .. Array Arguments .. + REAL(SiKi) ,intent(inout) :: AP( : ) !< AP is REAL array, dimension (N*(N+1)/2) + !! On entry, the upper or lower triangle of the symmetric matrix A, packed columnwise in a linear array. The j-th column of A + !! is stored in the array AP as follows: + !! if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j; + !! if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n. + !! See LAPACK_DPPTRF for further details. + !! On exit, if INFO = 0, the triangular factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T, in the same storage format as A. + + INTEGER(IntKi), intent( out) :: ErrStat !< Error level + CHARACTER(*), intent( out) :: ErrMsg !< Message describing error + CHARACTER(1), intent(in ) :: UPLO !< 'U': Upper triangle of A is stored; 'L': Lower triangle of A is stored. + + ! local variables + INTEGER :: INFO ! = 0: successful exit; < 0: if INFO = -i, the i-th argument had an illegal value; + ! > 0: if INFO = i, the leading minor of order i is not positive definite, and the factorization could not be completed + + + ErrStat = ErrID_None + ErrMsg = "" + + CALL SPPTRF (UPLO, N, AP, INFO) + + IF (INFO /= 0) THEN + ErrStat = ErrID_FATAL + WRITE( ErrMsg, * ) INFO + IF (INFO < 0) THEN + ErrMsg = "LAPACK_SPPTRF: illegal value in argument "//TRIM(ErrMsg)//"." + ELSE + ErrMsg = 'LAPACK_SPPTRF: Leading minor order '//TRIM(ErrMsg)//' of A is not positive definite, so Cholesky factorization could not be completed.' + END IF + END IF + + + RETURN END SUBROUTINE LAPACK_SPPTRF !======================================================================= !> Compute singular value decomposition (SVD) for a general matrix, A. From 1495a52447eaf8c4620cac8e88be01a8cd3bac82 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Mon, 2 Mar 2026 16:38:41 -0700 Subject: [PATCH 09/27] TurbSim: use xPOTRF instead of xPPTRF See PR #3123 --- modules/turbsim/src/TSsubs.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/modules/turbsim/src/TSsubs.f90 b/modules/turbsim/src/TSsubs.f90 index fd8753c5ab..4c0dbc5485 100644 --- a/modules/turbsim/src/TSsubs.f90 +++ b/modules/turbsim/src/TSsubs.f90 @@ -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 From 2e23291e155a1aa152f67ac21c00ebf53ad52956 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Mon, 2 Mar 2026 17:10:30 -0700 Subject: [PATCH 10/27] regtest: incorrect path check on py_ad reg testing --- reg_tests/executeAerodynPyRegressionCase.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/reg_tests/executeAerodynPyRegressionCase.py b/reg_tests/executeAerodynPyRegressionCase.py index 93a654fd90..66bcb3c233 100644 --- a/reg_tests/executeAerodynPyRegressionCase.py +++ b/reg_tests/executeAerodynPyRegressionCase.py @@ -90,7 +90,7 @@ # create the local output directory and initialize it with input files -if not os.path.isdir(dst): +if not os.path.isdir(testBuildDirectory): rtl.copyTree(inputsDirectory, testBuildDirectory, renameDict={'py_ad_driver.out':'py_ad_driver_ref.out'}) # , excludeExt=['.out','.outb']) From ef8544be13da7a32289699071565f9b378af84d9 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Tue, 3 Mar 2026 13:25:34 -0700 Subject: [PATCH 11/27] [BugFix] AD driver: hub acceleration value incorrect in RotMotionType==1 --- modules/aerodyn/src/AeroDyn_Driver_Subs.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/modules/aerodyn/src/AeroDyn_Driver_Subs.f90 b/modules/aerodyn/src/AeroDyn_Driver_Subs.f90 index 8f2483211d..2328a3530f 100644 --- a/modules/aerodyn/src/AeroDyn_Driver_Subs.f90 +++ b/modules/aerodyn/src/AeroDyn_Driver_Subs.f90 @@ -820,7 +820,7 @@ subroutine Set_Mesh_Motion(nt, dvr, ADI, FED, errStat, errMsg) call interpTimeValue(wt%hub%motion, time, wt%hub%iMotion, hubMotion) !print*,hubMotion wt%hub%rotSpeed = hubMotion(2) - wt%hub%rotAcc = hubMotion(2) + wt%hub%rotAcc = hubMotion(3) wt%hub%azimuth = MODULO(hubMotion(1)*R2D, 360.0_ReKi ) else if (wt%hub%motionType == idHubMotionUserFunction) then ! We call a user-defined function to determined the azimuth, speed (and potentially acceleration...) From 09fc1be906662fcca1479ecdd7f7a667534bd3e5 Mon Sep 17 00:00:00 2001 From: Roger Bergua Date: Tue, 3 Mar 2026 15:47:48 -0700 Subject: [PATCH 12/27] SeaState WaveTp logic --- modules/seastate/src/SeaState_Input.f90 | 13 +++++++++---- 1 file changed, 9 insertions(+), 4 deletions(-) diff --git a/modules/seastate/src/SeaState_Input.f90 b/modules/seastate/src/SeaState_Input.f90 index 1bc7fc80b8..5d7186363c 100644 --- a/modules/seastate/src/SeaState_Input.f90 +++ b/modules/seastate/src/SeaState_Input.f90 @@ -680,11 +680,16 @@ subroutine SeaStateInput_ProcessInitData( InitInp, p, InputFileData, ErrStat, Er ! WaveTp - Peak spectral period. - if ( InputFileData%Waves%WaveTp <= 0.0 ) then - call SetErrStat( ErrID_Fatal,'WaveTp must be greater than zero.',ErrStat,ErrMsg,RoutineName) - return - end if + if ( InputFileData%WaveMod == WaveMod_Regular .OR. & + InputFileData%WaveMod == WaveMod_RegularUsrPh .OR. & + InputFileData%WaveMod == WaveMod_JONSWAP ) then + if ( InputFileData%Waves%WaveTp <= 0.0 ) then + call SetErrStat( ErrID_Fatal,'WaveTp must be greater than zero.',ErrStat,ErrMsg,RoutineName) + return + end if + + end if ! WavePkShp - Peak shape parameter From bb92bb5ec1e35ac4550bafed96373a18ee3aa5c2 Mon Sep 17 00:00:00 2001 From: Roger Bergua Date: Tue, 3 Mar 2026 17:48:21 -0700 Subject: [PATCH 13/27] Small change WaveTp logic --- modules/seastate/src/SeaState_Input.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/modules/seastate/src/SeaState_Input.f90 b/modules/seastate/src/SeaState_Input.f90 index 5d7186363c..ded517105d 100644 --- a/modules/seastate/src/SeaState_Input.f90 +++ b/modules/seastate/src/SeaState_Input.f90 @@ -684,7 +684,7 @@ subroutine SeaStateInput_ProcessInitData( InitInp, p, InputFileData, ErrStat, Er InputFileData%WaveMod == WaveMod_RegularUsrPh .OR. & InputFileData%WaveMod == WaveMod_JONSWAP ) then - if ( InputFileData%Waves%WaveTp <= 0.0 ) then + if ( InputFileData%Waves%WaveTp <= 0.0_SiKi ) then call SetErrStat( ErrID_Fatal,'WaveTp must be greater than zero.',ErrStat,ErrMsg,RoutineName) return end if From 2941a282eea0e2f386d7dd45c6e11f70ed0d1b61 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Thu, 5 Mar 2026 13:06:10 -0700 Subject: [PATCH 14/27] SS c-bind: add SeaSt_C_GetFluidDynP routine to interface --- modules/seastate/src/SeaSt_WaveField.f90 | 98 +++++++++++++++++++++ modules/seastate/src/SeaState_C_Binding.f90 | 68 ++++++++++++++ 2 files changed, 166 insertions(+) diff --git a/modules/seastate/src/SeaSt_WaveField.f90 b/modules/seastate/src/SeaSt_WaveField.f90 index 04d66270d8..392576750a 100644 --- a/modules/seastate/src/SeaSt_WaveField.f90 +++ b/modules/seastate/src/SeaSt_WaveField.f90 @@ -17,6 +17,7 @@ MODULE SeaSt_WaveField PUBLIC WaveField_GetNodeWaveKin PUBLIC WaveField_GetNodeWaveVelAcc PUBLIC WaveField_GetWaveKin +PUBLIC WaveField_GetDynP CONTAINS @@ -310,6 +311,103 @@ logical function Failed() END SUBROUTINE WaveField_GetNodeWaveKin +!-------------------- Subroutine for dynamic pressure --------------------! +!NOTE: this is a subset of WaveField_GetNodeWaveKin +SUBROUTINE WaveField_GetDynP( WaveField, WaveField_m, Time, pos, forceNodeInWater, nodeInWater, FDynP, ErrStat, ErrMsg ) + type(SeaSt_WaveFieldType), intent(in ) :: WaveField + type(GridInterp_MiscVarType), intent(inout) :: WaveField_m + real(DbKi), intent(in ) :: Time + real(ReKi), intent(in ) :: pos(3) + logical, intent(in ) :: forceNodeInWater + real(SiKi), intent( out) :: FDynP + integer(IntKi), intent( out) :: nodeInWater + integer(IntKi), intent( out) :: ErrStat ! Error status of the operation + character(*), intent( out) :: ErrMsg ! Error message if errStat /= ErrID_None + + real(ReKi) :: posXY(2), posPrime(3), posXY0(3) + character(*), parameter :: RoutineName = 'WaveField_GetDynP' + integer(IntKi) :: errStat2 + character(ErrMsgLen) :: errMsg2 + + ! Temporary vars not kept + real(SiKi) :: WaveElev1 + real(SiKi) :: WaveElev2 + real(SiKi) :: WaveElev + + ErrStat = ErrID_None + ErrMsg = "" + + posXY = pos(1:2) + posXY0 = (/pos(1),pos(2),0.0_ReKi/) + + ! Wave elevation + WaveElev1 = WaveField_GetNodeWaveElev1( WaveField, WaveField_m, Time, pos, ErrStat2, ErrMsg2 ); if (Failed()) return; + WaveElev2 = WaveField_GetNodeWaveElev2( WaveField, WaveField_m, Time, pos, ErrStat2, ErrMsg2 ); if (Failed()) return; + WaveElev = WaveElev1 + WaveElev2 + + IF (WaveField%WaveStMod == 0) THEN ! No wave stretching + + IF ( pos(3) <= 0.0_ReKi) THEN ! Node is at or below the SWL + nodeInWater = 1_IntKi + ! Use location to obtain interpolated values of kinematics + CALL WaveField_Interp_Setup4D( Time+WaveField%WaveTimeShift, pos, WaveField%GridDepth, WaveField%VolGridParams, WaveField_m, ErrStat2, ErrMsg2 ); if (Failed()) return; + FDynP = GridInterp4D ( WaveField%WaveDynP, WaveField_m ) + ELSE ! Node is above the SWL + nodeInWater = 0_IntKi + FDynP = 0.0 + END IF + + ELSE ! Wave stretching enabled + IF ( (pos(3) <= WaveElev) .OR. forceNodeInWater ) THEN ! Node is submerged + nodeInWater = 1_IntKi + IF ( WaveField%WaveStMod < 3 ) THEN ! Vertical or extrapolated wave stretching + + IF ( pos(3) <= 0.0_SiKi) THEN ! Node is below the SWL - evaluate wave dynamics as usual + ! Use location to obtain interpolated values of kinematics + CALL WaveField_Interp_Setup4D( Time+WaveField%WaveTimeShift, pos, WaveField%GridDepth, WaveField%VolGridParams, WaveField_m, ErrStat2, ErrMsg2 ); if (Failed()) return; + FDynP = GridInterp4D ( WaveField%WaveDynP, WaveField_m ) + ELSE ! Node is above SWL - need wave stretching + + ! Vertical wave stretching + CALL WaveField_Interp_Setup4D( Time+WaveField%WaveTimeShift, posXY0, WaveField%GridDepth, WaveField%VolGridParams, WaveField_m, ErrStat2, ErrMsg2 ); if (Failed()) return; + FDynP = GridInterp4D ( WaveField%WaveDynP, WaveField_m ) + + ! Extrapoled wave stretching + IF (WaveField%WaveStMod == 2) THEN + CALL WaveField_Interp_Setup3D( Time+WaveField%WaveTimeShift, posXY, WaveField%SrfGridParams, WaveField_m, ErrStat2, ErrMsg2 ); if (Failed()) return; + FDynP = FDynP + GridInterp3D ( WaveField%PWaveDynP0, WaveField_m ) * pos(3) + END IF + + END IF ! Node is submerged + + ELSE ! Wheeler stretching - no need to check whether the node is above or below SWL + + ! Map the node z-position linearly from [-EffWtrDpth,m%WaveElev(j)] to [-EffWtrDpth,0] + posPrime = pos + posPrime(3) = WaveField%EffWtrDpth*(WaveField%EffWtrDpth+pos(3))/(WaveField%EffWtrDpth+WaveElev)-WaveField%EffWtrDpth + posPrime(3) = MIN( posPrime(3), 0.0_ReKi) ! Clamp z-position to zero. Needed when forceNodeInWater=.TRUE. + + ! Obtain the wave-field variables by interpolation with the mapped position. + CALL WaveField_Interp_Setup4D( Time+WaveField%WaveTimeShift, posPrime, WaveField%GridDepth, WaveField%VolGridParams, WaveField_m, ErrStat2, ErrMsg2 ); if (Failed()) return; + FDynP = GridInterp4D ( WaveField%WaveDynP, WaveField_m ) + END IF + + ELSE ! Node is out of water - zero-out all wave dynamics + + nodeInWater = 0_IntKi + FDynP = 0.0 + + END IF ! If node is in or out of water + END IF ! If wave stretching is on or off + +contains + logical function Failed() + call SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName ) + Failed = ErrStat >= AbortErrLev + end function +END SUBROUTINE WaveField_GetDynP + + !-------------------- Subroutine for wave field velocity only --------------------! SUBROUTINE WaveField_GetNodeWaveVelAcc( WaveField, WaveField_m, Time, pos, forceNodeInWater, nodeInWater, FV, FA, ErrStat, ErrMsg ) type(SeaSt_WaveFieldType), intent(in ) :: WaveField diff --git a/modules/seastate/src/SeaState_C_Binding.f90 b/modules/seastate/src/SeaState_C_Binding.f90 index 1a7a6d81d9..47298615c0 100644 --- a/modules/seastate/src/SeaState_C_Binding.f90 +++ b/modules/seastate/src/SeaState_C_Binding.f90 @@ -44,6 +44,7 @@ MODULE SeaState_C_Binding PUBLIC :: SeaSt_C_GetDens PUBLIC :: SeaSt_C_GetDpth PUBLIC :: SeaSt_C_GetMSL2SWL + PUBLIC :: SeaSt_C_GetDynPressure !------------------------------------------------------------------------------------ ! Debugging: DebugLevel -- passed at PreInit @@ -898,6 +899,73 @@ subroutine CheckWaveFieldPtr(callingRoutine,valid,ErrStat,ErrMsg) end subroutine +!> return the surface elevation and normal vector at a point. +subroutine SeaSt_C_GetDynPressure(Time_C, Pos_C, DynP_C, ErrStat_C,ErrMsg_C) BIND (C, NAME='SeaSt_C_GetDynPressure') +#ifndef IMPLICIT_DLLEXPORT +!DEC$ ATTRIBUTES DLLEXPORT :: SeaSt_C_GetDynPressure +!GCC$ ATTRIBUTES DLLEXPORT :: SeaSt_C_GetDynPressure +#endif + real(c_double), intent(in ) :: Time_C + real(c_float), intent(in ) :: Pos_c(3) + real(c_float), intent( out) :: DynP_C + integer(c_int), intent( out) :: ErrStat_C + character(kind=c_char), intent( out) :: ErrMsg_C(ErrMsgLen_C) + + real(DbKi) :: Time + real(ReKi) :: Pos(3) + real(SiKi) :: FDynP + integer :: ErrStat + character(ErrMsgLen) :: ErrMsg + character(*), parameter :: RoutineName = 'SeaSt_C_GetDynPressure' + logical :: valid + ! Temporary vars for calling GetNodeWaveKin -- remove if new interface developed + integer(IntKi) :: nodeInWater + + if (DebugLevel > 1) call ShowPassedData() + + ! convert position and time to fortran types + Time = real(Time_C, DbKi) + Pos = real(Pos_C(1:3), ReKi) + + ! verify there is actually wavefield data + call CheckWaveFieldPtr(RoutineName, valid, ErrStat, ErrMsg) + if (.not. valid) then + call Cleanup() + return + endif + + ! get wave elevation (total combined first and second order) + ! Notes: + ! - if position is outside the wave field boundary, it will simply return boundary edge value + ! - time must be positive or a fatal error occurs + ! - forceNodeInWater = .false. -- perhaps update this later to allow throug interface? + call WaveField_GetDynP( p%WaveField, m%WaveField_m, Time, pos, .false., nodeInWater, FDynP, ErrStat, ErrMsg ) + + + ! Store resulting elevation as C type + DynP_C = real(FDynP,c_float) + + if (DebugLevel > 1) call ShowReturnData() + call Cleanup() + return +contains + subroutine Cleanup() ! NOTE: we are ignoring any error reporting from here + CALL SetErrStat_F2C(ErrStat,ErrMsg,ErrStat_C,ErrMsg_C) + end subroutine Cleanup + subroutine ShowPassedData() + call WrScr("-----------------------------------------------------------") + call WrScr("Interface debugging: SeaSt_C_GetDynPressure") + call WrScr(" --------------------------------------------------------") + call WrScr(" Time_C -> "//trim(Num2LStr(Time_C))) + call WrScr(" Pos_C -> ("//trim(Num2LStr(Pos_C(1)))//","//trim(Num2LStr(Pos_C(2)))//")") + end subroutine ShowPassedData + subroutine ShowReturnData() + call WrScr(" DynP_C <- "//trim(Num2LStr(DynP_C))) + call WrScr("-----------------------------------------------------------") + end subroutine ShowReturnData +end subroutine SeaSt_C_GetDynPressure + + !FIXME: the following visualization writer should be merged into the library vtk.f90 ! this is a modified duplicate of the routine from FAST_Subs by the same name. !FIXME: this routine currently only grabs the closest timestep and does not do any interpolation From 4c603ad276712a539598d14405f005f2e697f490 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Thu, 5 Mar 2026 14:43:16 -0700 Subject: [PATCH 15/27] SS c-bind: improve debug message on GetDynP --- modules/seastate/src/SeaState_C_Binding.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/modules/seastate/src/SeaState_C_Binding.f90 b/modules/seastate/src/SeaState_C_Binding.f90 index 47298615c0..747b786d7d 100644 --- a/modules/seastate/src/SeaState_C_Binding.f90 +++ b/modules/seastate/src/SeaState_C_Binding.f90 @@ -957,7 +957,7 @@ subroutine ShowPassedData() call WrScr("Interface debugging: SeaSt_C_GetDynPressure") call WrScr(" --------------------------------------------------------") call WrScr(" Time_C -> "//trim(Num2LStr(Time_C))) - call WrScr(" Pos_C -> ("//trim(Num2LStr(Pos_C(1)))//","//trim(Num2LStr(Pos_C(2)))//")") + call WrScr(" Pos_C -> ("//trim(Num2LStr(Pos_C(1)))//","//trim(Num2LStr(Pos_C(2)))//","//trim(Num2LStr(Pos_C(3)))//")") end subroutine ShowPassedData subroutine ShowReturnData() call WrScr(" DynP_C <- "//trim(Num2LStr(DynP_C))) From f7454209daa0d6ae04f76391fa02920e043e74f6 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Mon, 9 Mar 2026 13:16:41 -0600 Subject: [PATCH 16/27] AWAE: restructure UpdateStates to match 5.0 before introducing changes in data structure --- modules/awae/src/AWAE.f90 | 180 ++++++++++++++++++++++++++------------ 1 file changed, 126 insertions(+), 54 deletions(-) diff --git a/modules/awae/src/AWAE.f90 b/modules/awae/src/AWAE.f90 index d2a77cd576..53df7dd679 100644 --- a/modules/awae/src/AWAE.f90 +++ b/modules/awae/src/AWAE.f90 @@ -1388,34 +1388,75 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM !! Output: Constraint states at t+dt type(AWAE_OtherStateType), intent(inout) :: OtherState !< Input: Other states at t; !! Output: Other states at t+dt - type(AWAE_MiscVarType), intent(inout) :: m !< Misc/optimization variables + type(AWAE_MiscVarType), target, intent(inout) :: m !< Misc/optimization variables integer(IntKi), intent( out) :: errStat !< Error status of the operation character(*), intent( out) :: errMsg !< Error message if errStat /= ErrID_None - type(AWAE_InputType) :: uInterp ! Interpolated/Extrapolated input + character(*), parameter :: RoutineName = 'AWAE_UpdateStates' integer(intKi) :: errStat2 ! temporary Error status character(ErrMsgLen) :: errMsg2 ! temporary Error message - character(*), parameter :: RoutineName = 'AWAE_UpdateStates' - integer(IntKi) :: n_high_low, nt, i_hl, i,j,k,c + type(AWAE_InputType) :: uInterp ! Interpolated/Extrapolated input + integer(IntKi) :: n_high_low, nt, i_hl + integer(IntKi) :: i,j,k,c + real(ReKi), pointer :: V_Grid(:,:,:,:) errStat = ErrID_None errMsg = "" - - ! Read the ambient wind data that is needed for t+dt, i.e., n+1 - + + ! If last time step, don't populate high-resolution grid if ( (n+1) == (p%NumDT-1) ) then n_high_low = 0 else n_high_low = p%n_high_low end if - if ( p%Mod_AmbWind == 1 ) then - ! read from file the ambient flow for the n+1 time step + !---------------------------------------------------------------------------- + ! Populate low resolution grids based on ambient wind source + !---------------------------------------------------------------------------- + + select case (p%Mod_AmbWind) + + ! File-based ambient wind + case (1) + + ! Read from file the ambient flow for the n+1 time step call ReadLowResWindFile(n+1, p, m%Vamb_Low, errStat2, errMsg2); if (Failed()) return; - !$OMP PARALLEL DO DEFAULT(Shared) PRIVATE(nt, i_hl, errStat2, errMsg2) !Private(nt,tm2,tm3) + ! InflowWind-based ambient wind (single or multiple instances) + case (2, 3) + +!FIXME: remove next 3 lines in 5.0.0 + ! Set the hub position and orientation to pass to IfW (IfW always calculates hub and disk avg vel) -- note that this is garbage data. + m%u_IfW_Low%HubPosition = (/ p%X0_low + 0.5*p%nX_low*p%dX_low, p%Y0_low + 0.5*p%nY_low*p%dY_low, p%Z0_low + 0.5*p%nZ_low*p%dZ_low /) + call Eye(m%u_IfW_Low%HubOrientation,ErrStat2,ErrMsg2); if (Failed()) return; + + ! Calculate the low-resolution grid inflow velocities + call InflowWind_CalcOutput(t+p%dt_low, m%u_IfW_Low, p%IfW(0), x%IfW(0), xd%IfW(0), z%IfW(0), OtherState%IfW(0), m%y_IfW_Low, m%IfW(0), errStat2, errMsg2) + if (Failed()) return + + ! Transfer velocities to low resolution grid + V_Grid(lbound(m%Vamb_low,1):ubound(m%Vamb_low,1),& + lbound(m%Vamb_low,2):ubound(m%Vamb_low,2),& + lbound(m%Vamb_low,3):ubound(m%Vamb_low,3),& + lbound(m%Vamb_low,4):ubound(m%Vamb_low,4)) => m%y_IfW_Low%VelocityUVW + m%Vamb_Low = V_Grid + + end select + + !---------------------------------------------------------------------------- + ! Populate high-resolution grid based on ambient wind source + !---------------------------------------------------------------------------- + + select case (p%Mod_AmbWind) + + ! File-based ambient wind + case (1) + + !$OMP PARALLEL DO DEFAULT(NONE) PRIVATE(nt, i_hl, errStat2, errMsg2) & + !$OMP SHARED(p, n_high_low, n, m, errStat, errMsg, AbortErrLev) do nt = 1,p%NumTurbines do i_hl=0, n_high_low + ! read from file the ambient flow for the current time step call ReadHighResWindFile(nt, (n+1)*p%n_high_low + i_hl, p, m%Vamb_high(nt)%data(:,:,:,:,i_hl), errStat2, errMsg2) if (ErrStat2 >= AbortErrLev) then @@ -1426,33 +1467,33 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM end do end do !$OMP END PARALLEL DO + if (errStat >= AbortErrLev) return - else ! p%Mod_AmbWind == 2 .or. 3 + ! Single InflowWind instance + case (2) + ! Loop through turbines + do nt = 1, p%NumTurbines + +!FIXME: remove next 4 lines in merge with 5.0 ! Set the hub position and orientation to pass to IfW (IfW always calculates hub and disk avg vel) - m%u_IfW_Low%HubPosition = (/ p%X0_low + 0.5*p%nX_low*p%dX_low, p%Y0_low + 0.5*p%nY_low*p%dY_low, p%Z0_low + 0.5*p%nZ_low*p%dZ_low /) - call Eye(m%u_IfW_Low%HubOrientation,ErrStat2,ErrMsg2); if (Failed()) return; - ! Set low-resolution inflow wind velocities - call InflowWind_CalcOutput(t+p%dt_low, m%u_IfW_Low, p%IfW(0), x%IfW(0), xd%IfW(0), z%IfW(0), OtherState%IfW(0), m%y_IfW_Low, m%IfW(0), errStat2, errMsg2); if (Failed()) return; - c = 1 - do k = 0,p%nZ_low-1 - do j = 0,p%nY_low-1 - do i = 0,p%nX_low-1 - m%Vamb_Low(:,i,j,k) = m%y_IfW_Low%VelocityUVW(:,c) - c = c+1 - end do - end do - end do - ! Set the high-resoultion inflow wind velocities for each turbine - if ( p%Mod_AmbWind == 2 ) then - do nt = 1,p%NumTurbines - m%u_IfW_High%PositionXYZ = p%Grid_high(:,:,nt) - ! Set the hub position and orientation to pass to IfW (IfW always calculates hub and disk avg vel) - m%u_IfW_High%HubPosition = (/ p%X0_high(nt) + 0.5*p%nX_high*p%dX_high(nt), p%Y0_high(nt) + 0.5*p%nY_high*p%dY_high(nt), p%Z0_high(nt) + 0.5*p%nZ_high*p%dZ_high(nt) /) - call Eye(m%u_IfW_High%HubOrientation,ErrStat2,ErrMsg2) - do i_hl=0, n_high_low - call InflowWind_CalcOutput(t+p%dt_low+i_hl*p%DT_high, m%u_IfW_High, p%IfW(0), x%IfW(0), xd%IfW(0), z%IfW(0), OtherState%IfW(0), m%y_IfW_High, m%IfW(0), errStat2, errMsg2); if (Failed()) return; + m%u_IfW_High%HubPosition = (/ p%X0_high(nt) + 0.5*p%nX_high*p%dX_high(nt), p%Y0_high(nt) + 0.5*p%nY_high*p%dY_high(nt), p%Z0_high(nt) + 0.5*p%nZ_high*p%dZ_high(nt) /) + call Eye(m%u_IfW_High%HubOrientation,ErrStat2,ErrMsg2) + +!FIXME: remove next 3 lines in merge with 5.0 + ! Set input position + m%u_IfW_High%PositionXYZ = p%Grid_high(:,:,nt) + + ! Loop through high resolution grids + do i_hl = 0, n_high_low + + ! Calculate wind velocities at grid locations from InflowWind + call InflowWind_CalcOutput(t+p%dt_low+i_hl*p%DT_high, m%u_IfW_High, p%IfW(0), x%IfW(0), xd%IfW(0), z%IfW(0), OtherState%IfW(0), m%y_IfW_High, m%IfW(0), errStat2, errMsg2) + if (Failed()) return + + ! Transfer velocities to high resolution grid +!FIXME: remove following 9 lines and uncomment block after during merge to 5.0 c = 1 do k = 0,p%nZ_high-1 do j = 0,p%nY_high-1 @@ -1462,25 +1503,46 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM end do end do end do - end do +! V_Grid(lbound(m%Vamb_high(nt)%data,1):ubound(m%Vamb_high(nt)%data,1),& +! lbound(m%Vamb_high(nt)%data,2):ubound(m%Vamb_high(nt)%data,2),& +! lbound(m%Vamb_high(nt)%data,3):ubound(m%Vamb_high(nt)%data,3),& +! lbound(m%Vamb_high(nt)%data,4):ubound(m%Vamb_high(nt)%data,4)) => m%y_IfW_High(nt)%VelocityUVW +! m%Vamb_high(nt)%data(:,:,:,:,i_hl) = V_Grid end do + end do - else !p%Mod_AmbWind == 3 - do nt = 1,p%NumTurbines - c = 1 - do k = 0,p%nZ_high-1 - do j = 0,p%nY_high-1 - do i = 0,p%nX_high-1 - m%u_IfW_High%PositionXYZ(:,c) = p%Grid_high(:,c,nt) - p%WT_Position(:,nt) - c = c+1 - end do + ! Multiple InflowWind instances (one per turbine) + case (3) + + ! Loop through turbines + do nt = 1, p%NumTurbines + +!FIXME: remove next 10 lines in merge with 5.0 + ! Set input velocity + c = 1 + do k = 0,p%nZ_high-1 + do j = 0,p%nY_high-1 + do i = 0,p%nX_high-1 + m%u_IfW_High%PositionXYZ(:,c) = p%Grid_high(:,c,nt) - p%WT_Position(:,nt) + c = c+1 end do end do - do i_hl=0, n_high_low - ! Set the hub position and orientation to pass to IfW (IfW always calculates hub and disk avg vel) - m%u_IfW_High%HubPosition = (/ p%X0_high(nt) + 0.5*p%nX_high*p%dX_high(nt), p%Y0_high(nt) + 0.5*p%nY_high*p%dY_high(nt), p%Z0_high(nt) + 0.5*p%nZ_high*p%dZ_high(nt) /) - p%WT_Position(:,nt) - call Eye(m%u_IfW_High%HubOrientation,ErrStat2,ErrMsg2) - call InflowWind_CalcOutput(t+p%dt_low+i_hl*p%DT_high, m%u_IfW_High, p%IfW(nt), x%IfW(nt), xd%IfW(nt), z%IfW(nt), OtherState%IfW(nt), m%y_IfW_High, m%IfW(nt), errStat2, errMsg2); if (Failed()) return; + end do + + ! Loop through high resolution grids + do i_hl = 0, n_high_low + +!FIXME: remove next 4 lines in merge with 5.0 + ! Set the hub position and orientation to pass to IfW (IfW always calculates hub and disk avg vel) + m%u_IfW_High%HubPosition = (/ p%X0_high(nt) + 0.5*p%nX_high*p%dX_high(nt), p%Y0_high(nt) + 0.5*p%nY_high*p%dY_high(nt), p%Z0_high(nt) + 0.5*p%nZ_high*p%dZ_high(nt) /) - p%WT_Position(:,nt) + call Eye(m%u_IfW_High%HubOrientation,ErrStat2,ErrMsg2) + + ! Calculate wind velocities at grid locations from InflowWind + call InflowWind_CalcOutput(t+p%dt_low+i_hl*p%DT_high, m%u_IfW_High, p%IfW(nt), x%IfW(nt), xd%IfW(nt), z%IfW(nt), OtherState%IfW(nt), m%y_IfW_High, m%IfW(nt), errStat2, errMsg2) + if (Failed()) return + + ! Transfer velocities to high resolution grid +!FIXME: remove following 9 lines and uncomment block after during merge to 5.0 c = 1 do k = 0,p%nZ_high-1 do j = 0,p%nY_high-1 @@ -1490,26 +1552,36 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM end do end do end do - end do +! V_Grid(lbound(m%Vamb_high(nt)%data,1):ubound(m%Vamb_high(nt)%data,1),& +! lbound(m%Vamb_high(nt)%data,2):ubound(m%Vamb_high(nt)%data,2),& +! lbound(m%Vamb_high(nt)%data,3):ubound(m%Vamb_high(nt)%data,3),& +! lbound(m%Vamb_high(nt)%data,4):ubound(m%Vamb_high(nt)%data,4)) => m%y_IfW_High(nt)%VelocityUVW +! m%Vamb_high(nt)%data(:,:,:,:,i_hl) = V_Grid end do - end if - end if + end do + + end select + + !---------------------------------------------------------------------------- + ! Propagate WAT tracer + !---------------------------------------------------------------------------- - ! WAT tracer propagation if (p%WAT_Enabled) then - ! find mean velocity of all turbine disks + + ! Find mean velocity of all turbine disks xd%Ufarm = 0.0_ReKi do nt=1,p%NumTurbines xd%Ufarm(1:3) = xd%Ufarm(1:3) + m%V_amb_low_disk(1:3,nt) enddo xd%Ufarm(1:3) = xd%Ufarm(1:3) / real(p%NumTurbines,ReKi) + ! add mean velocity * dt to the tracer for the position of the WAT box xd%WAT_B_Box(1:3) = xd%WAT_B_Box(1:3) + xd%Ufarm(1:3)*real(p%dt_low,ReKi) endif contains logical function Failed() - call SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName) + call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName) Failed = ErrStat >= AbortErrLev end function Failed end subroutine AWAE_UpdateStates From 6707322ecf454f6f62eb77c6cdafad332875b155 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Mon, 9 Mar 2026 13:38:21 -0600 Subject: [PATCH 17/27] AWAE.f90: further updatestates changes to more closely match 5.0 --- modules/awae/src/AWAE.f90 | 56 +++++++++++++++++++-------------------- 1 file changed, 28 insertions(+), 28 deletions(-) diff --git a/modules/awae/src/AWAE.f90 b/modules/awae/src/AWAE.f90 index 53df7dd679..0602479f5b 100644 --- a/modules/awae/src/AWAE.f90 +++ b/modules/awae/src/AWAE.f90 @@ -1425,7 +1425,7 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM ! InflowWind-based ambient wind (single or multiple instances) case (2, 3) -!FIXME: remove next 3 lines in 5.0.0 +!FIXME:merge5.0 remove next 3 lines ! Set the hub position and orientation to pass to IfW (IfW always calculates hub and disk avg vel) -- note that this is garbage data. m%u_IfW_Low%HubPosition = (/ p%X0_low + 0.5*p%nX_low*p%dX_low, p%Y0_low + 0.5*p%nY_low*p%dY_low, p%Z0_low + 0.5*p%nZ_low*p%dZ_low /) call Eye(m%u_IfW_Low%HubOrientation,ErrStat2,ErrMsg2); if (Failed()) return; @@ -1457,8 +1457,10 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM do nt = 1,p%NumTurbines do i_hl=0, n_high_low - ! read from file the ambient flow for the current time step + ! read from file the ambient flow for the current time step +!FIXME:merge5.0 replace next line with the following call ReadHighResWindFile(nt, (n+1)*p%n_high_low + i_hl, p, m%Vamb_high(nt)%data(:,:,:,:,i_hl), errStat2, errMsg2) +! call ReadHighResWindVTK(nt, n*p%n_high_low + i_hl, p, m%Vamb_high(nt)%data(:,:,:,:,i_hl), errStat2, errMsg2) if (ErrStat2 >= AbortErrLev) then !$OMP CRITICAL ! Needed to avoid data race on ErrStat and ErrMsg call SetErrStat( ErrStat2, ErrMsg2, errStat, errMsg, RoutineName ) @@ -1476,12 +1478,12 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM ! Loop through turbines do nt = 1, p%NumTurbines -!FIXME: remove next 4 lines in merge with 5.0 +!FIXME:merge5.0 remove next 4 lines ! Set the hub position and orientation to pass to IfW (IfW always calculates hub and disk avg vel) m%u_IfW_High%HubPosition = (/ p%X0_high(nt) + 0.5*p%nX_high*p%dX_high(nt), p%Y0_high(nt) + 0.5*p%nY_high*p%dY_high(nt), p%Z0_high(nt) + 0.5*p%nZ_high*p%dZ_high(nt) /) call Eye(m%u_IfW_High%HubOrientation,ErrStat2,ErrMsg2) -!FIXME: remove next 3 lines in merge with 5.0 +!FIXME:merge5.0 remove next 3 lines ! Set input position m%u_IfW_High%PositionXYZ = p%Grid_high(:,:,nt) @@ -1489,25 +1491,24 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM do i_hl = 0, n_high_low ! Calculate wind velocities at grid locations from InflowWind +!FIXME:merge5.0 replace next line with the next call InflowWind_CalcOutput(t+p%dt_low+i_hl*p%DT_high, m%u_IfW_High, p%IfW(0), x%IfW(0), xd%IfW(0), z%IfW(0), OtherState%IfW(0), m%y_IfW_High, m%IfW(0), errStat2, errMsg2) +! call IfW_FlowField_GetVelAcc(p%IfW(0)%FlowField, 1, t + i_hl*p%DT_high, & +! m%u_IfW_High(nt)%PositionXYZ, & +! m%y_IfW_High(nt)%VelocityUVW, AccUVW, errStat2, errMsg2) if (Failed()) return ! Transfer velocities to high resolution grid -!FIXME: remove following 9 lines and uncomment block after during merge to 5.0 - c = 1 - do k = 0,p%nZ_high-1 - do j = 0,p%nY_high-1 - do i = 0,p%nX_high-1 - m%Vamb_high(nt)%data(:,i,j,k,i_hl) = m%y_IfW_High%VelocityUVW(:,c) - c = c+1 - end do - end do - end do +!FIXME:merge5.0 remove following 4 lines and uncomment block after + V_Grid(lbound(m%Vamb_high(nt)%data,1):ubound(m%Vamb_high(nt)%data,1),& + lbound(m%Vamb_high(nt)%data,2):ubound(m%Vamb_high(nt)%data,2),& + lbound(m%Vamb_high(nt)%data,3):ubound(m%Vamb_high(nt)%data,3),& + lbound(m%Vamb_high(nt)%data,4):ubound(m%Vamb_high(nt)%data,4)) => m%y_IfW_High%VelocityUVW ! V_Grid(lbound(m%Vamb_high(nt)%data,1):ubound(m%Vamb_high(nt)%data,1),& ! lbound(m%Vamb_high(nt)%data,2):ubound(m%Vamb_high(nt)%data,2),& ! lbound(m%Vamb_high(nt)%data,3):ubound(m%Vamb_high(nt)%data,3),& ! lbound(m%Vamb_high(nt)%data,4):ubound(m%Vamb_high(nt)%data,4)) => m%y_IfW_High(nt)%VelocityUVW -! m%Vamb_high(nt)%data(:,:,:,:,i_hl) = V_Grid + m%Vamb_high(nt)%data(:,:,:,:,i_hl) = V_Grid end do end do @@ -1517,7 +1518,7 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM ! Loop through turbines do nt = 1, p%NumTurbines -!FIXME: remove next 10 lines in merge with 5.0 +!FIXME:merge5.0 remove next 10 lines ! Set input velocity c = 1 do k = 0,p%nZ_high-1 @@ -1532,31 +1533,30 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM ! Loop through high resolution grids do i_hl = 0, n_high_low -!FIXME: remove next 4 lines in merge with 5.0 +!FIXME:merge5.0 remove next 4 lines ! Set the hub position and orientation to pass to IfW (IfW always calculates hub and disk avg vel) m%u_IfW_High%HubPosition = (/ p%X0_high(nt) + 0.5*p%nX_high*p%dX_high(nt), p%Y0_high(nt) + 0.5*p%nY_high*p%dY_high(nt), p%Z0_high(nt) + 0.5*p%nZ_high*p%dZ_high(nt) /) - p%WT_Position(:,nt) call Eye(m%u_IfW_High%HubOrientation,ErrStat2,ErrMsg2) ! Calculate wind velocities at grid locations from InflowWind +!FIXME:merge5.0 replace next line with the next call InflowWind_CalcOutput(t+p%dt_low+i_hl*p%DT_high, m%u_IfW_High, p%IfW(nt), x%IfW(nt), xd%IfW(nt), z%IfW(nt), OtherState%IfW(nt), m%y_IfW_High, m%IfW(nt), errStat2, errMsg2) +! call IfW_FlowField_GetVelAcc(p%IfW(nt)%FlowField, 1, t + i_hl*p%DT_high, & +! m%u_IfW_High(nt)%PositionXYZ, & +! m%y_IfW_High(nt)%VelocityUVW, AccUVW, errStat2, errMsg2) if (Failed()) return ! Transfer velocities to high resolution grid -!FIXME: remove following 9 lines and uncomment block after during merge to 5.0 - c = 1 - do k = 0,p%nZ_high-1 - do j = 0,p%nY_high-1 - do i = 0,p%nX_high-1 - m%Vamb_high(nt)%data(:,i,j,k,i_hl) = m%y_IfW_High%VelocityUVW(:,c) - c = c+1 - end do - end do - end do +!FIXME:merge5.0 remove following 4 lines and uncomment block after + V_Grid(lbound(m%Vamb_high(nt)%data,1):ubound(m%Vamb_high(nt)%data,1),& + lbound(m%Vamb_high(nt)%data,2):ubound(m%Vamb_high(nt)%data,2),& + lbound(m%Vamb_high(nt)%data,3):ubound(m%Vamb_high(nt)%data,3),& + lbound(m%Vamb_high(nt)%data,4):ubound(m%Vamb_high(nt)%data,4)) => m%y_IfW_High%VelocityUVW ! V_Grid(lbound(m%Vamb_high(nt)%data,1):ubound(m%Vamb_high(nt)%data,1),& ! lbound(m%Vamb_high(nt)%data,2):ubound(m%Vamb_high(nt)%data,2),& ! lbound(m%Vamb_high(nt)%data,3):ubound(m%Vamb_high(nt)%data,3),& ! lbound(m%Vamb_high(nt)%data,4):ubound(m%Vamb_high(nt)%data,4)) => m%y_IfW_High(nt)%VelocityUVW -! m%Vamb_high(nt)%data(:,:,:,:,i_hl) = V_Grid + m%Vamb_high(nt)%data(:,:,:,:,i_hl) = V_Grid end do end do From 7eaeedd586f24dafe30558733bb70a62628d5154 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Mon, 9 Mar 2026 15:19:51 -0600 Subject: [PATCH 18/27] FF: add extra timestep at T=-DT_high to wind sent to OF instances --- glue-codes/fast-farm/src/FASTWrapper.f90 | 16 +++++++++------- glue-codes/fast-farm/src/FAST_Farm_Subs.f90 | 2 +- 2 files changed, 10 insertions(+), 8 deletions(-) diff --git a/glue-codes/fast-farm/src/FASTWrapper.f90 b/glue-codes/fast-farm/src/FASTWrapper.f90 index 3207a139b2..4ff6b281f2 100644 --- a/glue-codes/fast-farm/src/FASTWrapper.f90 +++ b/glue-codes/fast-farm/src/FASTWrapper.f90 @@ -44,7 +44,7 @@ MODULE FASTWrapper PUBLIC :: FWrap_t0 ! call to compute outputs at t0 [and initialize some more variables] PUBLIC :: FWrap_Increment ! call to update states to n+1 and compute outputs at n+1 - PUBLIC :: FWrap_SetInputs + PUBLIC :: FWrap_SetWindTStart PUBLIC :: FWrap_CalcOutput @@ -119,7 +119,7 @@ SUBROUTINE FWrap_Init( InitInp, u, p, x, xd, z, OtherState, y, m, Interval, Init ExternInitData%windGrid_n(1) = InitInp%nX_high ExternInitData%windGrid_n(2) = InitInp%nY_high ExternInitData%windGrid_n(3) = InitInp%nZ_high - ExternInitData%windGrid_n(4) = InitInp%n_high_low + ExternInitData%windGrid_n(4) = InitInp%n_high_low+1 ! include a step at t-dt_high ExternInitData%windGrid_delta(1) = InitInp%dX_high ExternInitData%windGrid_delta(2) = InitInp%dY_high @@ -392,7 +392,7 @@ SUBROUTINE FWrap_Increment( t, n, u, p, x, xd, z, OtherState, y, m, ErrStat, Err !ELSE ! ! set the inputs needed for FAST - !call FWrap_SetInputs(u, m, t) <<< moved up into FAST.Farm FARM_UpdateStates + !call FWrap_SetWindTStart(u, m, t) <<< moved up into FAST.Farm FARM_UpdateStates ! call FAST p%n_FAST_low times (p%n_FAST_low is simply the number of steps to make per wrapper call. It is affected by MooringMod) do n_ss = 1, p%n_FAST_low @@ -434,7 +434,7 @@ SUBROUTINE FWrap_t0( u, p, x, xd, z, OtherState, y, m, ErrStat, ErrMsg ) ErrMsg = '' ! set the inputs needed for FAST: - call FWrap_SetInputs(u, m, 0.0_DbKi) + call FWrap_SetWindTStart(u, m, 0.0_DbKi) ! compute the FAST t0 solution: call FAST_Solution0_T(m%Turbine, ErrStat2, ErrMsg2 ) @@ -685,16 +685,18 @@ SUBROUTINE FWrap_CalcOutput(p, u, y, m, ErrStat, ErrMsg) END SUBROUTINE FWrap_CalcOutput !---------------------------------------------------------------------------------------------------------------------------------- !> This subroutine sets the inputs needed before calling an instance of FAST -SUBROUTINE FWrap_SetInputs(u, m, t) +SUBROUTINE FWrap_SetWindTStart(u, m, t) TYPE(FWrap_InputType), INTENT(INOUT) :: u !< Inputs at t TYPE(FWrap_MiscVarType), INTENT(INOUT) :: m !< Misc variables for optimization (not copied in glue code) REAL(DbKi), INTENT(IN ) :: t !< current simulation time ! set the 4d-wind-inflow input array (a bit of a hack [simplification] so that we don't have large amounts of data copied in multiple data structures): - m%Turbine%IfW%p%FlowField%Grid4D%TimeStart = t + ! NOTE: the wind data starts at `t - DT_high` as one extra slice of wind data is added at start. If AeroDyn is updated to not require the `t-DT_high` + ! timestep, this can be changed + m%Turbine%IfW%p%FlowField%Grid4D%TimeStart = t - m%Turbine%IfW%p%FlowField%Grid4D%delta(4) -END SUBROUTINE FWrap_SetInputs +END SUBROUTINE FWrap_SetWindTStart !---------------------------------------------------------------------------------------------------------------------------------- END MODULE FASTWrapper !********************************************************************************************************************************** diff --git a/glue-codes/fast-farm/src/FAST_Farm_Subs.f90 b/glue-codes/fast-farm/src/FAST_Farm_Subs.f90 index 344241616b..5cd2381d1b 100644 --- a/glue-codes/fast-farm/src/FAST_Farm_Subs.f90 +++ b/glue-codes/fast-farm/src/FAST_Farm_Subs.f90 @@ -1204,7 +1204,7 @@ subroutine FARM_UpdateStates(t, n, farm, ErrStat, ErrMsg) ! set the inputs needed for FAST (these are slow-varying so can just be done once per farm time step) do nt = 1,farm%p%NumTurbines - call FWrap_SetInputs(farm%FWrap(nt)%u, farm%FWrap(nt)%m, t) + call FWrap_SetWindTStart(farm%FWrap(nt)%u, farm%FWrap(nt)%m, t) end do From a9fd45188be593bf42003b4675178a2aee81fe44 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Mon, 9 Mar 2026 15:26:27 -0600 Subject: [PATCH 19/27] AWAE: add wind timeslice at T=-DT_high to wind sent to OF instances - Increased arrays by one timeslice. Wind timeslice index=0 is now at `T_low-DT_high` - retrieval of wind now accounts for picking up one extra slice - if low res timestep `n==0`, set `T=-DT_high` wind equal to `T=0` wind - for lowres steps `n>0`, transfer last wind slice at `T=T_low_previous-DT_high` to `T=T_low_now-DT_high` - NOTE: it might also be possible to copy the `T=T_low_previous+n_high_low*DT_high` over to `T=T_low` since they are technically the same --- modules/awae/src/AWAE.f90 | 65 +++++++++++++++++++++--------- modules/awae/src/AWAE_Registry.txt | 1 + modules/awae/src/AWAE_Types.f90 | 4 ++ 3 files changed, 51 insertions(+), 19 deletions(-) diff --git a/modules/awae/src/AWAE.f90 b/modules/awae/src/AWAE.f90 index 0602479f5b..e304934f25 100644 --- a/modules/awae/src/AWAE.f90 +++ b/modules/awae/src/AWAE.f90 @@ -751,7 +751,7 @@ subroutine HighResGridCalcOutput(n, u, p, xd, y, m, errStat, errMsg) real(ReKi) :: WAT_k ! WAT scaling factor (averaged from overlapping wakes) real(SiKi) :: WAT_V(3) ! WAT velocity contribution real(ReKi) :: Pos_global(3) ! global position - real(ReKi), allocatable :: WAT_B_BoxHi(:,:) ! position of WAT box (global) for each intermediate steps, shape: 3 x n_high_low + real(ReKi), allocatable :: WAT_B_BoxHi(:,:) ! position of WAT box (global) for each intermediate steps, shape: 3 x n_high_low_p1 real(ReKi), allocatable :: wk_R_p2i(:,:,:)!< Orientations from plane to inertial for each wake, shape: 3x3xnWake real(ReKi), allocatable :: wk_V(:,:) !< Wake velocity from each overlapping wake, shape: 3xnWake real(ReKi), allocatable :: wk_WAT_k(:) !< WAT scaling factors for all wakes (for overlap) @@ -760,7 +760,7 @@ subroutine HighResGridCalcOutput(n, u, p, xd, y, m, errStat, errMsg) integer(IntKi) :: maxPln integer(IntKi) :: maxN_wake integer(IntKi) :: NumGrid_high !< number of points in high res grid grid - integer(IntKi) :: n_high_low + integer(IntKi) :: n_high_low_p1 integer(IntKi) :: WAT_iT,WAT_iY,WAT_iZ !< indexes for WAT point (Time interchangeable with X) integer(IntKi) :: errStat2 character(*), parameter :: RoutineName = 'HighResGridCalcOutput' @@ -770,10 +770,10 @@ subroutine HighResGridCalcOutput(n, u, p, xd, y, m, errStat, errMsg) maxPln = min(n,p%NumPlanes-2) ! We only need one high res file for that last simulation time - if ( (n/p%n_high_low) == (p%NumDT-1) ) then - n_high_low = 0 + if ( (n/p%n_high_low_p1) == (p%NumDT-1) ) then + n_high_low_p1 = 1 else - n_high_low = p%n_high_low + n_high_low_p1 = p%n_high_low_p1 end if maxN_wake = p%NumTurbines*( p%NumPlanes-1 ) @@ -786,10 +786,11 @@ subroutine HighResGridCalcOutput(n, u, p, xd, y, m, errStat, errMsg) ! Convect WAT Box tracer for each intermediate step ! Note: we substract because the high-res points are "before" current low res point if (p%WAT_Enabled) then - allocate ( WAT_B_BoxHi ( 3, 0:n_high_low), STAT=errStat2 ); if (errStat2 /= 0) call SetErrStat ( ErrID_Fatal, 'Could not allocate memory for WAT_B_BoxHi.', errStat, errMsg, RoutineName ) + allocate ( WAT_B_BoxHi ( 3, 0:n_high_low_p1), STAT=errStat2 ); if (errStat2 /= 0) call SetErrStat ( ErrID_Fatal, 'Could not allocate memory for WAT_B_BoxHi.', errStat, errMsg, RoutineName ) if (ErrStat >= AbortErrLev) return - do i_hl=0, n_high_low - WAT_B_BoxHi(1:3, i_hl) = xd%WAT_B_Box(1:3) - (n_high_low-i_hl) * xd%Ufarm(1:3) * real(p%DT_high,ReKi) + do i_hl=0, n_high_low_p1 + ! NOTE: i_hl=0 corresponds to T=-DT_high, so shift by one high res timestep + WAT_B_BoxHi(1:3, i_hl) = xd%WAT_B_Box(1:3) - (n_high_low_p1-i_hl+1) * xd%Ufarm(1:3) * real(p%DT_high,ReKi) enddo endif @@ -817,7 +818,7 @@ subroutine HighResGridCalcOutput(n, u, p, xd, y, m, errStat, errMsg) !$OMP& V_qs, & !$OMP& i_hl, Pos_global,& !$OMP& wk_WAT_k, WAT_k, WAT_iT, WAT_iY, WAT_iZ, WAT_V)& - !$OMP SHARED(NumGrid_High, m, u, p, y, xd, nt, maxPln, n_high_low, WAT_B_BoxHi, errStat, errMsg) + !$OMP SHARED(NumGrid_High, m, u, p, y, xd, nt, maxPln, n_high_low_p1, WAT_B_BoxHi, errStat, errMsg) ! Loop over all points of the high resolution ambient wind do iXYZ=1, NumGrid_high ! From flat index iXYZ to grid indices @@ -849,7 +850,7 @@ subroutine HighResGridCalcOutput(n, u, p, xd, y, m, errStat, errMsg) endif ! --- Store full velocity (Ambient + Wake QS + WAT) in grid - do i_hl=0, n_high_low + do i_hl=0, n_high_low_p1 ! Compute WAT velocity if (p%WAT_Enabled) then ! find location of grid point in the turbulent box, accounting for the convection of the box in between high res and low res @@ -929,7 +930,8 @@ subroutine AWAE_Init( InitInp, u, p, x, xd, z, OtherState, y, m, Interval, InitO p%NumRadii = InitInp%InputFileData%NumRadii p%NumTurbines = InitInp%InputFileData%NumTurbines p%WindFilePath = InitInp%InputFileData%WindFilePath ! TODO: Make sure this wasn't specified with the trailing folder separator. Note: on Windows a trailing / or \ causes no problem! GJH - p%n_high_low = InitInp%n_high_low + p%n_high_low = InitInp%n_high_low ! number of timesteps between low res steps + p%n_high_low_p1 = InitInp%n_high_low + 1 ! include a time slice at t_low-DT_high (for interpolation in AeroDyn -- this is a hack) p%NumDT = InitInp%NumDT p%NOutDisWindXY = InitInp%InputFileData%NOutDisWindXY p%NOutDisWindYZ = InitInp%InputFileData%NOutDisWindYZ @@ -991,7 +993,7 @@ subroutine AWAE_Init( InitInp, u, p, x, xd, z, OtherState, y, m, Interval, InitO ! This will establish certain parameters as well as all of the initialization outputs ! Sets: ! Parameters: nX_low, nY_low, nZ_low, nX_high, nY_high, nZ_high, Grid_low, - ! Grid_high, n_high_low, n_rp_max + ! Grid_high, n_high_low_p1, n_rp_max ! InitOutput: X0_high, Y0_high, Z0_high, dX_high, dY_high, dZ_high, nX_high, nY_high, nZ_high call AWAE_IO_InitGridInfo(InitInp, p, InitOut, errStat2, errMsg2); if(Failed()) return; @@ -1156,7 +1158,7 @@ subroutine AWAE_Init( InitInp, u, p, x, xd, z, OtherState, y, m, Interval, InitO allocate ( y%V_plane(3,0:p%NumPlanes-1,1:p%NumTurbines), STAT=ErrStat2 ); if (Failed0('y%V_plane.' )) return; allocate ( y%Vdist_High(1:p%NumTurbines), STAT=ErrStat2 ); if (Failed0('y%Vdist_High.')) return; do i = 1, p%NumTurbines - allocate ( y%Vdist_High(i)%data(3,0:p%nX_high-1,0:p%nY_high-1,0:p%nZ_high-1,0:p%n_high_low), STAT=ErrStat2 ); if (Failed0('y%Vdist_High%data.')) return; + allocate ( y%Vdist_High(i)%data(3,0:p%nX_high-1,0:p%nY_high-1,0:p%nZ_high-1,0:p%n_high_low_p1), STAT=ErrStat2 ); if (Failed0('y%Vdist_High%data.')) return; y%Vdist_High(i)%data = 0.0_Siki end do @@ -1200,7 +1202,7 @@ subroutine AWAE_Init( InitInp, u, p, x, xd, z, OtherState, y, m, Interval, InitO allocate ( m%Vamb_high(1:p%NumTurbines), STAT=ErrStat2 ); if (Failed0('Could not allocate memory for m%Vamb_high.')) return; do i = 1, p%NumTurbines - allocate ( m%Vamb_high(i)%data(3,0:p%nX_high-1,0:p%nY_high-1,0:p%nZ_high-1,0:p%n_high_low), STAT=ErrStat2 ); if (Failed0('m%Vamb_high%data.')) return; + allocate ( m%Vamb_high(i)%data(3,0:p%nX_high-1,0:p%nY_high-1,0:p%nZ_high-1,0:p%n_high_low_p1), STAT=ErrStat2 ); if (Failed0('m%Vamb_high%data.')) return; end do allocate ( m%parallelFlag( 0:p%NumPlanes-2,1:p%NumTurbines ), STAT=errStat2 ); if (Failed0('m%parallelFlag.')) return; @@ -1455,18 +1457,27 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM !$OMP PARALLEL DO DEFAULT(NONE) PRIVATE(nt, i_hl, errStat2, errMsg2) & !$OMP SHARED(p, n_high_low, n, m, errStat, errMsg, AbortErrLev) do nt = 1,p%NumTurbines + + ! Copy T=T_low_previous-DT_high (end-1 index in Vamb_high) into T=T_low_now-DT_high (0 index in Vamb_high) + if (n/=0_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,ubound(m%Vamb_high(nt)%data,5)-1) + do i_hl=0, n_high_low ! read from file the ambient flow for the current time step !FIXME:merge5.0 replace next line with the following - call ReadHighResWindFile(nt, (n+1)*p%n_high_low + i_hl, p, m%Vamb_high(nt)%data(:,:,:,:,i_hl), errStat2, errMsg2) -! call ReadHighResWindVTK(nt, n*p%n_high_low + i_hl, p, m%Vamb_high(nt)%data(:,:,:,:,i_hl), errStat2, errMsg2) + call ReadHighResWindFile(nt, (n+1)*p%n_high_low + i_hl, p, m%Vamb_high(nt)%data(:,:,:,:,i_hl+1), errStat2, errMsg2) +! call ReadHighResWindVTK(nt, n*p%n_high_low + i_hl, p, m%Vamb_high(nt)%data(:,:,:,:,i_hl+1), errStat2, errMsg2) if (ErrStat2 >= AbortErrLev) then !$OMP CRITICAL ! Needed to avoid data race on ErrStat and ErrMsg call SetErrStat( ErrStat2, ErrMsg2, errStat, errMsg, RoutineName ) !$OMP END CRITICAL endif end do + + ! Special handling at T=0 for time slice at -DT_high (0 index in Vamb_high) + ! -> Copy T=0 data into T=-DT_high for AD extrap/interp + if (n==0_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,1) + end do !$OMP END PARALLEL DO @@ -1487,6 +1498,9 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM ! Set input position m%u_IfW_High%PositionXYZ = p%Grid_high(:,:,nt) + ! Copy T=T_low_previous-DT_high (end-1 index in Vamb_high) into T=T_low_now-DT_high (0 index in Vamb_high) + if (n/=0_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,ubound(m%Vamb_high(nt)%data,5)-1) + ! Loop through high resolution grids do i_hl = 0, n_high_low @@ -1508,8 +1522,13 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM ! lbound(m%Vamb_high(nt)%data,2):ubound(m%Vamb_high(nt)%data,2),& ! lbound(m%Vamb_high(nt)%data,3):ubound(m%Vamb_high(nt)%data,3),& ! lbound(m%Vamb_high(nt)%data,4):ubound(m%Vamb_high(nt)%data,4)) => m%y_IfW_High(nt)%VelocityUVW - m%Vamb_high(nt)%data(:,:,:,:,i_hl) = V_Grid + m%Vamb_high(nt)%data(:,:,:,:,i_hl+1) = V_Grid end do + + ! Special handling at T=0 for time slice at -DT_high (0 index in Vamb_high) + ! -> Copy T=0 data into T=-DT_high for AD extrap/interp + if (n==0_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,1) + end do ! Multiple InflowWind instances (one per turbine) @@ -1530,8 +1549,11 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM end do end do + ! Copy T=T_low_previous-DT_high (end-1 index in Vamb_high) into T=T_low_now-DT_high (0 index in Vamb_high) + if (n/=0_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,ubound(m%Vamb_high(nt)%data,5)-1) + ! Loop through high resolution grids - do i_hl = 0, n_high_low + do i_hl = -1, n_high_low !FIXME:merge5.0 remove next 4 lines ! Set the hub position and orientation to pass to IfW (IfW always calculates hub and disk avg vel) @@ -1556,8 +1578,13 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM ! lbound(m%Vamb_high(nt)%data,2):ubound(m%Vamb_high(nt)%data,2),& ! lbound(m%Vamb_high(nt)%data,3):ubound(m%Vamb_high(nt)%data,3),& ! lbound(m%Vamb_high(nt)%data,4):ubound(m%Vamb_high(nt)%data,4)) => m%y_IfW_High(nt)%VelocityUVW - m%Vamb_high(nt)%data(:,:,:,:,i_hl) = V_Grid + m%Vamb_high(nt)%data(:,:,:,:,i_hl+1) = V_Grid end do + + ! Special handling at T=0 for time slice at -DT_high + ! -> Copy T=0 data into T=-DT_high for AD extrap/interp + if (n==0_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,1) + end do end select diff --git a/modules/awae/src/AWAE_Registry.txt b/modules/awae/src/AWAE_Registry.txt index e6afddb9d0..5627f31bf3 100644 --- a/modules/awae/src/AWAE_Registry.txt +++ b/modules/awae/src/AWAE_Registry.txt @@ -188,6 +188,7 @@ typedef ^ ParameterType ReKi Grid_low {:}{:} - - "XYZ typedef ^ ParameterType ReKi Grid_high {:}{:}{:} - - "XYZ components (global positions) of the spatial discretization of the high-resolution spatial domain for each turbine " m typedef ^ ParameterType ReKi WT_Position {:}{:} - - "X-Y-Z position of each wind turbine; index 1 = XYZ; index 2 = turbine number" meters typedef ^ ParameterType IntKi n_high_low - - - "Number of high-resolution time steps per low" - +typedef ^ ParameterType IntKi n_high_low_p1 - - - "Number of high-resolution time steps per low, plus one at t_low-dt_high" - typedef ^ ParameterType DbKi dt_low - - - "Low-resolution (FAST.Farm driver/glue code) time step" s typedef ^ ParameterType DbKi dt_high - - - "High-resolution (FAST) time step" s typedef ^ ParameterType IntKi NumDT - - - "Number of low-resolution (FAST.Farm driver/glue code) time steps" - diff --git a/modules/awae/src/AWAE_Types.f90 b/modules/awae/src/AWAE_Types.f90 index fbd5bde33a..3cbe84a334 100644 --- a/modules/awae/src/AWAE_Types.f90 +++ b/modules/awae/src/AWAE_Types.f90 @@ -211,6 +211,7 @@ MODULE AWAE_Types REAL(ReKi) , DIMENSION(:,:,:), ALLOCATABLE :: Grid_high !< XYZ components (global positions) of the spatial discretization of the high-resolution spatial domain for each turbine [m] REAL(ReKi) , DIMENSION(:,:), ALLOCATABLE :: WT_Position !< X-Y-Z position of each wind turbine; index 1 = XYZ; index 2 = turbine number [meters] INTEGER(IntKi) :: n_high_low = 0_IntKi !< Number of high-resolution time steps per low [-] + INTEGER(IntKi) :: n_high_low_p1 = 0_IntKi !< Number of high-resolution time steps per low, plus one at t_low-dt_high [-] REAL(DbKi) :: dt_low = 0.0_R8Ki !< Low-resolution (FAST.Farm driver/glue code) time step [s] REAL(DbKi) :: dt_high = 0.0_R8Ki !< High-resolution (FAST) time step [s] INTEGER(IntKi) :: NumDT = 0_IntKi !< Number of low-resolution (FAST.Farm driver/glue code) time steps [-] @@ -1999,6 +2000,7 @@ subroutine AWAE_CopyParam(SrcParamData, DstParamData, CtrlCode, ErrStat, ErrMsg) DstParamData%WT_Position = SrcParamData%WT_Position end if DstParamData%n_high_low = SrcParamData%n_high_low + DstParamData%n_high_low_p1 = SrcParamData%n_high_low_p1 DstParamData%dt_low = SrcParamData%dt_low DstParamData%dt_high = SrcParamData%dt_high DstParamData%NumDT = SrcParamData%NumDT @@ -2221,6 +2223,7 @@ subroutine AWAE_PackParam(RF, Indata) call RegPackAlloc(RF, InData%Grid_high) call RegPackAlloc(RF, InData%WT_Position) call RegPack(RF, InData%n_high_low) + call RegPack(RF, InData%n_high_low_p1) call RegPack(RF, InData%dt_low) call RegPack(RF, InData%dt_high) call RegPack(RF, InData%NumDT) @@ -2306,6 +2309,7 @@ subroutine AWAE_UnPackParam(RF, OutData) call RegUnpackAlloc(RF, OutData%Grid_high); if (RegCheckErr(RF, RoutineName)) return call RegUnpackAlloc(RF, OutData%WT_Position); if (RegCheckErr(RF, RoutineName)) return call RegUnpack(RF, OutData%n_high_low); if (RegCheckErr(RF, RoutineName)) return + call RegUnpack(RF, OutData%n_high_low_p1); if (RegCheckErr(RF, RoutineName)) return call RegUnpack(RF, OutData%dt_low); if (RegCheckErr(RF, RoutineName)) return call RegUnpack(RF, OutData%dt_high); if (RegCheckErr(RF, RoutineName)) return call RegUnpack(RF, OutData%NumDT); if (RegCheckErr(RF, RoutineName)) return From bf39a2f73cb1791b0630329df2d3b2677d1ed482 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Mon, 9 Mar 2026 16:46:59 -0600 Subject: [PATCH 20/27] AWAE: typo in loop limits on Mod_AmbWind==3 --- modules/awae/src/AWAE.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/modules/awae/src/AWAE.f90 b/modules/awae/src/AWAE.f90 index e304934f25..58f94e2ed4 100644 --- a/modules/awae/src/AWAE.f90 +++ b/modules/awae/src/AWAE.f90 @@ -1553,7 +1553,7 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM if (n/=0_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,ubound(m%Vamb_high(nt)%data,5)-1) ! Loop through high resolution grids - do i_hl = -1, n_high_low + do i_hl = 0, n_high_low !FIXME:merge5.0 remove next 4 lines ! Set the hub position and orientation to pass to IfW (IfW always calculates hub and disk avg vel) From 0efb7811eaf2f3452a3ef7e1c65a3cbfd7a05024 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Tue, 10 Mar 2026 00:13:43 -0600 Subject: [PATCH 21/27] AWAE: fix bug in calcoutput end time plane check --- modules/awae/src/AWAE.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/modules/awae/src/AWAE.f90 b/modules/awae/src/AWAE.f90 index 58f94e2ed4..582bec550c 100644 --- a/modules/awae/src/AWAE.f90 +++ b/modules/awae/src/AWAE.f90 @@ -770,7 +770,7 @@ subroutine HighResGridCalcOutput(n, u, p, xd, y, m, errStat, errMsg) maxPln = min(n,p%NumPlanes-2) ! We only need one high res file for that last simulation time - if ( (n/p%n_high_low_p1) == (p%NumDT-1) ) then + if ( (n/p%n_high_low) == (p%NumDT-1) ) then n_high_low_p1 = 1 else n_high_low_p1 = p%n_high_low_p1 From 756e48c8fa5045bbce3998dfe40f2df1a735648f Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Tue, 10 Mar 2026 00:34:35 -0600 Subject: [PATCH 22/27] AWAE: fix index bug in high res output vtk writing --- modules/awae/src/AWAE_IO.f90 | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/modules/awae/src/AWAE_IO.f90 b/modules/awae/src/AWAE_IO.f90 index 86ee0b77c5..17e79eed66 100644 --- a/modules/awae/src/AWAE_IO.f90 +++ b/modules/awae/src/AWAE_IO.f90 @@ -110,13 +110,14 @@ subroutine WriteDisWindFiles( n, WrDisSkp1, p, y, m, errStat, errMsg ) if (ErrStat >= AbortErrLev) return do nt= 1,p%NumTurbines - ! We are only writing out the first of the high res data for a given low res time step + ! We are only writing out the first of the high res data for a given low res time step + ! NOTE: y%Vdist_high(nt)%data(:,:,:,:,1) is at T=t_low, and index 0 is at T=t_low-DT_high FileName = trim(p%OutFileVTKRoot)//".HighT"//trim(num2lstr(nt))//".Dis."//trim(Tstr)//".vtk" call WrVTK_SP_header( FileName, "High resolution disturbed wind for time = "//trim(num2lstr(t_out))//" seconds.", Un, errStat2, errMsg2 ) call SetErrStat(errStat2, errMsg2, ErrStat, ErrMsg, RoutineName) if (ErrStat >= AbortErrLev) return - call WrVTK_SP_vectors3D( Un, "Velocity", (/p%nX_high,p%nY_high,p%nZ_high/), (/p%X0_high(nt),p%Y0_high(nt),p%Z0_high(nt)/), (/p%dX_high(nt),p%dY_high(nt),p%dZ_high(nt)/), y%Vdist_high(nt)%data(:,:,:,:,0), errStat2, errMsg2 ) + call WrVTK_SP_vectors3D( Un, "Velocity", (/p%nX_high,p%nY_high,p%nZ_high/), (/p%X0_high(nt),p%Y0_high(nt),p%Z0_high(nt)/), (/p%dX_high(nt),p%dY_high(nt),p%dZ_high(nt)/), y%Vdist_high(nt)%data(:,:,:,:,1), errStat2, errMsg2 ) call SetErrStat(ErrStat2, errMsg2, ErrStat, ErrMsg, RoutineName) if (ErrStat >= AbortErrLev) return From 9845567a4623e57827cae0484f1c8777e0c2877e Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Tue, 10 Mar 2026 01:20:29 -0600 Subject: [PATCH 23/27] AWAE: fix bug in time step indexing (starts n=-1) and WAT time indexing --- modules/awae/src/AWAE.f90 | 27 +++++++++++++-------------- 1 file changed, 13 insertions(+), 14 deletions(-) diff --git a/modules/awae/src/AWAE.f90 b/modules/awae/src/AWAE.f90 index 582bec550c..2eb80fff69 100644 --- a/modules/awae/src/AWAE.f90 +++ b/modules/awae/src/AWAE.f90 @@ -789,8 +789,7 @@ subroutine HighResGridCalcOutput(n, u, p, xd, y, m, errStat, errMsg) allocate ( WAT_B_BoxHi ( 3, 0:n_high_low_p1), STAT=errStat2 ); if (errStat2 /= 0) call SetErrStat ( ErrID_Fatal, 'Could not allocate memory for WAT_B_BoxHi.', errStat, errMsg, RoutineName ) if (ErrStat >= AbortErrLev) return do i_hl=0, n_high_low_p1 - ! NOTE: i_hl=0 corresponds to T=-DT_high, so shift by one high res timestep - WAT_B_BoxHi(1:3, i_hl) = xd%WAT_B_Box(1:3) - (n_high_low_p1-i_hl+1) * xd%Ufarm(1:3) * real(p%DT_high,ReKi) + WAT_B_BoxHi(1:3, i_hl) = xd%WAT_B_Box(1:3) - (n_high_low_p1-i_hl) * xd%Ufarm(1:3) * real(p%DT_high,ReKi) enddo endif @@ -1458,8 +1457,8 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM !$OMP SHARED(p, n_high_low, n, m, errStat, errMsg, AbortErrLev) do nt = 1,p%NumTurbines - ! Copy T=T_low_previous-DT_high (end-1 index in Vamb_high) into T=T_low_now-DT_high (0 index in Vamb_high) - if (n/=0_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,ubound(m%Vamb_high(nt)%data,5)-1) + ! Copy T=T_low_previous-DT_high (end-1 index in Vamb_high) into T=T_low_now-DT_high (0 index in Vamb_high). Note that n starts at -1 + if (n/=-1_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,ubound(m%Vamb_high(nt)%data,5)-1) do i_hl=0, n_high_low @@ -1474,9 +1473,9 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM endif end do - ! Special handling at T=0 for time slice at -DT_high (0 index in Vamb_high) + ! Special handling at T=0 for time slice at -DT_high (0 index in Vamb_high). Note that n starts at -1 ! -> Copy T=0 data into T=-DT_high for AD extrap/interp - if (n==0_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,1) + if (n==-1_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,1) end do !$OMP END PARALLEL DO @@ -1498,8 +1497,8 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM ! Set input position m%u_IfW_High%PositionXYZ = p%Grid_high(:,:,nt) - ! Copy T=T_low_previous-DT_high (end-1 index in Vamb_high) into T=T_low_now-DT_high (0 index in Vamb_high) - if (n/=0_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,ubound(m%Vamb_high(nt)%data,5)-1) + ! Copy T=T_low_previous-DT_high (end-1 index in Vamb_high) into T=T_low_now-DT_high (0 index in Vamb_high). Note that n starts at -1 + if (n/=-1_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,ubound(m%Vamb_high(nt)%data,5)-1) ! Loop through high resolution grids do i_hl = 0, n_high_low @@ -1525,9 +1524,9 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM m%Vamb_high(nt)%data(:,:,:,:,i_hl+1) = V_Grid end do - ! Special handling at T=0 for time slice at -DT_high (0 index in Vamb_high) + ! Special handling at T=0 for time slice at -DT_high (0 index in Vamb_high). Note that n starts at -1 ! -> Copy T=0 data into T=-DT_high for AD extrap/interp - if (n==0_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,1) + if (n==-1_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,1) end do @@ -1549,8 +1548,8 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM end do end do - ! Copy T=T_low_previous-DT_high (end-1 index in Vamb_high) into T=T_low_now-DT_high (0 index in Vamb_high) - if (n/=0_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,ubound(m%Vamb_high(nt)%data,5)-1) + ! Copy T=T_low_previous-DT_high (end-1 index in Vamb_high) into T=T_low_now-DT_high (0 index in Vamb_high). Note that n starts at -1 + if (n/=-1_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,ubound(m%Vamb_high(nt)%data,5)-1) ! Loop through high resolution grids do i_hl = 0, n_high_low @@ -1581,9 +1580,9 @@ subroutine AWAE_UpdateStates( t, n, u, p, x, xd, z, OtherState, m, errStat, errM m%Vamb_high(nt)%data(:,:,:,:,i_hl+1) = V_Grid end do - ! Special handling at T=0 for time slice at -DT_high + ! Special handling at T=0 for time slice at -DT_high. Note that n starts at -1 ! -> Copy T=0 data into T=-DT_high for AD extrap/interp - if (n==0_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,1) + if (n==-1_IntKi) m%Vamb_high(nt)%data(:,:,:,:,0) = m%Vamb_high(nt)%data(:,:,:,:,1) end do From 4ad9c3fd2f01fc20863f67bcf47d90753037c5f0 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Tue, 10 Mar 2026 10:33:09 -0600 Subject: [PATCH 24/27] release notes for 4.2.1 --- docs/changelogs/v4.2.1.md | 136 +++++++++++++++++++++++++++++++ docs/conf.py | 2 +- docs/source/user/api_change.rst | 7 ++ glue-codes/python/pyproject.toml | 2 +- openfast_io/pyproject.toml | 2 +- 5 files changed, 146 insertions(+), 3 deletions(-) create mode 100644 docs/changelogs/v4.2.1.md diff --git a/docs/changelogs/v4.2.1.md b/docs/changelogs/v4.2.1.md new file mode 100644 index 0000000000..1f0405a1cd --- /dev/null +++ b/docs/changelogs/v4.2.1.md @@ -0,0 +1,136 @@ +**Feature or improvement description** +Pull request to merge `dev` into `main` for release version 4.2.1 + +See the milestone and project pages for additional information + + https://github.com/OpenFAST/openfast/milestone/27 + +No test results change + +### Release checklist: +- [ ] Update the documentation version in docs/conf.py +- [ ] Update the versions in docs/source/user/api\_change.rst +- [ ] Update version info in openfast\_io/pyproject.toml (`openfast_io` package) +- [ ] Update version info in glue-codes/python/pyproject.toml (`pyOpenFAST` package for testing) +- [ ] Verify readthedocs builds correctly +- [ ] Create an annotated tag in OpenFAST during merge (mark as most recent if necessary) +- [ ] Create a merge commit in r-test and add a corresponding annotated tag +- [ ] Upload Docker image +- [ ] Compile executables for Windows builds + - [ ] `AeroDisk_Driver_x64.exe` + - [ ] `AeroDyn_Driver_x64.exe` + - [ ] `AeroDyn_Driver_x64_OpenMP.exe` + - [ ] `AeroDyn_Inflow_c_binding_x64.dll` + - [ ] `AeroDyn_Inflow_c_binding_x64_OpenMP.dll` + - [ ] `BeamDyn_Driver_x64.exe` + - [ ] `DISCON.dll (x64)` + - [ ] `DISCON_ITIBarge.dll (x64)` + - [ ] `DISCON_OC3Hywind.dll (x64)` + - [ ] `FAST.Farm_x64.exe` + - [ ] `FAST.Farm_x64_OMP.exe` + - [ ] `FAST_SFunc.mexw64` + - [ ] `HydroDynDriver_x64.exe` + - [ ] `HydroDyn_C_Binding_x64.dll` + - [ ] `IinflowWind_c_binding_x64.dll` + - [ ] `InflowWind_Driver_x64.exe` + - [ ] `InflowWind_Driver_x64_OpenMP.exe` + - [ ] `MoorDyn_Driver_x64.exe` + - [ ] `MoorDyn_c_binding_x64.dll` + - [ ] `OpenFAST-Simulink_x64.dll` + - [ ] `openfast_x64.exe` + - [ ] `SeaStateDriver_x64.exe` + - [ ] `SeaState_c_binding_x64.dll` + - [ ] `SimpleElastoDyn_x64.exe` + - [ ] `SubDyn_x64.exe` + - [ ] `Turbsim_x64.exe` + - [ ] `UnsteadyAero_x64.exe` + + + +# Release Overview +------------------ + + + +### Contribution Acknowledgements + +We are grateful for a first time code contributions from @IrisMeasure. + + + +# Changelog (from 4.2.1) +------------------------ + + +## Solvers + +### OpenFAST + +#3243 Add high res windslice at `T=t_low - DT_high` (Fix for FAST.Farm issue #3183) (@andrew-platt) + + + +## Modules + +## AeroDyn + +#### AeroDyn Driver / AeroDyn\_Inflow\_C\_Bindings interface + +#3214 [BugFix] ADI c-bind: comment out code that causes segfaults with interface debugging (@andrew-platt) + +#3227 [BugFix] AD: DTAero from input file was ignored (@andrew-platt) + +#3230 [BugFix] AD driver: hub acceleration value incorrect in RotMotionType==1 (@andrew-platt) + + +### NWTC-Library + +#3223 Switch from xPPTRF to xPOTRF to improve TurbSim speed on macOS (@IrisMeasure) + +#3229 Re-add the xPPTRF routines (rename xPOTRF routines from #3123) (@andrew-platt) + + +### SeaState + +#3235 backport of #3231 - SS WaveTp logic (@andrew-platt, @RBergua) + +#3236 Add SeaSt\_C\_GetFluidDynP routine to SeaState\_C\_Binding interface (@andrew-platt) + + + +## Testing and input file processing + +### Regression and Unit testing + +#3228 regtest: check if destination dirs exist before copyTree call (@andrew-platt) + + + +## Input file changes + +No input file changes since v4.2.0 + + +## Known issues +There are several issues that have not been addressed in this release due to time constraints, but will be addressed in future releases. These include: + +- No visualization of rectangular members from _HydroDyn_ or _SubDyn_ through the VTK output options +- Missing and broken features from several c-binding library interfaces: + - the _AeroDyn\_Inflow\_c-binding_ library interface does not allow for coupling to the tower. This will require an interface update. + - the _HydroDyn\_c-binding_ library interface does not currently support vizualization. This will require an interface update. + - the `InitNodePositions` input to _HydroDyn\_c-binding_ library interface does not currently work with any non-zero `x` or `y` coordinates (non-zero `z` is ok) + - the _MoorDyn\_c-binding_ library interface does not currently support vizualization. This will require an interface update. +- Documentation on the new _pyOpenFAST_ module is incomplete. Partial documentation exists on how to use it in regression testing, but no documentation or examples exist on using it to call c-bindings modules from Python. +- Documentation is incomplete for _HydroDyn_, _SubDyn_, and a few other modules. + + +# Precompiled Windows Binaries +The binary files in this release were built with the Visual Studio solution files distributed with OpenFAST (not using cmake), using + +- Intel Fortran Essentials 2025.3.0.333 +- Microsoft Visual Studio 2022 Version 17.14.23. +- MATLAB 2025.2.999 (R2025b) +- Executables with `_OpenMP` or `_OMP` in the name are built with OpenMP libraries and linked with dynamic libraries. + - You will need [this Intel Fortran redistributable package](https://registrationcenter-download.intel.com/akdlm/IRC_NAS/0dc56e76-d2c0-4bb8-9c83-c2ee3952b855/w_ifx_runtime_p_2025.2.1.1001.exe) installed to use these executables if you do not already have Intel Fortran OneAPI 2024 installed. See the installation instructions [here](https://software.intel.com/content/www/us/en/develop/articles/redistributable-libraries-for-intel-c-and-fortran-2022-compilers-for-windows.html). + +**The other OpenFAST executables DO NOT require these redistributable libraries to be installed. Instead, they were built with static libraries.** diff --git a/docs/conf.py b/docs/conf.py index 4733a3a914..41518da7f5 100644 --- a/docs/conf.py +++ b/docs/conf.py @@ -139,7 +139,7 @@ def runDoxygen(sourcfile, doxyfileIn, doxyfileOut): # The short X.Y version. version = f'4.2' # The full version, including alpha/beta/rc tags. -release = f'v4.2.0' +release = f'v4.2.1' # The language for content autogenerated by Sphinx. Refer to documentation # for a list of supported languages. diff --git a/docs/source/user/api_change.rst b/docs/source/user/api_change.rst index 1a11513ffe..ca51edbf05 100644 --- a/docs/source/user/api_change.rst +++ b/docs/source/user/api_change.rst @@ -10,6 +10,13 @@ The line number corresponds to the resulting line number after all changes are i Thus, be sure to implement each in order so that subsequent line numbers are correct. +OpenFAST v4.2.0 to OpenFAST v4.2.1 +---------------------------------- + +No input file changes were made. + + + OpenFAST v4.1.x to OpenFAST v4.2.0 ---------------------------------- diff --git a/glue-codes/python/pyproject.toml b/glue-codes/python/pyproject.toml index 9439622819..a56ca79d8d 100644 --- a/glue-codes/python/pyproject.toml +++ b/glue-codes/python/pyproject.toml @@ -4,7 +4,7 @@ build-backend = "setuptools.build_meta" [project] name = "pyOpenFAST" -version = "4.2.0" +version = "4.2.1" description = "Python interface to OpenFAST FAST Library and physics modules." readme = "README.md" requires-python = ">=3.9" diff --git a/openfast_io/pyproject.toml b/openfast_io/pyproject.toml index 75074c6760..965631b514 100644 --- a/openfast_io/pyproject.toml +++ b/openfast_io/pyproject.toml @@ -5,7 +5,7 @@ build-backend = "hatchling.build" [project] name = "openfast_io" # dynamic = ["version"] -version = "4.2.0" +version = "4.2.1" description = "Readers and writers for OpenFAST files." license = {file = "../LICENSE"} authors = [ From 79ce33b85abe1b3583055054c505328a3624aefe Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Tue, 10 Mar 2026 10:43:27 -0600 Subject: [PATCH 25/27] v4.2.1.md -- fix typo --- docs/changelogs/v4.2.1.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/changelogs/v4.2.1.md b/docs/changelogs/v4.2.1.md index 1f0405a1cd..0fcb9cf95a 100644 --- a/docs/changelogs/v4.2.1.md +++ b/docs/changelogs/v4.2.1.md @@ -85,7 +85,7 @@ We are grateful for a first time code contributions from @IrisMeasure. ### NWTC-Library -#3223 Switch from xPPTRF to xPOTRF to improve TurbSim speed on macOS (@IrisMeasure) +#3123 Switch from xPPTRF to xPOTRF to improve TurbSim speed on macOS (@IrisMeasure) #3229 Re-add the xPPTRF routines (rename xPOTRF routines from #3123) (@andrew-platt) From 71654a487eef4277ba171f854892fb99042bcac4 Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Tue, 10 Mar 2026 10:49:12 -0600 Subject: [PATCH 26/27] release notes: add missing PR --- docs/changelogs/v4.2.1.md | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/docs/changelogs/v4.2.1.md b/docs/changelogs/v4.2.1.md index 0fcb9cf95a..37b4fd4c4b 100644 --- a/docs/changelogs/v4.2.1.md +++ b/docs/changelogs/v4.2.1.md @@ -97,6 +97,11 @@ We are grateful for a first time code contributions from @IrisMeasure. #3236 Add SeaSt\_C\_GetFluidDynP routine to SeaState\_C\_Binding interface (@andrew-platt) +### TurbSim + +#3224 Fix TurbSim compile with OpenMP and IFX compiler (@deslaughter) + + ## Testing and input file processing From 5ca09e9cae65c417db09edf2fbf34159b3283cee Mon Sep 17 00:00:00 2001 From: andrew-platt Date: Tue, 10 Mar 2026 12:12:01 -0600 Subject: [PATCH 27/27] release notes typo --- docs/changelogs/v4.2.1.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/changelogs/v4.2.1.md b/docs/changelogs/v4.2.1.md index 37b4fd4c4b..917b60821f 100644 --- a/docs/changelogs/v4.2.1.md +++ b/docs/changelogs/v4.2.1.md @@ -64,7 +64,7 @@ We are grateful for a first time code contributions from @IrisMeasure. ## Solvers -### OpenFAST +### FAST.Farm #3243 Add high res windslice at `T=t_low - DT_high` (Fix for FAST.Farm issue #3183) (@andrew-platt)