# Output the angular velocities of a blade segment

Hi,

I check FASTv7 manual and I couldn’t find any description whether is it possible to get the angular velocity of a blade segment.
Is this included in the FAST code?

Regards,
Robert

Dear Robert,

While the angular velocities along the blade are calculated internally, they are not currently output from FAST v7 or FAST v8. It would take a minor modification to the FAST source code to include these outputs.

Best regards,

Dear Jason,

As I understood from the code, this value would be stored in the

variable, which is currently disabled in the code. If I simply uncomment the following lines, I would get the desired values?

[code]! AngVelHM(K,TipNode , = AngVelEH + QDT(DOF_BF(K,1))*PAngVelEM(K,TipNode,DOF_BF(K,1),0, & ! Currently
! + QDT(DOF_BF(K,2))*PAngVelEM(K,TipNode,DOF_BF(K,2),0, & ! unused
! + QDT(DOF_BE(K,1))*PAngVelEM(K,TipNode,DOF_BE(K,1),0, ! calculations

! NOTE: These are currently unused by the code, therefore, they need not
! be calculated. Thus, they are currently commented out. If it
! turns out that they are ever needed (i.e., if inertias of the
! computations:
! PAngVelEM(K,TipNode, :,1, = PAngVelEH(:,1, ! CALL CrossProd( PAngVelEM(K,TipNode,DOF_BF(K,1),1,:), AngVelEH, &
! PAngVelEM(K,TipNode,DOF_BF(K,1),0, )
! CALL CrossProd( PAngVelEM(K,TipNode,DOF_BF(K,2),1,:), AngVelEH, &
! PAngVelEM(K,TipNode,DOF_BF(K,2),0, )
! CALL CrossProd( PAngVelEM(K,TipNode,DOF_BE(K,1),1,:), AngVelEH, &
! PAngVelEM(K,TipNode,DOF_BE(K,1),0, )
[/code]

Regards,
Robert

Dear Robert,

Actually, what you want is PAngVelEM, where PAngVelEM(K,J,I,0, = the partial angular velocity of DOF I for element/node J of blade K. Summing these terms multiplied by their respective first derivative of the DOFs (QDT(I)) over all DOFs (I) will yield the total angular velocity of element/node J of blade K.

I hope that helps.

Best regards,

Hi,

Thank you for the help. This answered my question.

Regards,
Robert

Hello,
I’m currently also interested in the output of angular velocities of blade segments. These velocities have to be in the blade coordinate system and should be equal to the measurements a gyroscope mounted on the blade would output.

I’m using a workaround by extraction of the blade segments coordinate system from the vtp files and calculation of angular velocities from the change of orientation. This works somehow, but is computationally heavy and uncomfortable as I have to parse hundreds of files.

As I am not familiar with Fortran at all, I somehow went into troubles trying to implement the above mentioned code. So I would kindly ask if someone can tell me exactly what to change in the source code to enable output of angular velocities. And is the desired element/node J configurable via config file or hard coded?

Best regards
Eike

Dear Eike,

I’m not going to change the source code for you, but all you need to do is for node J of blade K you must loop through all of the structural DOFs and sum the partial angular velocity times the velocity of that DOF, i.e.:

```AngVelEM = 0.0 DO I = 1,p%DOFs%NActvDOF AngVelEM(K,J,:) = AngVelEM(K,J,:) + m%RtHS%PAngVelEM( K, J, p%DOFs%SrtPS(I), 0, : )*x%QDT( p%DOFs%SrtPS(I) ) ENDDO```
This code will give you the absolute angular velocity of node J of blade K in the inertial-frame coordinate system, i.e. AngVelEM(K,J,:). If you want this velocity projected into the local blade coordinate system, you can dot-product this vector with the local coordinate system direction, i.e.:

```AngVelEM_Lxb(K,J) = DOT_PRODUCT( AngVelEM(K,J,:), m%CoordSys%n1(K,J,:) ) AngVelEM_Lyb(K,J) = DOT_PRODUCT( AngVelEM(K,J,:), m%CoordSys%n2(K,J,:) ) AngVelEM_Lzb(K,J) = DOT_PRODUCT( AngVelEM(K,J,:), m%CoordSys%n3(K,J,:) )```
I hope that helps.

Best regards,

Dear Jason,
I try to adapt your code but went into some compiling errors. So, I reduced the code to just write a constant because I first want to get the new variable in the output files. Everything else will be the next step.

What I did so far:
ElastoDyn_Registry.txt

``typedef	^	ED_RtHndSide	ReKi	AngVelEM	{:}{:}{:}``

ElastoDyn.f90 @ SUBROUTINE CalculateAngularPosVelPAcc( p, x, CoordSys, RtHSdat )

``````DO K = 1,p%NumBl ! Loop through all blades
DO J = 1,p%BldNodes ! Loop through the blade nodes / elements
...
RtHSdat%AngVelEM(K,J,:) = 5.0 ! some constant
...
END DO !J = 1,p%BldNodes ! Loop through the blade nodes / elements
END DO !K = 1,p%NumBl``````

I’m now completely lost what to enter in ElastoDyn_IO.f90 to get the parameterizable output.

PS: Somehow nrel seems to be down. I cannot download the programmers instructions, which may also contain the necessary information. So please excuse my poor fortran skills…

Dear Eike,

What compiling errors are you getting? Did you ALLOCATE array AngVelEM?

I’m sorry, but I’m not sure what you mean by “parameterizable output”. I’m also not sure what mean by “NREL seems to be down”.

Best regards,

Dear Jason,

I really appreciate the effort you put in supporting people like me who try to extend/modify the source code with nearly no fortran skills
Please excuse the late answer, I was very busy on other projects.

How and where (ElastoDyn_Types.f90?) do I have to allocate AngVelEM?
I struggle with the Fortran syntax for arrays and all I did was some senseless try and error. So maybe it would be the best way to branch the OpenFast git project so you can see what I did?

Best regards,
Eike

Dear Eike,

I would put the code to calculate AngVelEM within SUBROUTINE ED_CalcOutput(). You shouldn’t need to update the ElastoDyn_Registry.txt file at all. You can then declare AngVelEM within SUBROUTINE ED_CalcOutput() as a fixed-size array, e.g.:

``REAL(ReKi) :: AngVelEM(p%NumBl,0:p%TipNode,NDims)``

I hope that helps.

Best regards

Dear Jason,

openFast now compiles without errors. I just get one warning

``Warning: ‘__builtin_memset’: specified size 18446744073709551587 exceeds maximum object size 9223372036854775807 [-Wstringop-overflow=]``

My code is

``````DO I = 1,p%DOFs%NActvDOF ! Calculate Angular Velocity
AngVelEM(K,J,:)   = AngVelEM(K,J,:)  + m%RtHS%PAngVelEM( K, J, p%DOFs%SrtPS(I), 0, : )*x%QDT( p%DOFs%SrtPS(I) )
ENDDO
AngVelEM_Lxb(K,J))= DOT_PRODUCT( AngVelEM(K,J,:), m%CoordSys%n1(K,J,:) )
AngVelEM_Lyb(K,J) = DOT_PRODUCT( AngVelEM(K,J,:), m%CoordSys%n2(K,J,:) )
AngVelEM_Lzb(K,J) = DOT_PRODUCT( AngVelEM(K,J,:), m%CoordSys%n3(K,J,:) )``````

and initialization

``````REAL(ReKi)                   :: AngVelEM(p%NumBl,0:p%TipNode,NDims)
REAL(ReKi)                   :: AngVelEM_Lxb(p%NumBl,0:p%TipNode)
REAL(ReKi)                   :: AngVelEM_Lyb(p%NumBl,0:p%TipNode)
REAL(ReKi)                   :: AngVelEM_Lzb(p%NumBl,0:p%TipNode)``````

So if I understand correctly, I need to add all variables I want to output in m%AllOuts.
If I do for example

``m%AllOuts(AngVelEM_Lxb(K,J)) = DOT_PRODUCT( AngVelEM(K,J,:), m%CoordSys%n1(K,J,:) )``

The compiler gives the following warning

``` m%AllOuts(AngVelEM_Lxb(K,J)) = DOT_PRODUCT( AngVelEM(K,J,:), m%CoordSys%n1(K,J,:) ) 1 Warning: Legacy Extension: REAL array index at (1)```
Can you please tell me how to fix it?

Best regards
Eike

Dear Eike,

The problem I see is that AngVelEM_Lxb(K,J) is declared as a REAL, but you are using it as an index in the m%AllOuts() array. Indexes must be INTEGERs. For now, I would not worry about including your new outputs in the FAST output file (by changing m%AllOuts()); instead, I would simply write your new outputs e.g. to the screen using PRINT statements.

Best regards,

Dear Jason,

thanks again. I’ll try out to print them to stdout. Nevertheless, a really like the binary output format and I parse it in matlab for postprocessing. So if AngVelEM* can somehow be included, it would be great benefit for me and my application Best regards