How to obtain the division units of the wind turbine structure?

[icnosvbovwb]

Dear @jjonkman
May I ask how the finite element model of a wind turbine is divided into elements? If it is matrix information, how can I obtain the overall matrix? Because I need to create an active control model, the linear quadratic regulator needs to obtain control parameters based on the finite element model of the structure. So I need to know where to obtain the division units of the wind turbine structure.

[jjonkman]

Dear @icnosvbovwb,

I'm not sure I fully understand your question, but only two of OpenFAST's many modules are based on a finite-element model; these modules include BeamDyn for blade structural dynamics and SubDyn for substructure structural dynamics. In these modules, the number of elements is specified by the user.

Given the modularization of OpenFAST, the equations of motion are not written in a single monolithic form. That said, you can extract the linear-equivalent form of the equations of motion in terms of linearized state-space matrices via the linearization functionality of OpenFAST.

Best regards,

[icnosvbovwb]

Dear @jjonkman
Yes, you are right，I want to extract the linear-equivalent form of the equations of motion in terms of linearized state-space matrices via the linearization functionality of OpenFAST. I went through some operational procedures because I am not familiar with linearization. Please help me if it is correct.
Q1: Is the following process a universal process for linearization? Is it right?
Q2：The first image is a state variable obtained through the linearization process, but there are only 29 of these variables. How can I add or subtract them as needed? For example, I only need the displacement, velocity, and acceleration in the tower, as well as the pitch and roll of the platform.
Q3：Where are these variables defined?

1）Open matlab-toolbox-main\Campbell\example folder_ Servo.csv file, inside
for pre-set work points. In order to generate a linearized model for the pitch section, the working points set are as follows:

2) Run the script runCampbell. m, and the program will determine the working condition of the working point. If the working point is in the MPPT stage, it will automatically adjust the electromagnetic torque based on the given pitch angle. If the working point is in the pitch stage, it will automatically adjust the pitch angle based on the given electromagnetic torque to ensure that the pre-set working point is a stable working point. The output of the working point, state matrix, and stable point after successful linearization appears in the form of. lin in the Matlab toolbox main_ ExampleData \ 5MW_ Land_ Lin_ Under the Trim folder.
3) Modify Main Flags in runCampbell. m as follows:

4) Copy fx_mbc3.m from the Matlab-toolbox-main\MBC\Source folder to the Matlab-toolbox-main\Campbell\example folder and rename it to MBC_mine.m. Because fx_mbc3 does not output control actions and output quantities for the working point, modifications need to be made based on the fx_mbc3 function. Change fx_getMats to fx_getMats_mine in line 41 of MBC_mine.m. Add the following code to line 253:

5) Copy fx_getMats.m from the Matlab-toolbox-main\MBC\Source folder to the Matlab toolbox main Campbell example folder and rename it to fx_getMats_mine.m. Add the following code snippet after line 143 to output the control input and system output for the work point:

6) Add the following code at the end of runCampbell. m, where windspeed_ Name is the wind speed name of the working point:

7) Run (runCampbell.m). In the workspace, MBC AvgA represents the state space matrix A that has undergone MBC transformation after azimuth averaging. The variables in the matData structure have not undergone MBC transformation, but contain information on the input and output of the working point as well as the state. Obtained MBC AvgA, MBC AvgB, MBC AvgC and MBC AvgD also needs to eliminate the rows and columns related to the rotor azimuth angle.
Open MBC DescStates, find the line where the ED Variable speed generator DOF is located (such as line 5). Then eliminate the row and column related to the fifth row in state x in the AB matrix. Meanwhile, the control variables for the 10th to 66th rows of leaf segments in the control input can be removed, and the relevant variables need to be removed from the CD matrix. The following functions were written for transformation:

8) Write a script to generate the state space matrix and system input/output structure for wind speeds of 14m/s

9) Copy StateSpace_14.mat to ROSCO-2.8.0\Matlab_Toolbox folder.
StateSpace_14.uop is the control input at the working point, and rows 4, 5, and 6 represent the pitch angle of 3 blades, respectively.
10) Load ('StateSpace_14. mat ') in runFAST. m
11) Modify the electromagnetic torque controller and pitch controller in Simulink based on the control variables in StateSpace_14.uop. The electromagnetic torque is set to a fixed 43094, and the pitch controller is set to a fixed 0.1472 (rad)

12) Modify the NRELOffshrBsline5MW_InflowWind.dat in the ROSCO-2.8.0\Test_Cases\NREL-5MW folder and set the HWindSpeed to 14.

[jjonkman]

Dear @icnosvbovwb,

Regarding (1), I'm not following everything in your process and don't know if you want me to comment specifically on something. Documentation on generating Campbell diagrams with the MATLAB scripts is provided here: https://github.com/OpenFAST/matlab-toolbox#campbell-diagram.

Regarding (2), the number of structural states that are included in the linearized models generated by OpenFAST are dictated by which degrees of freedom are enabled in the model. If you only want to see states associated with the tower-bending and platform pitch and roll degrees of freedom, disable the other degrees of freedom in ElastoDyn.

Regarding (3), these are the states (displacements and velocities) of the ElastoDyn degrees of freedom that you have enabled.

Best regards,

[icnosvbovwb]

Dear @jjonkman
Q1: I found through the linearization function that the output nx of onshore wind turbines is different from that of offshore wind turbines, with offshore wind turbines nx=138 and onshore wind turbines nx=32. This means that offshore wind turbines have an additional 94 outputs of hydrodynamic HD. What are these descriptions, similar to the first ExctnPtfmSg1? I don't know what these 94 outputs represent about HD?
Q2: I have processed the onshore wind turbine with equivalent linearization and compared it with the results of nonlinear model calculations, and the results (RootFxb1, RotSpeed) are relatively close. But now I need to compare the calculation results of linearization and nonlinearity for offshore wind turbines. I don't know what other factors need to be added besides wind speed processing? I guess it should be wave and current loads. For example, how to use waves, ocean currents, and wind as external loads to calculate a linearized model?

[jjonkman]

Dear @icnosvbovwb,

Regarding (1), the states labeled as HD ExctnPtfmSg1 etc, are hydrodynamic states associated with the linearized state-space form of wave-excitation from the potential-flow solution of HydroDyn. The following paper provides more information on this linearization capability of OpenFAST: https://www.nrel.gov/docs/fy19osti/71865.pdf.

Regarding (2), in addition to disturbance inputs from wind, the OpenFAST linearization functionality supports disturbance inputs associated with wave excitation for models with the potential-flow solution implemented in linearized state-space form. Again, the paper I linked above provides the theoretical background. Verification of this capability--comparing linear and nonlinear solutions--is presented in the following paper: https://www.nrel.gov/docs/fy20osti/74210.pdf.

Best regards,

[icnosvbovwb]

Dear @jjonkman
Q1: I want to enable linearization of the tuned-mass-damper (TMD) models and external user-specified controllers within ServoDyn.
Can you give me some advice? Do I need to set ServoDyn_StC file to enable TMD and then enable .fst file linearization function. Then run .fst file. It will get ABCD. Does these 4 matrixs include information about the TMD module?
Q2: Assuming my understanding and operation are correct, the information I generate includes six parameter information for TMD. I obtained four matrices ABCD by linearizing the TMD model of an onshore wind turbine, and then wanted to use these four matrices to actively control LQR, obtain the optimal control force, apply a linearized wind turbine, and form a closed-loop force. How can I transfer the calculated control force back to the linear model for calculation?

[jjonkman]

Dear @icnosvbovwb,

The linearization functionality of OpenFAST currently supports linearization with passive TMDs enabled within the Structural Control module, but the linearization functionality does not currently support any active TMD control options. So, the displacement and velocity of the TMDs will be included in the states (x) of the linear model; the motion of the body the TMD is connected to will show up in the inputs (u) of the linear model; and the mass, damping, and stiffness of the TMDs will appear within the matrices (ABCD) of the linear model; but you will not find external loads or related states showing up in the linear model. Including active control of TMDs within the OpenFAST linearization functionality would require further code development.

Best regards,

[icnosvbovwb]

Dear @jjonkman
Q1: I get the displacement and velocity of the TMDs included in the states (x) of the linear model. But I don't know how to develop these linearized displacements and velocities. Can you provide detailed suggestions？
Q2: Can you tell me if there is a old SIMULINK version of the active control interface in the versions before fast? It seems like it was modified by a doctor, but now the website can no longer be found. Is it fast v6 or fast V7?

[jjonkman]

Dear @icnosvbovwb,

I'm sorry, but I'm not sure I understand either of your questions. Please clarify what you are asking.

Best regards,

[icnosvbovwb]

Dear @jjonkman
Sorry, I didn't express my previous meaning clearly, what I want to ask now is that:

The left side of the first image shows the added AMD system, and the right side of the first image shows the particle spring system with AMD. Assuming m1 is the mass of the wind turbine and m2 is the mass of the TMD, if I enable StC_ DOF_ MODE=1, and Linearize=True, the resulting ABCD matrix contains TMD information. But now I need to design an active controller LQR and calculate an active control force. That is the control method of AMD. My question is: How to expand the linearized ABCD? How do B and D matrices expand to include the F(t) information circled in the first graph? Because without this F(t), the entire system only has TMD, so I can linearize it directly to obtain ABCD containing TMD information. Am I not sure if this understanding is correct? Because I linearized the wind turbine with TMD using Openfast to obtain ABCD, I want to obtain ABCD with active control information, that is, how to add the information in B and D.

[jjonkman]

Dear @icnosvbovwb,

It looks like you've derived the B matrix you want for force F for the simple 1 dimensional, 2 degree of freedom system. If you've linearized OpenFAST with only platform surge and TMD-x degrees of freedom, you would obtain A through an OpenFAST linearization. So, I think your challenge is to follow a similar process to derive B if A is expanded to include more OpenFAST degrees of freedom. Following the derivation of the TMD degrees of freedom as implemented in the Structural Control module would likely aid in this effort: https://openfast.readthedocs.io/en/main/source/user/servodyn-stc/StC_Theory.html#theoretical-background.

Best regards,

[icnosvbovwb]

Dear @jjonkman
I am currently expanding the linearized matrix of the onshore wind turbine. I have set up a TMD at the top of the onshore wind turbine tower ( NumTStC=1, StC_X_DOF=true, StC_Y_DOF=FALSE, StC_Z_DOF=FALSE, StC_P_Z=77.6, StC_P_X=StC_P_Y=0) and applied an active control force F(t), as shown in the first image of assuming the simple 1 dimensional, 2 degrees of freedom system (assuming m1 is the mass of the wind turbine and m2 is the mass of the TMD). I have made modifications to the B and D matrices.
Q1: By setting (Linearize=True, YawDOF=FALSE), it will get A(36x36) B(36X66) C(44X36) D(44X66). I only modified the B and D matrices, is it correct? as shown in the second image.
If it is correct, the modifications I made are: adding the 67th column to matrix B, modifying the fourth (1/mTMD) ('srvD Tower StC 1 local display state dX/dt m/s'), and modifying the 22nd (1/m rotor+tower+nacelle) (ED First time derivative of 1st tower for after bending mode DOF (internal DOF index=DOF-TFA1), m/s) in column 67 of matrix B. Also adding the 67th column to matrix D, modifying the 9th (1/Kx), (ED TTDspFA, (m)), and the 43rd (-1) XTowerHt) has been modified in column 67.
Q2: Is the 22nd in the 67th column of the B matrix the rotor mass+tower mass+nacelle mass? Is the 9th in the 67th column of the D matrix the displacement 1/kx? I don't know how to calculate the displacement caused by the active force at this position at the top of the tower. Is the 43th in the 67th column of the D matrix the -1XTowerHs?

[jjonkman]

Dear @icnosvbovwb,

I don't believe these changes are fully correct. Just a few comments regarding terms I don't believe you are correctly accounting for:

• With all of the turbine DOFs enabled, the effective mass will differ from M(rotor+tower+nacelle) because the system is not act as a rigid body- With all of the turbine DOFs enabled, the effective mass will differ from M(rotor+tower+nacelle) because the system is not act as a rigid body.
• The TMD DOFs in StC are relative to the motion of the body they are attached to (rather than expressed globally as in your analytical example).
• When the tower deflects due to bending, it also rotates, so, the x-based TMD will also move vertically by as small amount when the tower is deflected (so, the problem isn't 1D as in your analytical model)
• A force applied the tower will only deflect the tower based the inverse of the tower stiffness in a quasi-static condition; when dynamics are accounted for, forces induce accelerations (not deflections directly), which induce deflections only when integrated in time.

Best regards,

[icnosvbovwb]

Dear @jjonkman
Can you provide me with specific guidance on how to proceed? Because I don't quite understand the principle of TMD, my understanding of TMD is to use the Newmark algorithm to calculate TMD (bying M/K/C), so I may only think it is a unidirectional vibration. However, the TMD calculation of onshore wind turbines may include calculations in three directions, so could you please instruct me how to operate it?

[jjonkman]

Dear @icnosvbovwb,

I'm sorry, but I'm not really sure I understand your question. I doubt I can provide more guidance than I already provided given that NREL is not currently funded to work on TMDs and add the linearization functionality you desire.

Best regards,

[icnosvbovwb]

Dear @jjonkman



Q1: 5MW Baseline wind turbine Tower-top Mass is 349389.844kg. Is Tower-top Mass equal to Rotor Mass(109389.852kg) + nacelle mass(240000kg) ?
Q2: I don't know where to find the R matrix circled in the picture?
Q3: Are both YawBrFxn and YawBrFxp forces acting on the top of the tower? What is the difference between them?

[jjonkman]

Dear @icnosvbovwb,

Here are my answers to your questions:

1. Yes, that is correct.
2. R_N/G is the 3x3 direction cosine matrix to transform a vector in global coordinates (G) into the local coordinates of the body where the TMD is placed (N = nacelle for a nacelle-based TMD). Given that the dynamics of the body where the TMD is placed is solved outside of the StC submodel of ServoDyn, this matrix is an input to ServoDyn from the calling program (derived from the output of ElastoDyn, BeamDyn, or SubDyn depending on what's enabled and where the TMD is placed).
3. These ElastoDyn outputs are both reaction forces at the yaw bearing, just expressed in different coordinate systems. YawBrFxn is the shear force along the x axis of the nacelle coordinate system and YawBrFxp is the shear force along the x axis of the tower-top coordinate system (the nacelle- and tower-top coordinate systems are aligned when the nacelle-yaw angle is zero).

Best regards,

[icnosvbovwb]

Dear @jjonkman

How can I specifically obtain this R matrix? I want to extract it as a constant and place it in the Simulink module, and then multiply it by the three forces and moments of XYZ to output. Because I only added an active force to the control force in the X direction of TMD and modified the B and D matrices.

[jjonkman]

Dear @icnosvbovwb,

I'm not sure I really understand your question. But the direction cosine matrix (R) that I described in my prior post comes from a sequence of rotations from the platform to the nacelle (or wherever you are placing the TMD), including platform rotation (roll, pitch yaw), tower-top deflection (fore-aft, side-side), and nacelle yaw. Of course, the transformation matrix is stored within the OpenFAST source code and passed directly from the glue code through SerovDyn to the StC submodel.

Best regards,

[icnosvbovwb]

Dear @jjonkman
The image shows a linearized and nonlinearized calculation model for a 5MW onshore wind turbine. I have successfully constructed a linearized model through MATLAB/Simulink, and the RootFxb1 operation results are shown in the second figure (14m/s to 22m/s step wind), and the TTDspFA operation results are shown in the third figure (14m/s to 15m/s step wind).
My goal is to derive the active control force from the linearized onshore wind turbine model based on the above two degrees of freedom models, and add this active force to the nonlinear onshore wind turbine model.
Based on your previous suggestions, I have successfully modified the source code and set up StC_ DOF_ MODE=4, which means that the nonlinear wind turbine has two NRELOffshrBsline5MW_ ServoDyn_ StC. dat. I set StC_ DOF_ MODE=1 in ServoDyn_ StC1.dat, set StC_ DOF_ MODE=4 in ServoDyn_ StC2.dat, which applies the active force and TMD control force together to the top position of the tower.
The current minor issue is: I add six outlists (TStC1Fxl, TStC1Fyl, TStC1Fzl, TStC1Mxl, TStC1Myl, TStC1Mzl) in NRELOffshrBsline5MW_Onshore_ServoDyn.dat, and linearized them to generate four matrices ABCD. I modified the B and D matrices (Marix B Last line and Marix D Last Line). I need to multiply the TMD (XYZ forces and XYZ moments) in the D matrix by cosine rotation conversion, and then transfer the force in the last column to the interface of the nonlinear wind turbine. I'm not sure how to derive this cosine matrix?
Can you provide me with detailed assistance on how to extract this R matrix, as I need this cosine matrix to convert the force and moment of the active control force.

[jjonkman]

Dear @icnosvbovwb,

I'm not following everything you are saying, but just one clarification: When you use StC_DOF_MODE = 4 for an StC located at the tower-top where another StC is placed, the StC load from StC_DOF_MODE` = 4 will be applied to the tower top, but not to the other StC. So, this does not exactly match the square red box in your first image.

It sounds like your question is about how to find the direction cosine matrix of the tower-top within the StC module during a nonlinear time-domain simulation. This direction cosine matrix can be found via input array u%Mesh(i_pt)%Orientation(:,:,1), which you can find already used within SUBROUTINE StC_CalcOutput() of source file StructCntrl.f90 when PrescribedForcesCoord = 2 (local coordinates).

Best regards,