From 4f1c8450907803dd9e53b853513f3444af953fcf Mon Sep 17 00:00:00 2001 From: Roger Bergua Date: Wed, 3 Jun 2026 12:50:23 -0600 Subject: [PATCH 1/9] Update beams self-weight (only sensor outputs) --- modules/subdyn/src/SubDyn_Output.f90 | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/modules/subdyn/src/SubDyn_Output.f90 b/modules/subdyn/src/SubDyn_Output.f90 index 970b6de36..a5d1602ad 100644 --- a/modules/subdyn/src/SubDyn_Output.f90 +++ b/modules/subdyn/src/SubDyn_Output.f90 @@ -247,10 +247,9 @@ SUBROUTINE ConfigOutputNode_MKF_ID(pLst, iElem, iiNode, iStore, NodeID2) CALL ElemM(p%ElemProps(iElem), pLst%Me(:,:,iiNode,iStore)) CALL ElemK(p%ElemProps(iElem), pLst%Ke(:,:,iiNode,iStore)) CALL ElemF(p%ElemProps(iElem), Init%g, pLst%Fg(:,iiNode,iStore), FCe) - ! NOTE: Removing this force contribution for now - ! The output of subdyn will just be the "Kx" part for now - !pLst%Fg(:,iiNode,iStore) = pLst%Fg(:,iiNode,iStore) + FCe(1:12) ! Adding cable element force - pLst%Fg(:,iiNode,iStore) = FCe(1:12) ! Adding cable element force + ! Apply superposition correction: include both the element gravity force vector + ! and the cable pretension nodal force when evaluating internal member loads. + pLst%Fg(:,iiNode,iStore) = pLst%Fg(:,iiNode,iStore) + FCe(1:12) ! gravity force + cable element force END SUBROUTINE ConfigOutputNode_MKF_ID From 65d51f86992cc05cc27a6efc980dbd7cf20e87de Mon Sep 17 00:00:00 2001 From: Roger Bergua Date: Thu, 4 Jun 2026 08:40:21 -0600 Subject: [PATCH 2/9] Fix (SubDyn): remove incorrect comment about Fg --- modules/subdyn/src/SubDyn_Output.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/modules/subdyn/src/SubDyn_Output.f90 b/modules/subdyn/src/SubDyn_Output.f90 index a5d1602ad..f66e1f95f 100644 --- a/modules/subdyn/src/SubDyn_Output.f90 +++ b/modules/subdyn/src/SubDyn_Output.f90 @@ -532,7 +532,7 @@ SUBROUTINE CALC_NODE_FORCES(DIRCOS,Me,Ke,Udotdot,Y2 ,Fg, FirstOrSecond, FM_nod, FM_glb = matmul(Me,Udotdot) ! GLOBAL REFERENCE FF_glb = matmul(Ke,Y2) ! GLOBAL REFERENCE - FF_glb = FF_glb - Fg ! GLOBAL REFERENCE ! NOTE: Fg is now 0, only the "Kx" part in Fk + FF_glb = FF_glb - Fg ! GLOBAL REFERENCE DO L=1,4 ! Transforming coordinates 3 at a time FM_elm((L-1)*3+1:L*3) = matmul(DIRCOS, FM_glb( (L-1)*3+1:L*3 ) ) FF_elm((L-1)*3+1:L*3) = matmul(DIRCOS, FF_glb( (L-1)*3+1:L*3 ) ) From 65cbd8aa77cac966f90af855a0f0f7bb29ca85df Mon Sep 17 00:00:00 2001 From: Roger Bergua Date: Thu, 4 Jun 2026 09:15:48 -0600 Subject: [PATCH 3/9] SubDyn: update CALC_NODE_FORCES description and comments --- modules/subdyn/src/SubDyn_Output.f90 | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/modules/subdyn/src/SubDyn_Output.f90 b/modules/subdyn/src/SubDyn_Output.f90 index f66e1f95f..61e80cea9 100644 --- a/modules/subdyn/src/SubDyn_Output.f90 +++ b/modules/subdyn/src/SubDyn_Output.f90 @@ -515,15 +515,15 @@ subroutine ElementForce(pLst, iiNode, JJ, FM_elm, FK_elm, sgn, DIRCOS, bUseInput end subroutine ElementForce !==================================================================================================== - !> Calculates static and dynamic forces for a given element, using K and M of the element, and - !output quantities Udotdot and Y2 containing the - !and K2 indicating wheter the 1st (1) or 2nd (2) node is to be picked + !> Calculates static and dynamic nodal forces for a given element using its stiffness and mass matrices. + ! Udotdot contains nodal accelerations and Y2 contains nodal displacements. + ! FirstOrSecond selects whether the node of interest is the first (1) or second (2) node of the element. !---------------------------------------------------------------------------------------------------- SUBROUTINE CALC_NODE_FORCES(DIRCOS,Me,Ke,Udotdot,Y2 ,Fg, FirstOrSecond, FM_nod, FK_nod) Real(FEKi), DIMENSION (3,3), INTENT(IN) :: DIRCOS !direction cosice matrix (global to local) (3x3) Real(FEKi), DIMENSION (12,12), INTENT(IN) :: Me,Ke !element M and K matrices (12x12) in GLOBAL REFERENCE (DIRCOS^T K DIRCOS) - Real(ReKi), DIMENSION (12), INTENT(IN) :: Udotdot, Y2 !acceleration and velocities, gravity forces - Real(FEKi), DIMENSION (12), INTENT(IN) :: Fg !acceleration and velocities, gravity forces + Real(ReKi), DIMENSION (12), INTENT(IN) :: Udotdot, Y2 ! nodal accelerations and nodal displacements + Real(FEKi), DIMENSION (12), INTENT(IN) :: Fg ! constant element load vector (gravity + cable) Integer(IntKi), INTENT(IN) :: FirstOrSecond !1 or 2 depending on node of interest REAL(ReKi), DIMENSION (6), INTENT(OUT) :: FM_nod, FK_nod !output static and dynamic forces and moments !Locals @@ -546,7 +546,7 @@ END SUBROUTINE SDOut_MapOutputs !==================================================================================================== SUBROUTINE SDOut_CloseSum( UnSum, ErrStat, ErrMsg ) - INTEGER, INTENT( IN ) :: UnSum ! the unit number for the SubDyn summary file + INTEGER, INTENT( IN ) :: UnSum ! the unit number for the SubDyn summary file INTEGER, INTENT( OUT ) :: ErrStat ! returns a non-zero value when an error occurs CHARACTER(*), INTENT( OUT ) :: ErrMsg ! Error message if ErrStat /= ErrID_None ! Local variables From 4ec470daa009685d68ca83d7b50f45e7da989184 Mon Sep 17 00:00:00 2001 From: Roger Bergua Date: Fri, 5 Jun 2026 17:34:43 -0600 Subject: [PATCH 4/9] SubDyn: self-weight reconstruction in floating systems --- modules/subdyn/src/SubDyn_Output.f90 | 17 ++++++++++++++--- 1 file changed, 14 insertions(+), 3 deletions(-) diff --git a/modules/subdyn/src/SubDyn_Output.f90 b/modules/subdyn/src/SubDyn_Output.f90 index 61e80cea9..407bbaca8 100644 --- a/modules/subdyn/src/SubDyn_Output.f90 +++ b/modules/subdyn/src/SubDyn_Output.f90 @@ -248,7 +248,10 @@ SUBROUTINE ConfigOutputNode_MKF_ID(pLst, iElem, iiNode, iStore, NodeID2) CALL ElemK(p%ElemProps(iElem), pLst%Ke(:,:,iiNode,iStore)) CALL ElemF(p%ElemProps(iElem), Init%g, pLst%Fg(:,iiNode,iStore), FCe) ! Apply superposition correction: include both the element gravity force vector - ! and the cable pretension nodal force when evaluating internal member loads. + ! and the initial cable pretension nodal force when evaluating internal member loads. + ! Note: pLst%Fg stores only element gravity + cable (excludes lumped masses, which are + ! external point forces, not internal element loads). For floating systems, these forces + ! are rotated in ElementForce using RRg2b to account for platform rigid body motion. pLst%Fg(:,iiNode,iStore) = pLst%Fg(:,iiNode,iStore) + FCe(1:12) ! gravity force + cable element force END SUBROUTINE ConfigOutputNode_MKF_ID @@ -498,6 +501,7 @@ subroutine ElementForce(pLst, iiNode, JJ, FM_elm, FK_elm, sgn, DIRCOS, bUseInput integer(IntKi) :: FirstOrSecond !< 1 or 2 if first node or second node integer(IntKi), dimension(2) :: ElemNodes ! Node IDs for element under consideration (may not be consecutive numbers) real(ReKi) , dimension(12) :: X_e, Xdd_e ! Displacement and acceleration for an element + real(R8Ki) , dimension(12) :: Fg_e ! Gravity + initial cable pretension load vector (rotated for floating, static for fixed) integer(IntKi), dimension(2), parameter :: NodeNumber_To_Sign = (/-1, +1/) iElem = pLst%ElmIDs(iiNode,JJ) ! element number @@ -508,10 +512,17 @@ subroutine ElementForce(pLst, iiNode, JJ, FM_elm, FK_elm, sgn, DIRCOS, bUseInput X_e(7:12) = m%U_full_elast (p%NodesDOF(ElemNodes(2))%List(1:6)) ! No additional transformation required Xdd_e(1:6) = matmul(RRg2b,m%U_full_dotdot(p%NodesDOF(ElemNodes(1))%List(1:6))) ! Transform acceleration to be back in the Guyan frame Xdd_e(7:12) = matmul(RRg2b,m%U_full_dotdot(p%NodesDOF(ElemNodes(2))%List(1:6))) + ! For floating, gravity forces need rotation, but use only element gravity (pLst%Fg), not concentrated masses + if (p%Floating) then + Fg_e(1:6) = matmul(RRg2b, real(pLst%Fg(1:6,iiNode,JJ), R8Ki)) + Fg_e(7:12) = matmul(RRg2b, real(pLst%Fg(7:12,iiNode,JJ), R8Ki)) + else + Fg_e = real(pLst%Fg(:,iiNode,JJ), R8Ki) ! Fixed-bottom: use static gravity + initial cable pretension (no rotation needed) + endif if (.not. bUseInputDirCos) then DIRCOS=transpose(p%ElemProps(iElem)%DirCos)! global to local endif - CALL CALC_NODE_FORCES( DIRCOS, pLst%Me(:,:,iiNode,JJ),pLst%Ke(:,:,iiNode,JJ), Xdd_e, X_e, pLst%Fg(:,iiNode,JJ), FirstOrSecond, FM_elm, FK_elm) + CALL CALC_NODE_FORCES( DIRCOS, pLst%Me(:,:,iiNode,JJ),pLst%Ke(:,:,iiNode,JJ), Xdd_e, X_e, Fg_e, FirstOrSecond, FM_elm, FK_elm) end subroutine ElementForce !==================================================================================================== @@ -523,7 +534,7 @@ SUBROUTINE CALC_NODE_FORCES(DIRCOS,Me,Ke,Udotdot,Y2 ,Fg, FirstOrSecond, FM_nod, Real(FEKi), DIMENSION (3,3), INTENT(IN) :: DIRCOS !direction cosice matrix (global to local) (3x3) Real(FEKi), DIMENSION (12,12), INTENT(IN) :: Me,Ke !element M and K matrices (12x12) in GLOBAL REFERENCE (DIRCOS^T K DIRCOS) Real(ReKi), DIMENSION (12), INTENT(IN) :: Udotdot, Y2 ! nodal accelerations and nodal displacements - Real(FEKi), DIMENSION (12), INTENT(IN) :: Fg ! constant element load vector (gravity + cable) + Real(FEKi), DIMENSION (12), INTENT(IN) :: Fg ! element load vector from gravity and initial cable pretension Integer(IntKi), INTENT(IN) :: FirstOrSecond !1 or 2 depending on node of interest REAL(ReKi), DIMENSION (6), INTENT(OUT) :: FM_nod, FK_nod !output static and dynamic forces and moments !Locals From 5160c031a668484ab7dc24da7bbfbc74e64e3afc Mon Sep 17 00:00:00 2001 From: Roger Bergua Date: Fri, 5 Jun 2026 22:09:34 -0600 Subject: [PATCH 5/9] SubDyn floating: equivalent nodal moments based on the actual direction cosine --- modules/subdyn/src/SubDyn_Output.f90 | 33 ++++++++++++++++++---------- 1 file changed, 22 insertions(+), 11 deletions(-) diff --git a/modules/subdyn/src/SubDyn_Output.f90 b/modules/subdyn/src/SubDyn_Output.f90 index 407bbaca8..5da5442f0 100644 --- a/modules/subdyn/src/SubDyn_Output.f90 +++ b/modules/subdyn/src/SubDyn_Output.f90 @@ -247,11 +247,9 @@ SUBROUTINE ConfigOutputNode_MKF_ID(pLst, iElem, iiNode, iStore, NodeID2) CALL ElemM(p%ElemProps(iElem), pLst%Me(:,:,iiNode,iStore)) CALL ElemK(p%ElemProps(iElem), pLst%Ke(:,:,iiNode,iStore)) CALL ElemF(p%ElemProps(iElem), Init%g, pLst%Fg(:,iiNode,iStore), FCe) - ! Apply superposition correction: include both the element gravity force vector + ! Apply superposition correction: include both the element gravity force vector ! and the initial cable pretension nodal force when evaluating internal member loads. - ! Note: pLst%Fg stores only element gravity + cable (excludes lumped masses, which are - ! external point forces, not internal element loads). For floating systems, these forces - ! are rotated in ElementForce using RRg2b to account for platform rigid body motion. + ! Note: for floating systems, pLst%Fg loads will be rotated in ElementForce into the current body/Guyan frame. pLst%Fg(:,iiNode,iStore) = pLst%Fg(:,iiNode,iStore) + FCe(1:12) ! gravity force + cable element force END SUBROUTINE ConfigOutputNode_MKF_ID @@ -501,7 +499,9 @@ subroutine ElementForce(pLst, iiNode, JJ, FM_elm, FK_elm, sgn, DIRCOS, bUseInput integer(IntKi) :: FirstOrSecond !< 1 or 2 if first node or second node integer(IntKi), dimension(2) :: ElemNodes ! Node IDs for element under consideration (may not be consecutive numbers) real(ReKi) , dimension(12) :: X_e, Xdd_e ! Displacement and acceleration for an element - real(R8Ki) , dimension(12) :: Fg_e ! Gravity + initial cable pretension load vector (rotated for floating, static for fixed) + real(FEKi) , dimension(12) :: Fg_e ! Gravity + initial cable pretension load vector (reorient for floating, constant for fixed-bottom) + real(FEKi) , dimension(3,3) :: CurDirCos ! Current element direction cosine matrix in the floating body frame + real(R8Ki) , dimension(3,3) :: Rb2g ! Rotation matrix from body/Guyan to global frame (transpose of Rg2b) integer(IntKi), dimension(2), parameter :: NodeNumber_To_Sign = (/-1, +1/) iElem = pLst%ElmIDs(iiNode,JJ) ! element number @@ -512,12 +512,22 @@ subroutine ElementForce(pLst, iiNode, JJ, FM_elm, FK_elm, sgn, DIRCOS, bUseInput X_e(7:12) = m%U_full_elast (p%NodesDOF(ElemNodes(2))%List(1:6)) ! No additional transformation required Xdd_e(1:6) = matmul(RRg2b,m%U_full_dotdot(p%NodesDOF(ElemNodes(1))%List(1:6))) ! Transform acceleration to be back in the Guyan frame Xdd_e(7:12) = matmul(RRg2b,m%U_full_dotdot(p%NodesDOF(ElemNodes(2))%List(1:6))) - ! For floating, gravity forces need rotation, but use only element gravity (pLst%Fg), not concentrated masses + ! Load the computed at initialization gravity + initial cable pretension load vector. + ! For floating systems, project the load vector into the current body/Guyan frame for Fx/Fy/Fz. + ! For the self-weight bending moments, recompute them based on the current element direction cosine matrix. + Fg_e = real(pLst%Fg(:,iiNode,JJ), R8Ki) if (p%Floating) then - Fg_e(1:6) = matmul(RRg2b, real(pLst%Fg(1:6,iiNode,JJ), R8Ki)) - Fg_e(7:12) = matmul(RRg2b, real(pLst%Fg(7:12,iiNode,JJ), R8Ki)) - else - Fg_e = real(pLst%Fg(:,iiNode,JJ), R8Ki) ! Fixed-bottom: use static gravity + initial cable pretension (no rotation needed) + Fg_e(1:6) = matmul(RRg2b, Fg_e(1:6)) ! project node-1 forces/moments into body frame + Fg_e(7:12) = matmul(RRg2b, Fg_e(7:12)) ! project node-2 forces/moments into body frame + + Rb2g = transpose(Rg2b) ! current body-to-global rotation needed for moment projection + CurDirCos = matmul(Rb2g, p%ElemProps(iElem)%DirCos) + + ! Element self-weight bending moment contribution at the nodes: + Fg_e(4) = -p%ElemProps(iElem)%Length**2 * p%ElemProps(iElem)%Rho * p%ElemProps(iElem)%Area * p%g / 12.0_FEKi * CurDirCos(2,3) + Fg_e(5) = p%ElemProps(iElem)%Length**2 * p%ElemProps(iElem)%Rho * p%ElemProps(iElem)%Area * p%g / 12.0_FEKi * CurDirCos(1,3) + Fg_e(10) = -Fg_e(4) + Fg_e(11) = -Fg_e(5) endif if (.not. bUseInputDirCos) then DIRCOS=transpose(p%ElemProps(iElem)%DirCos)! global to local @@ -528,13 +538,14 @@ end subroutine ElementForce !==================================================================================================== !> Calculates static and dynamic nodal forces for a given element using its stiffness and mass matrices. ! Udotdot contains nodal accelerations and Y2 contains nodal displacements. + ! Fg is the element gravity + initial cable pretension load vector. ! FirstOrSecond selects whether the node of interest is the first (1) or second (2) node of the element. !---------------------------------------------------------------------------------------------------- SUBROUTINE CALC_NODE_FORCES(DIRCOS,Me,Ke,Udotdot,Y2 ,Fg, FirstOrSecond, FM_nod, FK_nod) Real(FEKi), DIMENSION (3,3), INTENT(IN) :: DIRCOS !direction cosice matrix (global to local) (3x3) Real(FEKi), DIMENSION (12,12), INTENT(IN) :: Me,Ke !element M and K matrices (12x12) in GLOBAL REFERENCE (DIRCOS^T K DIRCOS) Real(ReKi), DIMENSION (12), INTENT(IN) :: Udotdot, Y2 ! nodal accelerations and nodal displacements - Real(FEKi), DIMENSION (12), INTENT(IN) :: Fg ! element load vector from gravity and initial cable pretension + Real(FEKi), DIMENSION (12), INTENT(IN) :: Fg ! element load vector from gravity and initial cable pretension (orientation dependent for floating, constant for fixed-bottom) Integer(IntKi), INTENT(IN) :: FirstOrSecond !1 or 2 depending on node of interest REAL(ReKi), DIMENSION (6), INTENT(OUT) :: FM_nod, FK_nod !output static and dynamic forces and moments !Locals From c87830862d2aa008ec1f2d6e43141bac6468642b Mon Sep 17 00:00:00 2001 From: Roger Bergua Date: Sat, 6 Jun 2026 10:04:38 -0600 Subject: [PATCH 6/9] SubDyn floating: refactored self-weight output calculation --- modules/subdyn/src/SubDyn_Output.f90 | 17 +++++++++-------- 1 file changed, 9 insertions(+), 8 deletions(-) diff --git a/modules/subdyn/src/SubDyn_Output.f90 b/modules/subdyn/src/SubDyn_Output.f90 index 5da5442f0..5a0842ae8 100644 --- a/modules/subdyn/src/SubDyn_Output.f90 +++ b/modules/subdyn/src/SubDyn_Output.f90 @@ -501,7 +501,6 @@ subroutine ElementForce(pLst, iiNode, JJ, FM_elm, FK_elm, sgn, DIRCOS, bUseInput real(ReKi) , dimension(12) :: X_e, Xdd_e ! Displacement and acceleration for an element real(FEKi) , dimension(12) :: Fg_e ! Gravity + initial cable pretension load vector (reorient for floating, constant for fixed-bottom) real(FEKi) , dimension(3,3) :: CurDirCos ! Current element direction cosine matrix in the floating body frame - real(R8Ki) , dimension(3,3) :: Rb2g ! Rotation matrix from body/Guyan to global frame (transpose of Rg2b) integer(IntKi), dimension(2), parameter :: NodeNumber_To_Sign = (/-1, +1/) iElem = pLst%ElmIDs(iiNode,JJ) ! element number @@ -513,21 +512,23 @@ subroutine ElementForce(pLst, iiNode, JJ, FM_elm, FK_elm, sgn, DIRCOS, bUseInput Xdd_e(1:6) = matmul(RRg2b,m%U_full_dotdot(p%NodesDOF(ElemNodes(1))%List(1:6))) ! Transform acceleration to be back in the Guyan frame Xdd_e(7:12) = matmul(RRg2b,m%U_full_dotdot(p%NodesDOF(ElemNodes(2))%List(1:6))) ! Load the computed at initialization gravity + initial cable pretension load vector. - ! For floating systems, project the load vector into the current body/Guyan frame for Fx/Fy/Fz. + ! For floating systems, update the force components (Fx/Fy/Fz) from the initial body/Guyan frame to the current one. ! For the self-weight bending moments, recompute them based on the current element direction cosine matrix. Fg_e = real(pLst%Fg(:,iiNode,JJ), R8Ki) if (p%Floating) then - Fg_e(1:6) = matmul(RRg2b, Fg_e(1:6)) ! project node-1 forces/moments into body frame - Fg_e(7:12) = matmul(RRg2b, Fg_e(7:12)) ! project node-2 forces/moments into body frame + ! Rotate only the force (Fx,Fy,Fz) components into the body frame + Fg_e(1:3) = matmul(Rg2b, Fg_e(1:3)) + Fg_e(7:9) = matmul(Rg2b, Fg_e(7:9)) - Rb2g = transpose(Rg2b) ! current body-to-global rotation needed for moment projection - CurDirCos = matmul(Rb2g, p%ElemProps(iElem)%DirCos) - - ! Element self-weight bending moment contribution at the nodes: + ! Recompute bending moment components from current element orientation + ! CurDirCos = Rb2g * DirCos0 (DirCos0 is the element direction cosine matrix at initialization) + CurDirCos = matmul(transpose(Rg2b), p%ElemProps(iElem)%DirCos) Fg_e(4) = -p%ElemProps(iElem)%Length**2 * p%ElemProps(iElem)%Rho * p%ElemProps(iElem)%Area * p%g / 12.0_FEKi * CurDirCos(2,3) Fg_e(5) = p%ElemProps(iElem)%Length**2 * p%ElemProps(iElem)%Rho * p%ElemProps(iElem)%Area * p%g / 12.0_FEKi * CurDirCos(1,3) + Fg_e(6) = 0.0_FEKi ! no torsional self-weight moment Fg_e(10) = -Fg_e(4) Fg_e(11) = -Fg_e(5) + Fg_e(12) = 0.0_FEKi endif if (.not. bUseInputDirCos) then DIRCOS=transpose(p%ElemProps(iElem)%DirCos)! global to local From 5d55d8a74e07c46a3027d44886583a1af86f4b1d Mon Sep 17 00:00:00 2001 From: Roger Bergua Date: Sat, 6 Jun 2026 10:35:53 -0600 Subject: [PATCH 7/9] SubDyn_Output.f90: clean up trailing whitespace --- modules/subdyn/src/SubDyn_Output.f90 | 232 +++++++++++++-------------- 1 file changed, 116 insertions(+), 116 deletions(-) diff --git a/modules/subdyn/src/SubDyn_Output.f90 b/modules/subdyn/src/SubDyn_Output.f90 index 5a0842ae8..5ab01acc5 100644 --- a/modules/subdyn/src/SubDyn_Output.f90 +++ b/modules/subdyn/src/SubDyn_Output.f90 @@ -63,24 +63,24 @@ SUBROUTINE SDOut_Init( Init, y, p, misc, InitOut, WtrDpth, ErrStat, ErrMsg ) INTEGER(IntKi) :: iElem ! Index of element in Element List INTEGER(IntKi) :: iNode ! Index of node in Node list INTEGER(IntKi) :: iiElem ! Loop counter on element index - INTEGER(IntKi) :: nElemPerNode, nNodesPerElem ! Number of elements connecting to a node, Number of nodes per elem + INTEGER(IntKi) :: nElemPerNode, nNodesPerElem ! Number of elements connecting to a node, Number of nodes per elem type(MeshAuxDataType), pointer :: pLst !< Alias to shorten notation and highlight code similarities real(ReKi), allocatable :: T_TIreact(:,:) ! Transpose of TIreact, temporary - ErrStat = 0 + ErrStat = 0 ErrMsg="" p%OutAllDims=12*p%NMembers*2 !size of AllOut Member Joint forces - ! Check that the variables in OutList are valid + ! Check that the variables in OutList are valid CALL SDOut_ChkOutLst( Init%SSOutList, p, ErrStat2, ErrMsg2 ); if(Failed()) return ! --- Allocation (size 0 if not outputs) - !IF ( ALLOCATED( p%OutParam ) .AND. p%NumOuts > 0 ) THEN ! Output has been requested + !IF ( ALLOCATED( p%OutParam ) .AND. p%NumOuts > 0 ) THEN ! Output has been requested ! Allocate SDWrOuput which is used to store a time step's worth of output channels, prior to writing to a file. CALL AllocAry(misc%SDWrOutput , p%NumOuts + p%OutAllInt*p%OutAllDims, 'SDWrOutupt' , ErrStat2, ErrMsg2) ; if(Failed()) return - ! Allocate WriteOuput + ! Allocate WriteOuput CALL AllocAry(y%WriteOutput , p%NumOuts + p%OutAllInt*p%OutAllDims, 'WriteOutput', ErrStat2, ErrMsg2); if(Failed()) return - allocate(misc%AllOuts(0:MaxOutPts + p%OutAllInt*p%OutAllDims)) ! Need to start at 0... + allocate(misc%AllOuts(0:MaxOutPts + p%OutAllInt*p%OutAllDims)) ! Need to start at 0... ! Header, and Units, copy of data already available in the OutParam data structure ! TODO TODO TODO remove copy CALL AllocAry(InitOut%WriteOutputHdr, p%NumOuts + p%OutAllint*p%OutAllDims, 'WriteOutputHdr', ErrStat2, ErrMsg2); if(Failed()) return CALL AllocAry(InitOut%WriteOutputUnt, p%NumOuts + p%OutAllint*p%OutAllDims, 'WriteOutputUnt', ErrStat2, ErrMsg2); if(Failed()) return @@ -90,9 +90,9 @@ SUBROUTINE SDOut_Init( Init, y, p, misc, InitOut, WtrDpth, ErrStat, ErrMsg ) y%WriteOutput = 0 DO I = 1,p%NumOuts+p%OutAllint*p%OutAllDims InitOut%WriteOutputHdr(I) = TRIM( p%OutParam(I)%Name ) - InitOut%WriteOutputUnt(I) = TRIM( p%OutParam(I)%Units ) - END DO - + InitOut%WriteOutputUnt(I) = TRIM( p%OutParam(I)%Units ) + END DO + !_________________________________ OUTPUT FOR REQUESTED MEMBERS _______________________________ DO I=1,p%NMOutputs pLst => p%MOutLst(I) ! Alias to shorten notations @@ -120,13 +120,13 @@ SUBROUTINE SDOut_Init( Init, y, p, misc, InitOut, WtrDpth, ErrStat, ErrMsg ) IF (K2 == 2) EXIT ! we found both elements already, error... K2=K2+1 call ConfigOutputNode_MKF_ID(pLst, iElem, iiNode=J, iStore=K2, NodeID2=iNode) - END IF + END IF ENDDO ! iiElem, nElemPerNode ENDDO !J, Noutcnt ENDDO !I, NMOutputs - + !_________________________________ OUTPUT FOR ALL MEMBERS __________________________________ - IF (p%OutAll) THEN !I need to store all member end forces and moments + IF (p%OutAll) THEN !I need to store all member end forces and moments ! MOutLst2: nodal output info by members, for all members, First and Last Node ALLOCATE ( p%MOutLst2(p%NMembers), STAT = ErrStat2 ); ErrMsg2 = 'Error allocating p%MOutLst2 array in SDOut_Init'; if(Failed()) return @@ -134,7 +134,7 @@ SUBROUTINE SDOut_Init( Init, y, p, misc, InitOut, WtrDpth, ErrStat, ErrMsg ) DO iMember=1,p%NMembers pLst => p%MOutLst2(iMember) ! Alias pLst%MemberID = Init%Members(iMember,1) - nNodesPerElem = count(Init%MemberNodes(iMember,:) >0 ) + nNodesPerElem = count(Init%MemberNodes(iMember,:) >0 ) CALL AllocAry(pLst%NodeIDs, nNodesPerElem, 'MOutLst2(I)%NodeIDs', ErrStat2, ErrMsg2); if(Failed()) return CALL AllocAry(pLst%ElmIDs, 2, 1, 'MOutLst2(I)%ElmIDs' , ErrStat2, ErrMsg2); if(Failed()) return CALL AllocAry(pLst%ElmNds, 2, 1, 'MOutLst2(I)%ElmNds' , ErrStat2, ErrMsg2); if(Failed()) return @@ -163,7 +163,7 @@ SUBROUTINE SDOut_Init( Init, y, p, misc, InitOut, WtrDpth, ErrStat, ErrMsg ) p%OutReact = .FALSE. DO I=1,p%NumOuts if ( ANY( p%OutParam(I)%Indx == ReactSS) ) THEN ! bjj: removed check of first 5 characters being "React" because (1) cases matter and (2) we can also ask for "-React*" or "mREACT" - p%OutReact =.TRUE. + p%OutReact =.TRUE. EXIT ENDIF ENDDO @@ -181,8 +181,8 @@ SUBROUTINE SDOut_Init( Init, y, p, misc, InitOut, WtrDpth, ErrStat, ErrMsg ) CALL AllocAry(pLst%Ke, 12, 12 , 1, nElemPerNode, ' p%MOutLst3(I)%Ke' , ErrStat2, ErrMsg2); if(Failed()) return CALL AllocAry(pLst%Fg, 12 , 1, nElemPerNode, ' p%MOutLst3(I)%Fg' , ErrStat2, ErrMsg2); if(Failed()) return DO iiElem = 1, nElemPerNode - iElem = Init%NodesConnE(iNode, iiElem+1) ! iiElem-th Element Number in the set of elements attached to the selected node - call ConfigOutputNode_MKF_ID(pLst, iElem, iiNode=1, iStore=iiElem, NodeID2=iNode) + iElem = Init%NodesConnE(iNode, iiElem+1) ! iiElem-th Element Number in the set of elements attached to the selected node + call ConfigOutputNode_MKF_ID(pLst, iElem, iiNode=1, iStore=iiElem, NodeID2=iNode) ENDDO ENDDO ! Compute p%TIreact, rigid transf. matrix from reaction DOFs to base structure point (0,0,-WD) @@ -196,13 +196,13 @@ SUBROUTINE SDOut_Init( Init, y, p, misc, InitOut, WtrDpth, ErrStat, ErrMsg ) CONTAINS LOGICAL FUNCTION Failed() - call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, 'SDOut_Init') + call SetErrStat(ErrStat2, ErrMsg2, ErrStat, ErrMsg, 'SDOut_Init') Failed = ErrStat >= AbortErrLev END FUNCTION Failed !> Returns true if an element is connected to node iNode, and along member iMember LOGICAL FUNCTION ThisElementIsAlongMember(iElem, iNode, iMember) - integer(IntKi), intent(in) :: iElem !< Element index + integer(IntKi), intent(in) :: iElem !< Element index integer(IntKi), intent(in) :: iNode !< Node index integer(IntKi), intent(in) :: iMember !< Member index integer(IntKi), dimension(2) :: ElemNodes ! Node IDs for element under consideration (may not be consecutive numbers) @@ -215,7 +215,7 @@ LOGICAL FUNCTION ThisElementIsAlongMember(iElem, iNode, iMember) iOtherNode=ElemNodes(1) else ThisElementIsAlongMember=.false. ! Not along member since nodes don't match - return + return endif ! Being along the member means the second node of the element is in the node list of the member ThisElementIsAlongMember= ANY(Init%MemberNodes(iMember,:) == iOtherNode) @@ -224,10 +224,10 @@ LOGICAL FUNCTION ThisElementIsAlongMember(iElem, iNode, iMember) !> Set different "data" for a given output node, and possibly store more than one "data" per node: !! The "data" is: !! - Mass, stiffness matrices and constant element force (gravity and cable) - !! - A flag whether the node is the 1st or second node of an element + !! - A flag whether the node is the 1st or second node of an element !! The "data" is stored at the index (iiNode,iStore): !! - iiNode: node index within the list of nodes that are to be used for output for this member - !! - iStore: index over the number of "data" stored per node. E.g. Member1 and 2 connecting to a node + !! - iStore: index over the number of "data" stored per node. E.g. Member1 and 2 connecting to a node SUBROUTINE ConfigOutputNode_MKF_ID(pLst, iElem, iiNode, iStore, NodeID2) type(MeshAuxDataType), intent(inout) :: pLst !< Info for one member output integer(IntKi) , intent(in) :: iElem !< Element index to which the node belong @@ -238,8 +238,8 @@ SUBROUTINE ConfigOutputNode_MKF_ID(pLst, iElem, iiNode, iStore, NodeID2) REAL(FEKi) :: FCe(12) ! Pretension force from cable element pLst%ElmIDs(iiNode,iStore) = iElem ! This array has for each joint requested the elements' ID to get results for ElemNodes = p%Elems(iElem,2:3) ! 1st and 2nd node of the k-th element - if (ElemNodes(2) == NodeID2) then - pLst%ElmNds(iiNode,iStore) = 2 ! store whether first or second node of element + if (ElemNodes(2) == NodeID2) then + pLst%ElmNds(iiNode,iStore) = 2 ! store whether first or second node of element else pLst%ElmNds(iiNode,iStore) = 1 ! store whether first or second node of element endif @@ -259,7 +259,7 @@ END SUBROUTINE SDOut_Init !> Writes the data stored in the y variable to the correct indexed postions in WriteOutput !! This is called by SD_CalcOutput() at each time step. !! This routine does fill Allouts -!! note that this routine assumes m%u_TP and m%udotdot_TP have been set before calling +!! note that this routine assumes m%u_TP and m%udotdot_TP have been set before calling !! this routine (which is done in SD_CalcOutput() and SD CalcContStateDeriv) SUBROUTINE SDOut_MapOutputs(u,p,x, y, m, AllOuts, ErrStat, ErrMsg ) type(SD_InputType), intent( in ) :: u ! SubDyn module's input data @@ -287,7 +287,7 @@ SUBROUTINE SDOut_MapOutputs(u,p,x, y, m, AllOuts, ErrStat, ErrMsg ) INTEGER(IntKi) :: ErrStat2 ! Error status of the operation CHARACTER(ErrMsgLen) :: ErrMsg2 ! Error message if ErrStat /= ErrID_None - ErrStat = ErrID_None + ErrStat = ErrID_None ErrMsg = "" if ( p%Floating ) then @@ -307,18 +307,18 @@ SUBROUTINE SDOut_MapOutputs(u,p,x, y, m, AllOuts, ErrStat, ErrMsg ) RRb2g = transpose(RRg2b) AllOuts = 0.0_ReKi ! initialize for those outputs that aren't valid (and thus aren't set in this routine) - + ! -------------------------------------------------------------------------------- ! --- Requested member-outputs (Node kinematics and loads) ! -------------------------------------------------------------------------------- ! p%MOutLst has the mapping for the member, node, elements per node, to be used - ! MXNYZZZ will need to connects to p%MOutLst(X)%ElmIDs(Y,1:2) if it is a force or accel; else to u%UFL(p%MOutLst(X)%NodeIDs(Y)) + ! MXNYZZZ will need to connects to p%MOutLst(X)%ElmIDs(Y,1:2) if it is a force or accel; else to u%UFL(p%MOutLst(X)%NodeIDs(Y)) if (p%NumOuts > 0) then !bjj: some of these fields aren't allocated when NumOuts==0 ! Loop over member-outputs requested DO iMemberOutput=1,p%NMOutputs - pLst=>p%MOutLst(iMemberOutput) ! List for a given member-output - DO iiNode=1,pLst%NOutCnt !Iterate on requested nodes for that member - ! --- Forces (potentially averaged on 2 elements) + pLst=>p%MOutLst(iMemberOutput) ! List for a given member-output + DO iiNode=1,pLst%NOutCnt !Iterate on requested nodes for that member + ! --- Forces (potentially averaged on 2 elements) call ElementForce(pLst, iiNode, 1, FM_elm, FK_elm, sgn, DIRCOS, .false.) FM_elm2=sgn*FM_elm FK_elm2=sgn*FK_elm @@ -351,15 +351,15 @@ SUBROUTINE SDOut_MapOutputs(u,p,x, y, m, AllOuts, ErrStat, ErrMsg ) ENDDO ! iiNode, Loop on requested nodes for that member ENDDO ! iMemberOutput, Loop on member outputs END IF - + ! -------------------------------------------------------------------------------- - ! --- All nodal loads from stiffness and mass matrix + ! --- All nodal loads from stiffness and mass matrix ! -------------------------------------------------------------------------------- ! "MaaaJbFKxe, MaaaJbMKxe MaaaJbFMxe, MaaaJbMMxe for member aaa and node b." - IF (p%OutAll) THEN + IF (p%OutAll) THEN DO iMemberOutput=1,p%NMembers !Cycle on all members pLst=>p%MOutLst2(iMemberOutput) - DO iiNode=1,2 !Iterate on requested nodes for that member (first and last) + DO iiNode=1,2 !Iterate on requested nodes for that member (first and last) call ElementForce(pLst, iiNode, 1, FM_elm, FK_elm, sgn, DIRCOS, .false.) ! Store in All Outs L = MaxOutPts+(iMemberOutput-1)*24+(iiNode-1)*12+1 @@ -368,7 +368,7 @@ SUBROUTINE SDOut_MapOutputs(u,p,x, y, m, AllOuts, ErrStat, ErrMsg ) ENDDO !iiNode, nodes 1 and 2 ENDDO ! iMemberOutput, Loop on members ENDIF - + ! -------------------------------------------------------------------------------- ! --- Interface kinematics and loads (TP/platform reference point) ! -------------------------------------------------------------------------------- @@ -384,7 +384,7 @@ SUBROUTINE SDOut_MapOutputs(u,p,x, y, m, AllOuts, ErrStat, ErrMsg ) AllOuts(IntfSS(1:6,iTP)) = - (/y%Y1Mesh(iTP)%Force(:,1), y%Y1Mesh(iTP)%Moment(:,1)/) !-y%Y1 !Note this is the force that the TP applies to the Jacket, opposite to what the GLue Code needs thus "-" sign end do - ! Interface translations and rotations in SS coordinate system + ! Interface translations and rotations in SS coordinate system ! "IntfTDXss, IntfTDYss, IntfTDZss, IntfRDXss, IntfRDYss IntfRDZss" do iTP = 1,nTP AllOuts(IntfTRss(1:3,iTP)) = u%TPMesh(iTP)%TranslationDisp(:,1) @@ -431,13 +431,13 @@ SUBROUTINE SDOut_MapOutputs(u,p,x, y, m, AllOuts, ErrStat, ErrMsg ) ! --- Modal parameters "SSqmXX, SSqmdotXX, SSqmddXX" amplitude, speed and acceleration ! -------------------------------------------------------------------------------- maxOutModes = min(p%nDOFM,99) ! We only have space for the first 99 values - IF ( maxOutModes > 0 ) THEN + IF ( maxOutModes > 0 ) THEN !BJJ: TODO: is there a check to see if we requested these channels but didn't request the modes? (i.e., retain 2 modes but asked for 75th mode?) AllOuts(SSqm01 :SSqm01 +maxOutModes-1) = x%qm (1:maxOutModes) AllOuts(SSqmd01 :SSqmd01 +maxOutModes-1) = x%qmdot (1:maxOutModes) AllOuts(SSqmdd01:SSqmdd01+maxOutModes-1) = m%qmdotdot(1:maxOutModes) END IF - + ! --------------------------------------------------------------------------------} ! --- Base reaction loads ! --------------------------------------------------------------------------------{ @@ -533,7 +533,7 @@ subroutine ElementForce(pLst, iiNode, JJ, FM_elm, FK_elm, sgn, DIRCOS, bUseInput if (.not. bUseInputDirCos) then DIRCOS=transpose(p%ElemProps(iElem)%DirCos)! global to local endif - CALL CALC_NODE_FORCES( DIRCOS, pLst%Me(:,:,iiNode,JJ),pLst%Ke(:,:,iiNode,JJ), Xdd_e, X_e, Fg_e, FirstOrSecond, FM_elm, FK_elm) + CALL CALC_NODE_FORCES( DIRCOS, pLst%Me(:,:,iiNode,JJ),pLst%Ke(:,:,iiNode,JJ), Xdd_e, X_e, Fg_e, FirstOrSecond, FM_elm, FK_elm) end subroutine ElementForce !==================================================================================================== @@ -551,30 +551,30 @@ SUBROUTINE CALC_NODE_FORCES(DIRCOS,Me,Ke,Udotdot,Y2 ,Fg, FirstOrSecond, FM_nod, REAL(ReKi), DIMENSION (6), INTENT(OUT) :: FM_nod, FK_nod !output static and dynamic forces and moments !Locals INTEGER(IntKi) :: L !counter - REAL(DbKi), DIMENSION(12) :: FM_glb, FF_glb, FM_elm, FF_elm ! temporary storage + REAL(DbKi), DIMENSION(12) :: FM_glb, FF_glb, FM_elm, FF_elm ! temporary storage FM_glb = matmul(Me,Udotdot) ! GLOBAL REFERENCE FF_glb = matmul(Ke,Y2) ! GLOBAL REFERENCE FF_glb = FF_glb - Fg ! GLOBAL REFERENCE DO L=1,4 ! Transforming coordinates 3 at a time FM_elm((L-1)*3+1:L*3) = matmul(DIRCOS, FM_glb( (L-1)*3+1:L*3 ) ) - FF_elm((L-1)*3+1:L*3) = matmul(DIRCOS, FF_glb( (L-1)*3+1:L*3 ) ) + FF_elm((L-1)*3+1:L*3) = matmul(DIRCOS, FF_glb( (L-1)*3+1:L*3 ) ) ENDDO - FM_nod = FM_elm(6*(FirstOrSecond-1)+1:FirstOrSecond*6) ! k2=1, 1:6, k2=2 7:12 - FK_nod = FF_elm(6*(FirstOrSecond-1)+1:FirstOrSecond*6) + FM_nod = FM_elm(6*(FirstOrSecond-1)+1:FirstOrSecond*6) ! k2=1, 1:6, k2=2 7:12 + FK_nod = FF_elm(6*(FirstOrSecond-1)+1:FirstOrSecond*6) - END SUBROUTINE CALC_NODE_FORCES + END SUBROUTINE CALC_NODE_FORCES END SUBROUTINE SDOut_MapOutputs !==================================================================================================== SUBROUTINE SDOut_CloseSum( UnSum, ErrStat, ErrMsg ) INTEGER, INTENT( IN ) :: UnSum ! the unit number for the SubDyn summary file - INTEGER, INTENT( OUT ) :: ErrStat ! returns a non-zero value when an error occurs + INTEGER, INTENT( OUT ) :: ErrStat ! returns a non-zero value when an error occurs CHARACTER(*), INTENT( OUT ) :: ErrMsg ! Error message if ErrStat /= ErrID_None ! Local variables - INTEGER :: Stat ! status from I/) operation - ErrStat = ErrID_None + INTEGER :: Stat ! status from I/) operation + ErrStat = ErrID_None ErrMsg = "" ! Write any closing information in the summary file IF ( UnSum > 0 ) THEN @@ -590,31 +590,31 @@ SUBROUTINE SDOut_CloseSum( UnSum, ErrStat, ErrMsg ) ErrMsg = TRIM(ErrMsg)//' Problem closing summary file.' END IF IF ( ErrStat /= ErrID_None ) ErrMsg = 'SDOut_CloseSum'//TRIM(ErrMsg) - END IF -END SUBROUTINE SDOut_CloseSum + END IF +END SUBROUTINE SDOut_CloseSum !==================================================================================================== SUBROUTINE SDOut_OpenSum( UnSum, SummaryName, SD_Prog, ErrStat, ErrMsg ) - INTEGER, INTENT( OUT ) :: UnSum ! the unit number for the SubDyn summary file + INTEGER, INTENT( OUT ) :: UnSum ! the unit number for the SubDyn summary file CHARACTER(*), INTENT( IN ) :: SummaryName ! the name of the SubDyn summary file TYPE(ProgDesc), INTENT( IN ) :: SD_Prog ! the name/version/date of the program - INTEGER, INTENT( OUT ) :: ErrStat ! returns a non-zero value when an error occurs + INTEGER, INTENT( OUT ) :: ErrStat ! returns a non-zero value when an error occurs CHARACTER(*), INTENT( OUT ) :: ErrMsg ! Error message if ErrStat /= ErrID_None integer :: ErrStat2 - ErrStat = ErrID_None - ErrMsg = "" + ErrStat = ErrID_None + ErrMsg = "" CALL GetNewUnit( UnSum ) - CALL OpenFOutFile ( UnSum, SummaryName, ErrStat, ErrMsg ) + CALL OpenFOutFile ( UnSum, SummaryName, ErrStat, ErrMsg ) IF ( ErrStat >= AbortErrLev ) THEN ErrMsg = 'Failed to open SubDyn summary file: '//TRIM(ErrMsg) RETURN END IF - + ! Write the summary file header WRITE (UnSum,'(/,A/)', IOSTAT=ErrStat2) '#This summary file was generated by '//TRIM( SD_Prog%Name )//& ' '//TRIM( SD_Prog%Ver )//' on '//CurDate()//' at '//CurTime()//'.' -END SUBROUTINE SDOut_OpenSum +END SUBROUTINE SDOut_OpenSum !==================================================================================================== SUBROUTINE SDOut_OpenOutput( ProgVer, OutRootName, p, InitOut, ErrStat, ErrMsg ) @@ -624,17 +624,17 @@ SUBROUTINE SDOut_OpenOutput( ProgVer, OutRootName, p, InitOut, ErrStat, ErrMsg ! Passed variables TYPE(ProgDesc), INTENT( IN ) :: ProgVer CHARACTER(*), INTENT( IN ) :: OutRootName ! Root name for the output file - TYPE(SD_ParameterType), INTENT( INOUT ) :: p + TYPE(SD_ParameterType), INTENT( INOUT ) :: p TYPE(SD_InitOutPutType ), INTENT( IN ) :: InitOut ! - INTEGER, INTENT( OUT ) :: ErrStat ! a non-zero value indicates an error occurred + INTEGER, INTENT( OUT ) :: ErrStat ! a non-zero value indicates an error occurred CHARACTER(*), INTENT( OUT ) :: ErrMsg ! Error message if ErrStat /= ErrID_None ! Local variables - INTEGER :: I ! Generic loop counter + INTEGER :: I ! Generic loop counter CHARACTER(1024) :: OutFileName ! The name of the output file including the full path. CHARACTER(200) :: Frmt ! a string to hold a format statement - INTEGER :: ErrStat2 + INTEGER :: ErrStat2 - ErrStat = ErrID_None + ErrStat = ErrID_None ErrMsg = "" ! Initialize to -1 to indicate that the output file unit is not valid @@ -645,28 +645,28 @@ SUBROUTINE SDOut_OpenOutput( ProgVer, OutRootName, p, InitOut, ErrStat, ErrMsg call WrScr('SubDyn: no outputs were requested, so separate output file will not be generated.') return end if - + ! Open the file for output OutFileName = TRIM(OutRootName)//'.out' call GetNewUnit( p%UnJckF ) - call OpenFOutFile ( p%UnJckF, OutFileName, ErrStat, ErrMsg ) + call OpenFOutFile ( p%UnJckF, OutFileName, ErrStat, ErrMsg ) if (ErrStat >= AbortErrLev) then ErrMsg = ' Error opening SubDyn-level output file: '//TRIM(ErrMsg) return end if - + ! Write the output file header write(p%UnJckF,'(/,A/)', IOSTAT=ErrStat2) 'These predictions were generated by '//TRIM(GETNVD(ProgVer))//& ' on '//CurDate()//' at '//CurTime()//'.' - + write(p%UnJckF, '(//)') ! add 3 lines to make file format consistant with FAST v8 (headers on line 7; units on line 8) [this allows easier post-processing] - + ! Write the names of the output parameters: Frmt = '(A8,'//TRIM(Int2LStr(p%NumOuts+p%OutAllInt*p%OutAllDims))//'(:,A,'//TRIM( p%OutSFmt )//'))' write(p%UnJckF,Frmt, IOSTAT=ErrStat2) TRIM( 'Time' ), ( p%Delim, TRIM( InitOut%WriteOutputHdr(I) ), I=1,p%NumOuts+p%OutAllInt*p%OutAllDims ) - - ! Write the units of the output parameters: + + ! Write the units of the output parameters: write(p%UnJckF,Frmt, IOSTAT=ErrStat2) TRIM( 's'), ( p%Delim, TRIM( InitOut%WriteOutputUnt(I) ), I=1,p%NumOuts+p%OutAllInt*p%OutAllDims ) END SUBROUTINE SDOut_OpenOutput @@ -702,19 +702,19 @@ SUBROUTINE SDOut_WriteOutputNames( UnJckF, p, ErrStat, ErrMsg ) INTEGER, INTENT( IN ) :: UnJckF ! file unit for the output file TYPE(SD_ParameterType), INTENT( IN ) :: p ! SubDyn module's parameter data - INTEGER, INTENT( OUT ) :: ErrStat ! returns a non-zero value when an error occurs + INTEGER, INTENT( OUT ) :: ErrStat ! returns a non-zero value when an error occurs CHARACTER(*), INTENT( OUT ) :: ErrMsg ! Error message if ErrStat /= ErrID_None - + CHARACTER(200) :: Frmt ! a string to hold a format statement INTEGER :: I ! Generic loop counter - - ErrStat = ErrID_None + + ErrStat = ErrID_None ErrMsg = "" - + Frmt = '(A8,'//TRIM(Int2LStr(p%NumOuts+p%OutAllInt*p%OutAllDims))//'(:,A,'//TRIM( p%OutSFmt )//'))' WRITE(UnJckF,Frmt) TRIM( p%OutParam(0)%Name ), ( p%Delim, TRIM( p%OutParam(I)%Name ), I=1,p%NumOuts+p%OutAllInt*p%OutAllDims ) - + END SUBROUTINE SDOut_WriteOutputNames !==================================================================================================== @@ -722,35 +722,35 @@ END SUBROUTINE SDOut_WriteOutputNames SUBROUTINE SDOut_WriteOutputUnits( UnJckF, p, ErrStat, ErrMsg ) INTEGER, INTENT( IN ) :: UnJckF ! file unit for the output file TYPE(SD_ParameterType), INTENT( IN ) :: p ! SubDyn module's parameter data - INTEGER, INTENT( OUT ) :: ErrStat ! returns a non-zero value when an error occurs + INTEGER, INTENT( OUT ) :: ErrStat ! returns a non-zero value when an error occurs CHARACTER(*), INTENT( OUT ) :: ErrMsg ! Error message if ErrStat /= ErrID_None CHARACTER(200) :: Frmt ! a string to hold a format statement INTEGER :: I ! Generic loop counter - ErrStat = ErrID_None + ErrStat = ErrID_None ErrMsg = "" - + Frmt = '(A8,'//TRIM(Int2LStr(p%NumOuts+p%OutAllInt*p%OutAllDims))//'(:,A,'//TRIM( p%OutSFmt )//'))' WRITE(UnJckF,Frmt) TRIM( p%OutParam(0)%Units ), ( p%Delim, TRIM( p%OutParam(I)%Units ), I=1,p%NumOuts+p%OutAllInt*p%OutAllDims ) - + END SUBROUTINE SDOut_WriteOutputUnits !==================================================================================================== SUBROUTINE SDOut_WriteOutputs( UnJckF, Time, SDWrOutput, p, ErrStat, ErrMsg ) ! This subroutine writes the data stored in WriteOutputs (and indexed in OutParam) to the file ! opened in SDOut_Init() -!---------------------------------------------------------------------------------------------------- +!---------------------------------------------------------------------------------------------------- INTEGER, INTENT( IN ) :: UnJckF ! file unit for the output file REAL(DbKi), INTENT( IN ) :: Time ! Time for this output REAL(ReKi), INTENT( IN ) :: SDWrOutput(:) ! SubDyn module's output data TYPE(SD_ParameterType), INTENT( IN ) :: p ! SubDyn module's parameter data - INTEGER, INTENT( OUT ) :: ErrStat ! returns a non-zero value when an error occurs + INTEGER, INTENT( OUT ) :: ErrStat ! returns a non-zero value when an error occurs CHARACTER(*), INTENT( OUT ) :: ErrMsg ! Error message if ErrStat /= ErrID_None ! Local variables INTEGER :: I ! Generic loop counter CHARACTER(200) :: Frmt ! a string to hold a format statement - ErrStat = ErrID_None + ErrStat = ErrID_None ErrMsg = "" ! If output file is not open, return @@ -769,13 +769,13 @@ END SUBROUTINE SDOut_WriteOutputs !==================================================================================================== SUBROUTINE SDOut_ChkOutLst( OutList, p, ErrStat, ErrMsg ) ! This routine checks the names of inputted output channels, checks to see if any of them are ill- -! conditioned (returning an error if so), and assigns the OutputDataType settings (i.e, the index, -! name, and units of the output channels). +! conditioned (returning an error if so), and assigns the OutputDataType settings (i.e, the index, +! name, and units of the output channels). ! NOTE OutParam is populated here -!---------------------------------------------------------------------------------------------------- +!---------------------------------------------------------------------------------------------------- TYPE(SD_ParameterType), INTENT( INOUT ) :: p ! SubDyn module parameter data - CHARACTER(ChanLen), INTENT( IN ) :: OutList (:) ! An array holding the names of the requested output channels. - INTEGER, INTENT( OUT ) :: ErrStat ! a non-zero value indicates an error occurred + CHARACTER(ChanLen), INTENT( IN ) :: OutList (:) ! An array holding the names of the requested output channels. + INTEGER, INTENT( OUT ) :: ErrStat ! a non-zero value indicates an error occurred CHARACTER(*), INTENT( OUT ) :: ErrMsg ! Error message if ErrStat /= ErrID_None ! Local variables. INTEGER :: I,J,K ! Generic loop-counting index. @@ -785,9 +785,9 @@ SUBROUTINE SDOut_ChkOutLst( OutList, p, ErrStat, ErrMsg ) CHARACTER(ChanLen), DIMENSION(24) :: ToTUnits,ToTNames,ToTNames0 LOGICAL :: InvalidOutput(0:MaxOutPts) ! This array determines if the output channel is valid for this configuration LOGICAL :: CheckOutListAgain - ErrStat = ErrID_None + ErrStat = ErrID_None ErrMsg = "" - + InvalidOutput = .FALSE. ! mark invalid output channels: @@ -796,7 +796,7 @@ SUBROUTINE SDOut_ChkOutLst( OutList, p, ErrStat, ErrMsg ) InvalidOutput(SSqmd01 +k-1) = .true. InvalidOutput(SSqmdd01+k-1) = .true. END DO - + DO I=1,99 !I know el # and whether it is 1st node or second node if (I <= p%NMOutputs) then @@ -804,8 +804,8 @@ SUBROUTINE SDOut_ChkOutLst( OutList, p, ErrStat, ErrMsg ) else INDX = 1 end if - - DO J=INDX,9 !Iterate on requested nodes for that member + + DO J=INDX,9 !Iterate on requested nodes for that member !Forces and moments InvalidOutput(MNfmKe (:,J,I)) = .true. !static forces and moments (6) Local Ref InvalidOutput(MNfmMe (:,J,I)) = .true. !dynamic forces and moments (6) Local Ref @@ -832,32 +832,32 @@ SUBROUTINE SDOut_ChkOutLst( OutList, p, ErrStat, ErrMsg ) InvalidOutput(IntfTRAss(:,I)) = .true. END DO END IF - + !------------------------------------------------------------------------------------------------- ! ALLOCATE the OutParam array - !------------------------------------------------------------------------------------------------- + !------------------------------------------------------------------------------------------------- ALLOCATE ( p%OutParam(1:p%NumOuts+p%OutAllInt*p%OutAllDims) , STAT=ErrStat ) IF ( ErrStat /= 0 ) THEN ErrMsg = ' Error allocating memory for the OutParam array.' ErrStat = ErrID_Fatal RETURN END IF - - + + !------------------------------------------------------------------------------------------------- ! Set index, name, and units for the output channels ! If a selected output channel is not available in this module, set error flag and return. - !------------------------------------------------------------------------------------------------- + !------------------------------------------------------------------------------------------------- !!!p%OutParam(0)%Name = 'Time' ! OutData(0) is the time channel by default. !!!p%OutParam(0)%Units = '(sec)' ! !!!p%OutParam(0)%Indx = Time !!!p%OutParam(0)%SignM = 1 - + DO I = 1,p%NumOuts - - p%OutParam(I)%Name = OutList(I) + + p%OutParam(I)%Name = OutList(I) OutListTmp = OutList(I) - + CALL Conv2UC( OutListTmp ) ! Convert OutListTmp to upper case ! Interface output backward compatibility @@ -868,9 +868,9 @@ SUBROUTINE SDOut_ChkOutLst( OutList, p, ErrStat, ErrMsg ) ! Reverse the sign (+/-) of the output channel if the user prefixed the ! channel name with a '-', '_', 'm', or 'M' character indicating "minus". - + CheckOutListAgain = .FALSE. - + IF ( INDEX( '-_', OutListTmp(1:1) ) > 0 ) THEN p%OutParam(I)%SignM = -1 ! ex, '-TipDxc1' causes the sign of TipDxc1 to be switched. OutListTmp = OutListTmp(2:) @@ -880,53 +880,53 @@ SUBROUTINE SDOut_ChkOutLst( OutList, p, ErrStat, ErrMsg ) ELSE p%OutParam(I)%SignM = 1 END IF - + if ( INDEX( 'mM', OutListTmp(1:1) ) > 0 .and. INDEX( '0123456789', OutListTmp(2:2) ) > 0 .and. INDEX( 'nN', OutListTmp(3:3) ) > 0 ) then ! an old-style output without the leading zero on the member number OutListTmp = OutListTmp(1:1)//'0'//OutListTmp(2:) CheckOutListAgain = .FALSE. end if Indx = IndexCharAry( OutListTmp(1:OutStrLenM1), ValidParamAry ) - - IF ( CheckOutListAgain .AND. Indx < 1 ) THEN ! Let's assume that "M" really meant "minus" and then test again + + IF ( CheckOutListAgain .AND. Indx < 1 ) THEN ! Let's assume that "M" really meant "minus" and then test again p%OutParam(I)%SignM = -1 ! ex, 'MTipDxc1' causes the sign of TipDxc1 to be switched. OutListTmp = OutListTmp(2:) - - Indx = IndexCharAry( OutListTmp(1:10), ValidParamAry ) + + Indx = IndexCharAry( OutListTmp(1:10), ValidParamAry ) END IF - + IF ( Indx > 0 ) THEN p%OutParam(I)%Indx = ParamIndxAry(Indx) IF ( InvalidOutput( ParamIndxAry(Indx) ) ) THEN - p%OutParam(I)%Units = 'INVALID' - p%OutParam(I)%SignM = 0 + p%OutParam(I)%Units = 'INVALID' + p%OutParam(I)%SignM = 0 ELSE p%OutParam(I)%Units = ParamUnitsAry(Indx) END IF ELSE ErrMsg = p%OutParam(I)%Name//' is not an available output channel.' ErrStat = ErrID_Fatal - p%OutParam(I)%Units = 'INVALID' + p%OutParam(I)%Units = 'INVALID' p%OutParam(I)%Indx = 0 p%OutParam(I)%SignM = 0 ! this will print all zeros END IF - + END DO - + IF (p%OutAll) THEN !Finish populating the OutParam with all the joint forces and moments ToTNames0=RESHAPE(SPREAD( (/"FKxe", "FKye", "FKze", "MKxe", "MKye", "MKze", "FMxe", "FMye", "FMze", "MMxe", "MMye", "MMze"/), 2, 2), (/24/) ) ToTUnits=RESHAPE(SPREAD( (/"(N) ","(N) ","(N) ", "(N*m)","(N*m)","(N*m)", "(N) ","(N) ","(N) ", "(N*m)","(N*m)","(N*m)"/), 2, 2), (/24/) ) DO I=1,p%NMembers DO K=1,2 - DO J=1,12 + DO J=1,12 TotNames(J+(K-1)*12)=TRIM("M"//Int2Lstr(I))//TRIM("J"//Int2Lstr(K))//TRIM(TotNames0(J)) - ENDDO + ENDDO ENDDO p%OutParam(p%NumOuts+(I-1)*12*2+1:p%NumOuts+I*12*2)%Name = ToTNames p%OutParam(p%NumOuts+(I-1)*12*2+1:p%NumOuts+I*12*2)%Units = ToTUnits ENDDO p%OutParam(p%NumOuts+1:p%NumOuts+p%OutAllDims)%SignM = 1 - p%OutParam(p%NumOuts+1:p%NumOuts+p%OutAllDims)%Indx= MaxOutPts+(/(J, J=1, p%OutAllDims)/) + p%OutParam(p%NumOuts+1:p%NumOuts+p%OutAllDims)%Indx= MaxOutPts+(/(J, J=1, p%OutAllDims)/) ENDIF END SUBROUTINE SDOut_ChkOutLst From 4e611c0df9dca849988ac1d1fb0bf334c92153b6 Mon Sep 17 00:00:00 2001 From: Roger Bergua Date: Mon, 8 Jun 2026 18:58:06 -0600 Subject: [PATCH 8/9] New SD regression test self-weight CTestList --- reg_tests/CTestList.cmake | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/reg_tests/CTestList.cmake b/reg_tests/CTestList.cmake index efe0587ca..4d2def7a8 100644 --- a/reg_tests/CTestList.cmake +++ b/reg_tests/CTestList.cmake @@ -56,11 +56,11 @@ function(regression TEST_SCRIPT EXECUTABLE SOURCE_DIRECTORY BUILD_DIRECTORY STEA if(STEADYSTATE_FLAG STREQUAL " ") set(STEADYSTATE_FLAG "") endif() - + if(OTHER_FLAGS STREQUAL " ") set(OTHER_FLAGS "") endif() - + add_test( ${TESTNAME} ${Python_EXECUTABLE} ${TEST_SCRIPT} @@ -111,7 +111,7 @@ function(of_fastlib_regression TESTNAME LABEL) regression(${TEST_SCRIPT} ${OPENFAST_EXECUTABLE} ${SOURCE_DIRECTORY} ${BUILD_DIRECTORY} " " "${TESTNAME}_fastlib" "${LABEL}" " " ${TESTNAME}) endfunction(of_fastlib_regression) -# openfast aeroacoustic +# openfast aeroacoustic function(of_regression_aeroacoustic TESTNAME LABEL) set(TEST_SCRIPT "${CMAKE_CURRENT_LIST_DIR}/executeOpenfastAeroAcousticRegressionCase.py") set(OPENFAST_EXECUTABLE "${CTEST_OPENFAST_EXECUTABLE}") @@ -514,7 +514,7 @@ sd_regression("SD_PendulumDamp" "subdyn;offshore") sd_regression("SD_Rigid" "subdyn;offshore") sd_regression("SD_SparHanging" "subdyn;offshore") sd_regression("SD_AnsysComp1_PinBeam" "subdyn;offshore") # TODO Issue #855 -sd_regression("SD_AnsysComp2_Cable" "subdyn;offshore") +sd_regression("SD_AnsysComp2_Cable" "subdyn;offshore") sd_regression("SD_AnsysComp3_PinBeamCable" "subdyn;offshore") # TODO Issue #855 sd_regression("SD_Spring_Case1" "subdyn;offshore") sd_regression("SD_Spring_Case2" "subdyn;offshore") @@ -523,6 +523,7 @@ sd_regression("SD_Revolute_Joint" "subdyn;offshore") sd_regression("SD_2Beam_Spring" "subdyn;offshore") sd_regression("SD_2Beam_Cantilever" "subdyn;offshore") sd_regression("SD_CantileverBeam_Rectangular" "subdyn;offshore") +sd_regression("SD_SelfWeight_FloatingSystem" "subdyn;offshore") # TODO test below are bugs, should be added when fixed # sd_regression("SD_Force" "subdyn;offshore") # sd_regression("SD_AnsysComp4_UniversalCableRigid" "subdyn;offshore") From 0aff0398ed6441b3b3aec6afd75560c6657ec972 Mon Sep 17 00:00:00 2001 From: Roger Bergua Date: Mon, 8 Jun 2026 19:22:21 -0600 Subject: [PATCH 9/9] Add SubDyn self-weight regression test --- reg_tests/r-test | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/reg_tests/r-test b/reg_tests/r-test index 25466634a..d87852239 160000 --- a/reg_tests/r-test +++ b/reg_tests/r-test @@ -1 +1 @@ -Subproject commit 25466634afabc3e697f3e7ca0897e1b84b3045b5 +Subproject commit d87852239e5f5fa79ff896c491b113fd64b0ede2