Hello,
I am working with the linearized state vector from the state-space system. After the full linearization at a operating point I get the state vector and the averaged Systemmatrices.
Can I interpret the values of the state vector as factors to the mode shapes described with a sixth order polynomial in the ElastoDyn file?
Thanks for your help.
Best regards!
Dear @Judith.Hess,
Yes, that is correct. The bending degrees of freedom of the blades and tower are factors that scale the user-specified mode shapes (shape functions). That is, the actual transverse blade and tower displacements are linear combinations of the degree of freedom multiplied by the shape functions.
Best regards,
1 Like
Dear @Jason.Jonkman,
thanks for your reply.
I have a further question concerning the stability of the system after the linearization.
I ran the runCampbell.m skript to generate a linearized system at the OP (WindSpeed_[m/s] , RotorSpeed_[rpm] , PitchAngle_[deg] , GeneratorTorque_[Nm]:15, 12.1, 10.45, 43094). I saved the averaged Systemmatrices from mbc_data in the function getCampbellData.m.
(I did so by adding the following lines to getCampbellData.m:
AvgA=mbc_data.AvgA;
AvgB=mbc_data.AvgB;
AvgC=mbc_data.AvgC;
AvgD=mbc_data.AvgD;
filename = [âSystemmatrices_OP_â, num2str(WindSpeed(iOP)),â.matâ];
save(filename, âAvgAâ, âAvgBâ, âAvgCâ, âAvgDâ))
When I have a look at the stability properties of that system loaded to matlab, it turns out that the first two modes (I think ordered by frequency, with a very low frequeny of 3.59e-03) from the modalanalysis damp(sys) are unstable with damping ratio = -1.00 & 1.00 (sys=ss(AvgA,AvgB,AvgC,AvgD)).
I also had a look at the Campbell_DataSummary.xlsx as an output from runCampbell.m. Comparing the damping ratios with the damping ratios as the output from damp(sys) it looks like damp() shows two more modes, the two unstable ones. I can´t explain this.
How can I deal with this problem or must there be a mistake in my initialization of the linear System?
Furthermore I have a question regarding the paper âA digital-twin solution for floating offshore wind turbines validated using a full-scale prototypeâ (WESD - A digital-twin solution for floating offshore wind turbines validated using a full-scale prototype
).
I want do a very similar approach like stated in the paper- (measuring top acceleration with additional SCADA Data to then use the dual kalman filter or augmented kalman filter to estimate the states). I donât understand how I can find the aerodynamic torque in the input vector [u]. Or asked in a diffrent way: if the aerodynamic torque is augmented in the state, how to extract the Jacobian out of the [B]- matrix, as I donât see it in the input vector [u].
Best regards,
Dear @Judith.Hess,
Iâm not following everything you are saying, but I would guess the modes showing up with +/- unity damping ratio are from a rigid-body mode, e.g., from generator rotation.
Iâm not sure I understand your question about our digital twin preprint.
Best regards,
Dear @Jason.Jonkman,
thanks for the quick reply. Iâll try to formulate my questions more clearly.
-
The system (with the matrices A,B,C,D of the state space model) that I get after linearization is unstable because it contains +/- unity damping ratio of a mode. If this comes from a rigid-body motion, I wonder how I have to adjust the linearization so that this does not happen. What would be an alternative to deal with the obtained system matrices to get a stable system? The goal is to use this system for a state estimation with a Kalman filter, for which it must be stable.
In the .fst file for the linearization the following two parameters are still zero, as i used the linearization example 5MW_Land_BD_Linear from git as a starting point, but want to use CalcSteady to linearize at a steady point.
0 Twr_Kdmp - Damping factor for the tower [used only if CalcSteady=True] (N/(m/s))
0 Bld_Kdmp - Damping factor for the blades [used only if CalcSteady=True] (N/(m/s))
I dontât know if what values are needed here and if the instability has something to do with the leeking damping factor here. -
In the first formulation of my problem, I have shown the first nine rows of the input vector of the .lin file, the remaining entries of the input vector refer to the blades. These only contain the four parameters of the operations point specified for the linearization process, that is WindSpeed_[m/s] , RotorSpeed_[rpm] , PitchAngle_[deg] , GeneratorTorque_[Nm].
I refer to the following equations in the above mentioned paper âA digital-twin solution for floating offshore wind turbines validated using a full-scale prototypeâ:
It is described that âThe Jacobians with respect to the loads are extracted from the B and D matricesâ.
The B and D matrices are related to the input vector. If the required loads are not included in this input vector, e.g. aerodynamic torque, these can also not be taken from the B and D matrix. How can I adapt my linearization, and thus the input vector, so that the Jacobians required for the augmented Kalman filter are contained in the system matrices B and D?
I hope I was able to formulate my questions more cleary.
Thank you so much for your help.
Best regards,
Judith Hess
Dear @Judith.Hess,
Regarding (1), does your model have the generator degree-of-freedom enabled (GenDOF = TRUE)? The generator DOF state will result in a rigid-body mode. You can simply eliminate this state (by deleting the associated rows and columns from the A and B matrices) to eliminate the rigid-body mode. (See related forum topics for more information.)
Setting inputs Twr_Kdmp and Bld_Kdmp can be used used to help speed up the calculation of the steady-state solution, but once the operating point has been found, these inputs will not directly effect the linearized state-space matrices extracted from OpenFAST.
Regarding (2), Iâll ask an NREL colleague to respond.
Best regards,
Hi @Judith.Hess,
-
A small note, make sure you output the linear models with enough precision, setting
OutFmtto âES17.8E2â or even more significant digits. The precision in the .lin file can greatly affect the eigenvalue analysis. You can try to track down how the eigenvalue analysis is done with the Matlab toolbox and haw itâs done in your own approach to see why you observe differences. -
For the Aerodynamic torque, I used the openfast output channel
RtFldMxh_[N-m]. To get the full âBâ and âDâ matrices, you need to setLinInputs=2'LinOutputs=2` in the OpenFAST input file. Youâll get more inputs and outputs from which you can extract the relevant elements of the B and D matrices. You need care to convert the units of your matrices if the output channels are not in SI units. I introduced some tools to ease this in Python but not for Matlab.
Iâm enclosing a list of relevant OpenFAST outputs I used at the end of this post.
' x ' : ' PtfmSurge_[m] ' ,
' y ' : ' PtfmSway_[m] ' ,
' z ' : ' {PtfmHeave_[m]} ' ,
' phi_x ' : ' {PtfmRoll_[deg]} * np.pi/180 ' , # SI [deg] -> [rad]
' phi_y ' : ' {PtfmPitch_[deg]} * np.pi/180 ', # SI [deg] -> [rad]
' phi_z ' : ' {PtfmYaw_[deg]} * np.pi/180 ' , # SI [deg] -> [rad]
' psi ' : ' {Azimuth_[deg]} * np.pi/180 ' , # SI [deg] -> [rad]
' q_FA1 ' : ' Q_TFA1_[m] ' ,
' q_SS1 ' : ' Q_TSS1_[m] ' ,
' dpsi ' : ' {RotSpeed_[rpm]} * 2*np.pi/60 ' , # SI [rpm] -> [rad/s]
' dq_FA1 ' : ' QD_TFA1_[m/s] ' ,
' dq_SS1 ' : ' QD_TSS1_[m/s] ' ,
' dx ' : ' QD_Sg_[m/s] ' ,
' dy ' : ' QD_Sw_[m/s] ' ,
' dz ' : ' QD_Hv_[m/s] ' ,
' dphi_x ' : ' QD_R_[rad/s] ' ,
' dphi_y ' : ' QD_P_[rad/s] ',
' dphi_z ' : ' QD_Y_[rad/s] ' ,
' ddpsi ' : 'QD2_GeAz_[rad/s^2]' ,
' ddq_FA1' : 'QD2_TFA1_[m/s^2] ' ,
' ddq_SS1' : 'QD2_TSS1_[m/s^2] ' ,
' ddx ' : 'QD2_Sg_[m/s^2] ' ,
' ddy ' : 'QD2_Sw_[m/s^2] ' ,
' ddz ' : 'QD2_Hv_[m/s^2] ' ,
' ddphi_x' : 'QD2_R_[rad/s^2] ' ,
' ddphi_y' : 'QD2_P_[rad/s^2] ',
' ddphi_z' : 'QD2_Y_[rad/s^2] ' ,
' Thrust ' : ' RtFldFxh_[N] ' ,
' Qaero ' : ' RtFldMxh_[N-m] ' ,
' Qgen ' : ' {GenTq_[kN-m]} *1000 '+'*{}'.format(nGear) , # [kNm] -> [Nm]
' WS ' : ' RtVAvgxh_[m/s] ' ,
' pitch ' : ' {BldPitch1_[deg]} * np.pi/180 ' , # SI [deg]->[rad]
' NcIMUAx ' : ' NcIMUTAxs_[m/s^2] ',
' NcIMUAy ' : ' NcIMUTAys_[m/s^2] ',
' NcIMUAz ' : ' NcIMUTAzs_[m/s^2] ',
' NcIMUVx ' : ' NcIMUTVxs_[m/s] ',
' NcIMUVy ' : ' NcIMUTVys_[m/s] ',
' NcIMUVz ' : ' NcIMUTVzs_[m/s] ',
I hope that helps,
Emmanuel
Dear @Emmanuel.Branlard and @Jason.Jonkman,
Thank you for your reply.
-
Changing OutFmt did not solve the problem, but turning off GenDOF as recommended by Jason Jonkman did. I am not sure if this is the right approach for my linearization purpose as I want to use the linearized matrices for Kalman filtering for a turbulent wind field, similar to your work in the paper where the generator is needed for the wind turbine. Still I want to linearize without accounting for aerodynamics, so WindSpeed=0 & RotorSpeed=0 as an operation point.
In this case, does it make sense to turn off GenDOF or is it better to delete the associated rows and columns of the matrices after the linearization process.
Sorry for the continuing questions, the functionality of the generator is for me, as a civil engineer, not easy to understand.
I have looked at the eiganalysis() function from the Matlab toolbox.
in eiganalysis(),line 55:
positiveImagEvals = find( imag(origEvals) > 0);
The two eigenvalues that are associated with the unstable damping +/- 1.00 are eliminated here, as their imaginary part is 0. This is the reason why those values are not to find in Campbell_DataSummary.xlsx. -
Thanks for the advise of changing the LinInputs and LinOutputs to get the full matrices. Unfortunately I canât relate the description of the rows of the following input vector in the ws00.0.1.lin, now with the size 1695, with the needed inputs for the kalman filtering approach. This is the reason why I canât identify the corresponding matrix elements, the Jacobian for the thrust and the aerodynamic torque.
I now know how to generate the output values of an openfast simulation, which I then need as inputs for the kalman filtering approach and that some of the needed inputs need to be estimated/calculated out of the given ones. Nevertheless am still not sure how to build up my B and D matrix for this.
I would very much appreciate your further help.
Best regards,
Judith
Dear @Judith.Hess,
Instead of disabling GenDOF, you can keep GenDOF enabled and eliminate the generator azimuth state from the linear model (but keeping the generator speed state). This has been discussed in several topics on this forum, e.g. see: FAST linearization V7.
Best regards,
Dear @Jason.Jonkman,
thanks for your reply.
I had a look at the previous topic of the forum.
Looking at my state vector, I canât find the generator azimuth state that I want to eliminate in the A and C matrix after the linearization.
Here is my state vector:
I am not sure if my settings for my intended linearization are right. I want to linearize without accounting for aerodynamics, so windspeed=0. Unfortunately. I couldnât find topics in the forum for this case.
My current settings are:
I am using the excel file provided in the matlab toolbox for the generation of the Campbell Diagramm defining the operating point, and the runCampell.m file to generate the linearization.
in .fst file:
1 CompElast
1 CompInflow
2 CompAero
1 CompServo
true Linearize
true CalcSteady
2 TrimCase
0.0001 TrimTol
0.001 TrimGain
0 Twr_Kdmp
0 Bld_Kdmp
in ServoDyn:
in AeroDyn:
in ElastoDyn:
Is this a valid approach for the linearization?
@Emmanuel.Branlard:
I canât find the Jacobian for the thrust and the aerodynamic torque in the B&D matric, as those two inputs are not in the input vector. What I understood you derive the input values seperately in your paper and you donât take them from the OpenFAST Simulation directly.
Can you explain me this in more detail?
Thanks for your help!
Best regards,
Dear @Judith.Hess,
Item 6 in your list is the generator azimuth state, so, you should eliminate column 6 from your A and C matrices and row 6 from your A and B matrices.
To eliminate aerodynamics, it is better to disable InflowWind and AeroDyn altogether by setting CompInflow = CompAero = 0.
Best regards,
Dear @Jason.Jonkman,
thanks for your help.
I eliminated the columns and rows from the A,B,C matrices.
Unfortunately, there is still the +1,0 damping in after the frequency analysis.
Do you have an explanation why this is the case?
Do you know, how to identify the Jacobian for the thrust and the aerodynamic torque in the B&D matrix, as those two inputs are not in the input vector, even in the full version. I understood the approach to derive the input values seperately in your paper and that you donât take them from the OpenFAST Simulation directly.
Can you explain me this in more detail?
Best regards,
Dear @Judith.Hess,
Can you confirm that the +/- damping associated with the rigid-body mode has now been has been eliminated? A damping > 1 indicates an overly damped mode, which often comes from aerodynamic damping in the flapwise direction.
@Emmanuel.Branlard would be better answering your question about the B and D matrices used in the digital twin paper.
Best regards,
Dear @Jason.Jonkman,
this is the result for the modal analysis that I get after eliminating the 6th state.
Only the -1.00 damping rate was eliminated here.
I changed CompInflow & CompAero to 0 in the .fst file.
Unfortunately, I cannot confirm, that the +1.00 damping is eliminated.
Best regards,
Dear @Judith.Hess,
I would guess the mode with the +1 damping corresponds to the generator speed state, which no longer has the corresponding displacement state. What is you specific concern?
Best regards,
Hi,
In my message above, I listed the mapping between the variables I used and the channels of OpenFAST (at least thatâs what I think I used).
For the aerodynamic thrust and torque, I mentioned RtFldFxh_[N] and RtFldMxh_[N-m] (you can replace Fld by Aero, this is an OpenFAST version issue). These are output channels of AeroDyn, so you need to list them in your outlist of AeroDyn, or use LinInput=2, LinOutput=2. Then you can scroll through the list in the .lin file to find the right index.
For the thrust, I think that I used the ElastoDyn loads at the nacelle (which the python-toolbox nicknames as NacFxN1_[N]). It is also a good option. Typically ElastoDyn loads include gravity and inertial contributions, whereas AeroDyn loads contain only the aerodynamics. So it depends what you are after.
The python-toolbox attributes a shorter names to the input and outputs: dfs = FASTLinearizationFile('file.lin').toDataFrame(); print(dfs['B']) ; print(dfs['B'].columns); You can then easily access coefficients as follows: dfs['B']['PtfmSurge_[m]', 'NacFxN1_[N]'] without having to worry about the index. It could be nice to implement something similar in Matlab. If you do please share it with us via a pull request to the matlab-toolbox.
I hope that helps,
Emmanuel
1 Like
Dear @Emmanuel.Branlard,
To clarify, my aim is not to extract the linearized system of a floating wind turbine, as previously mentioned. Rather, Iâm interested in a simpler approach similar to the methodology outlined in your earlier paper titled âA digital twin based on OpenFAST linearizations for real-time load and fatigue estimation of land-based turbines.â (doi:10.1088/1742-6596/1618/2/022030) I am sorry for the confusion.
My query pertains to whether the matrices required for the Kalman filter are directly available in this reduced version from OpenFAST. Specifically, Iâm uncertain if the settings in OpenFAST need to be adjusted for the linearization to accommodate a state vector containing only 5 states (How would this work?), or if the matrices subsequently need manual reduction for the Kalman filter.
I am referring to these lines in the paper:
My goal is it to get the linear state-space matrices A,B,C,D from OpenFAST corresponding to these states, inputs & ouputs, preferably also including the second mode as part of the states.
Furthermore, in the paper you compare the linearized system from YAMS and OpenFAST. You are saying that you canât get the monolitic system of equation in a straightforward way, but directly the matrices of the state-space form. For the comparison it is necessary that the dimensions of the matrices are the same, so I hope that it is possible to get them in the right dimension directly from OpenFAST
Given that thrust and aerodynamic torque are not included in the input vector, even in the full version with lininput=2, I have been contemplating whether itâs feasible to reconstruct the Jacobians for these values in the B matrix for the Kalman filter. My assumption was that since these values can be included in the output vector, the Jacobians might be composed from the C and D matrices of the output equation. However, due to the involvement of numerous components from the input vector that play a role here, Iâm uncertain if this approach is appropriate. You said, that you are using the loads at the nacelle for the thrust. Are those loads these values of the input vector? And therefore the corresponding values from the B and D matrix the Jacobians for the augmented approach?
Is there a similar approach for aerodynamic torque? Does the aerodynamic torque appear in the input vector resulting from the systems compression to the reduced amount of states?
Setting LinInputs=2 and LinOutputs=2 significantly increases the state space system, whereas I am particularly interested in the highly reduced version, as mentioned in the paper, with 4 states before the state vector is augmented.
I would greatly appreciate any further guidance regarding the implementation of such a reduced linearized state space system for a Kalman filter in the context of wind energy systems.
Thank you very much for your time and assistance.
Best regards,
Hi Judith,
I think I understand the approach you are trying to use. You want to get a reduced order model similar to what we did in these publications. To do so, I extracted a subset of the ABCD matrices from OpenFAST (few lines and columns). And since I used an augmented approach, sometimes I had to do sume surgery to take some values from the large B matrix of OpenFAST, and put it into the A matrix of the reduced order model. You might need to do similar surgery using values from the D matrix.
I think you are on the right track. The best is to scroll through all the inputs and outputs vector, to figure out which subset of terms you need to use. For inscance, the term dfs['B']['QD_FA1_[m/s]', 'NacFxN1_[N]'] , would give the jacobian of the tower top fore-aft degree of freedom velocity with respect to the nacelle tower top force. Which you might put in the A matrix of your reduced order system (if you use an augmented appraoch), or in the B matrix (if you chose to have the thrust as an input).
For the aerodynamic torque and generator torque, you can look for: RtAreoMxh_[N-m] and GenTq_[kN-m]. These are the short hand notations from the python toolbox.
I hope that helps.
Emmanuel