Verification of 10MW On Shore Linearization with OpenFAST v.2.4.0


I’m currently working on the linearization of 10MW OnShore RWT model. The parameters fixed for the linearization are rotor speed equal to 9.6 rpm, blade pitch equal to 9.292° for the three blades and steady wind with speed equal to 14 m/s, while the rotor azimuth angles selected for the linearizations are 0°:22.5°:337.5°(hence 16 angles were considered). The inputs considered are Blade1(B1) individual pitch, B2 individual pitch, B3 individual pitch and collective pitch, while the outputs considered are FlapWise bending moment at B1, B2 and B3 roots.

The model is defined using InflowWind, AeroDyn15, ElastoDyn and BeamDyn, in particular the WT was modelled without BeamDyn first; the purpose was to compare the results obtained by linearization performed on a flexible WT(FlapDOF1 EdgeDOF TwFADOF1 TwSSDOF1 activated) featured by isotropic blades(hence NO Bend-Twist Coupling) modelled by ElastoDyn(CompElast=1) first and then modelled by ElastoDyn+BeamDyn for blades(CompElast=2).

The desired result from the comparison was the coherence between the TFs obtained from linearizations performed on models in which different tools were implied to model same structural characteristics, hence input-output transfer functions obtained by a model defined through ElastoDyn must be quite similar to the ones obtained by a model defined through ElastoDyn+BeamDyn.

The desired results was not reached, in particular problems related to numerical instability(I think) are highlighted using BeamDyn; it turns out that using different OrderElement(I used 4 with time step=0.001s,5 and 6 with time step=0.0005s) number the resulting TFs are different. I attached a folder containing images related to TFs obtained using ElastoDyn only and ElastoDyn+BeamDyn with different OrderElement considering for all the cases rotor azimuth angle equal to 0°(the TFs are considered in a frequency range of [0.12pi rad/s : 2*pi rad/s] since that is my range of interest). I attached a folder which includes the models I used to perform linearization.

I’m asking help to understand the cause of TFs’ variation due to OrderElement changes and I would like to define the reason of the differences between the TFs computed using ElastoDyn only and using ElastoDyn+BeamDyn which are(should be) defining the same structural model.

Thanks in advance.

Simone Covre
Models.7z (19 KB)
Images.7z (1.06 MB)

Dear Simone,

I don’t look at transfer functions often, but that sounds like something that would be difficult to debug problems using. I would first confirm some simpler things, like:

  • Do the time series match well between ElastoDyn and BeamDyn?
  • Do simpler linearization post-processing, e.g. performing an Eigenanalysis on the state matrix “A” result in similar natural frequencies and damping ratios between ElastoDyn and BeamDyn?

When performing a linearization analysis, please keep a couple basic steps in mind, such as:

  • Ensure the solution is in steady state before linearizing (plot the time series to verify if need be).
  • Use a high precision for the linearization output, e.g. OutFmt = “ES17.9E3”.

Best regards,

Dear Jason,

About your useful advices:

  1. I think the time series match well between ElastoDyn and BeamDyn, the small gap is propably due to the fact that torsional DOF is added in BeamDyn model; I attach a folder which includes images related to the matching, in particular both matching and steady state could be verified looking at these.

  2. OutFmt = “ES17.9E3” was very useful since solved almost totally the problem related to the differences between resulting TFs obtained with different OrderElement. In the folder previusly mentioned there are images related to Bode plot where the improvement is very clear comparing with previous post’s results.
    I have one more doubt about this issue: the doubt is related to the azimuth average configuration, I think that the obtained azimuth-average is consistent only for the ElastoDyn model, since for what concern BeamDyn model the average has a strange behaviour that I couldn’t completely explain.

  3. I used the matlab function eiganalysis.m provided by the latest release of the Matlab Toolbox in order to perform the suggested post processing elaboration; I’m a little bit confused about the obtained results, in particular the number of natural frequencies are obviusly much higher with BeamDyn(89 for BeamDyn, 8 for ElastoDyn) and the damping ratios evaluated in the frequency range of ElastoDyn natural frequencies are not similar between Elastodyn model and BeamDyn model.

-The TFs between ElastoDyn model and BeamDyn model are still different, but the similarity between the two is improved.
-The TFs derived from azimuth average of BeamDyn model shows strange behaviour according to my opinion.
-The difference between TFs obtained by BeamDyn model with different OrderElement is nearly completely vanished, but the reason related to
diffent magnitude in corrispondence of spikes observed from BeamDyn model featured by different OrderElement is not clear.

Update images for TFs are included in the attached folder.

Thanks for your important support.

Simone Covre
ImagesReply1.7z (2.45 MB)

Dear Simone,

OK, it definitely looks like you’ve made progress since your prior post!

Is the green curve the TF calculated from the azimuth averaged matrices? Have you verified that the matrices are being averaged correctly? Are you applying MBC3 before azimuth averaging?

Regarding the damping, have you set the damping in BeamDyn so as to match ElastoDyn at least for the first bending modes? With stiffness-proportional damping used in BeamDyn, the damping will increase with natural frequency. So, it is likely only possible to match the structural damping well for the first mode, if the damping between ElastoDyn and BeamDyn is set consistently.

Best regards,

Dear Jason,

The green curve is calculated from the azimuth averaged matrices, in particular to compute those matrices I used the matlab function fx_getMats.m from Matlab Toolbox; turns out that this function works well with ElastoDyn model, while BeamDyn model requires fx_mbc3.m as you suggested.
I tried fx_mbc3.m also with ElastoDyn model but the bode diagram obtained from the azimuth averaged matrices after Multi-Blade Coordinate transformation wasn’t good, hence I would ask you: Multi-Blade Coordinate transformation is required only with BeamDyn?

I modified the damping section in BeamDyn Blade File in order to match ElastoDyn 1st Bending modes, in particular I followed the explanation you gave in another post related to BeamDyn Damping; I had infos only for flapwise and edgewise damping, hence I compute only mu44 and mu55, I set the other damping coefficient equal to 0.01s since mu44=0.01s and mu55=0.016s.
After this modification there is much closer matching between damping ratios for first modes, I attached a .xlsx file which includes the natural frequencies and related damping ratios for both ElastoDyn and BeamDyn cases.

Despite all the improvements, the TFs given by BeamDyn model is different from TFs given by ElastoDyn model.
In my opinion the most relevant difference between the two is given by the peaks, in fact the shape is quite similar between the two except for the behaviour in correspondence of peaks. Do you think it’s a problem related to damping? If you think so, do you have any idea about how it could be fixed?

Thanks for your help.

Simone Covre
ExcelReply2.xlsx (13.4 KB)

Dear Simone,

In ElastoDyn, the blade states are in a rotating frame of reference and the tower states are a fixed frame of reference. You’ve set RotStates = TRUE in BeamDyn, so, the blade states in BeamDyn are also in a rotating frame of reference.

Thus, for the Eigenanalysis, you should apply MBC3 and azimuth average both the ElastoDyn and ElastoDyn+BeamDyn results before performing the Eigenanalysis. Are the natural frequencies and damping you are showing based on this method?

I’m less clear how you should post-process the transfer functions between blade pitch and blade flapwise moment. Given the mixture of rotating and nonrotating states, I generally recommend applying MBC3 and azimuth averaging before subsequent post-processing. Of course, this will convert the individual blade states/inputs/outputs into collective, asymmetric sine, and asymmetric cosine terms. If you eliminate the tower modes–keeping everything in the rotating frame of reference–perhaps it makes sense not to apply MBC3.

How do some of the key elements of the state, input, and output matrices compare before and after the azimuth averaging? Making plots e.g. of A(i,j) as a function of azimuth angle and the azimuth-averaged A(i,j) may help clarify what the azimuth averaging is doing. Again, I would generally recommend doing this after applying MBC3.

Best regards,

Dear Jason,

I updated the results of the eigenanalysis; I considered 2 models entirely flexible (flexible tower+flexible rotor and flexible tower+flexible rotor modelled by BeamDyn) and 2 models featured by rigid tower in order to obtain continuos states belonging to the rotating frame. According to my opinion the results are a little bit confusing since it seems like I’m not setting BeamDyn Blade damping consistently; I should probably start from the natural frequencies obtained from this linearization to define first bending mode frequencies, in order to obtain closer matching. Do you think that the maching is satisfying?

I used tower rigid models for the other part of the processing since I would like to have all the states in the rotational reference frame.
I tested the azimuth averaging on key elements as you suggested and the results is perfect, hence the azimuth averaging capabilities are on point with BeamDyn model too, so I can’t figure out the reason why the TFs of the system related to the averaged matrices doesn’t match the azimuth TFs distribution just in case the rotor is modelled by BeamDyn.

I considered the TFs of the model with and without MBC (since I’m using a model featured by rigid tower the analysis makes sense without MBC too as we noticed from previous post) and the main problem is still related to the model associated to azimuth averaged matrices in which the flexible rotor is model by BeamDyn. I left two images to re-show this fact.

Recap: I think that currently the two main problems are related to BeamDyn Blade damping and to the unclear reason of the mismatching between TFs related to azimuth averaged model featured by flexible rotor modelled by BeamDyn and TFs related to azimuth distribution models with flexible rotor modelled by BeamDyn.

Thanks for your assistance and sorry for the huge time required for my reply.

Simone Covre
ImagesReply3.7z (785 KB)
ExcelReply3.xlsx (24.8 KB)

Dear Simone,

Regarding the blade damping, I’m not sure what to make of the differences you are showing. To better isolate the structural damping, I would first recommend comparing the eigensolution without aerodynamic loading/damping. Does the structural damping alone compare well between the ED only and ED+BD models, at least for the first flapwise and edgewise modes? For additional insight you could generate the Campbell diagram, plotting the natural frequencies and damping ratios as a function of rotational speed without aerodynamic loads. I would generally expect that you’d get good agreement between the ED only and ED+BD models for lowest bending eigenmodes of the blades and tower.

Regarding the transfer functions and azimuth averaging, I’m not sure I understand what you plotting in the MBC_AzimuthAverageWorking figures.

Best regards,

Dear Jason,

I am sorry for the lack of clarity of previous data. I followed your advice and I perfomed linearization procedure without aerodynamic damping/loads (CompAero=0) in order to highlight the matching between ED and ED+BD models between 1st flapwise and 1st edgewise bending mode frequencies and damping ratio.

I want to post the results without going further with an additional insight that requires Campbell diagram, because I think there is a problem related to the 1st flapwise bending mode damping ratio that must be solved before doing any other step.
From the results I consider that matching between the two models is achieved for 1st edgewise and flapwise frequencies and for 1st edgewise mode damping ratio, while for what concern 1st flapwise mode damping ratio a strange behaviour is highlighted.

1st flap mode damping ratio seems to depend on rotor azimuth angle and I can’t understand why. BeamDyn blade damping is tuned starting from the first mode natural frequencies obtained throught free decay analysis of the same model implied for the linearization(obviusly only ED model was considered, since ED+BD model doesn’t support out-of-plane and in-plane tip deflection).

You gave me this question with related suggestion two posts ago and MBC_AzimuthAverageWorking figures were supposed to show some key elements values of matrices A,B and C for different rotor azimuth angles and on the same plot a line y=AzimuthAverageKeyElementValue is plotted. In particular in each figure related to ED model I used different markers for different blades(B1,B2,B3) and the input and output for the matrix element considered are explained in the legend. The results I found are clearer(highlghted) in ED model case.
If the images are still not clear I’ll re-post with a clearer explanation in a post entirely dedicated to azimuth averaging.

Thanks a lot for your support.

Simone Covre
NoAeroDyn_LinearizationData.7z (2.04 MB)
ModelsReply4.7z (20.1 KB)
ExcelReply4.xlsx (16.3 KB)

Dear Simone,

Regarding the azimuth variation of the damping in ED+BD that you don’t see in ED or in the first edgewise mode of BD, I’m curious if this is related to some impact of gravity or shaft tilt and the geometric nonlinearities BD. Do you get more consistency if you perform the ED+BD analysis with ShftTilt = 0 deg and Gravity = 0 m/s^2 in ED? Since you are referring only to flapwise and edgewise modes, and not collective flap, asymmetric flap, etc., presumably this eigenanalysis is based on post-processing that does not include application of the Coleman transformation (MBC3); is that correct?

Regarding the azimuth averaging; thanks for clarifying. I think I understand what is plotted, at least for the ED azimuth variation results. What I don’t see in these plots are the actual azimuthally averaged values, which is what you are trying to verify and presumably would be the mean value of the azimuthally varying result and independent of azimuth angle. And I have a harder time interpreting the BD results because the figures are separated by blade, and some of the legends are not showing up correctly.

Best regards,

Dear Jason,

I performed the simulations that you suggested and the results are still not consistent, in particular the damping ratio related to 1st flapwise collective mode depends on rotor azimuth angle. You won’t find any attachment because the results are the same presented in previous file .xlsx.

I’m sorry about previous post, I forget to mention that the results are obtained after MBC transformation, hence should have specified that the results were related to 1st edgewise collective mode and 1st flapwise collective mode.

I’m sure that MBC transformation is not the source of inconsistency, since I performed the eiganalysis in rotating frame too and ED+BD model still shows inconsistency for what concern 1st flapwise mode(in this case we don’t talk about collective or asymmetric since we didn’t perform MBC transformation). We can also keep out from the possible inconsistency causes the gravity effect related to the rotation after last simulations.

Currently I have two main questions, first one regarding 1st mode flapwise collective damping ratio inconsistency, second one regarding the azimuth averaging, since considering the state space model derived from azimuth averaged matrices, after the MBC transformation and the application of eiganalysis we obtain consistent frequency of 1st mode flapwise and edgewise collective, consistent damping ratio related to 1st mode edgewise collective and unconsistent damping ratio related to 1st flapwise collective mode. In particular the previously mentioned inconsistent damping ratio derived from azimuth averaged matrices is way higher than the damping ratios obtain considering rotor azimuth angle distribution.
All the issue mentioned in previous lines are reported and highlighted in the file .xlsx attached in previous post.

Thanks for your support.

Simone Covre

Dear Simone,

I’m not sure. Here are few other items I’m curious about to understand this issue better:

  • Are the results converged with OutFmt = “ES17.9E3”, or have you tried even more precision, e.g. “ES20.12E3”?
  • Does rhoinf have any effect on the damping levels you are seeing (rhoinf = 1 would result in no numerical damping)?
  • What level of damping do you see (if any) when you eliminate all structural damping (mu = 0)?
  • Do you see the same azimuth-dependence on damping when linearizing with zero rotational speed? (this would require separate linearization simulations for each azimuth angle)

Best regards,

Dear Jason,

I would like to share with you the results currently obtained with the initial 3 points of the list of the recommended simulations I should perform; I think it’s a clever way to proceed because we could have the chance to set a more correct focus in last point simulations(null rotational speed simulations).

Increase of output format precision didn’t lead to better results, I think that convergence was reached with format “ES17.9E3” too.

The same simulation was carried out with and without numerical damping and the result is the same, hence numerical damping gives a negligible contribution to the overall damping computed by eiganalysis.

I dug mainly into this issue because I saw from the results that, including aerodynamic loads/damping(AD),numerical damping(BD) and keeping out structural damping(BD) , damping ratio’s dependency on azimuth angle disappeared. In particular, I performed 3 significative simulations, first one using output format “ES20.12E3” and including aerodinamic and numerical damping, second one using output format “ES17.9E3” and including aerodinamic and numerical damping, third one using output format “ES17.9E3” and including aerodinamic damping; all the simulation previusly mentioned shown off no azimuth angle dependency on damping ratios.

I think it is also important to mention the results obtained using eiganalysis on the azimuth averaged matrices, for both fixed frame(MBC applied) and rotational frame. First it’s worthwhile saying that next lines are related to BD+ED model, since ED model keeps on working perfectly.

Let’s start from rotational frame, hence where MBC transormation was not applied to matrices: it can be noticed that the results obtained by eiganalysis.m applied on matrices in rotational frame are significantly not consistent, except in case aerodynamic loads/damping are neglected (CompAero=0); in particular in rotatinal frame 6 natural frequencies below 1.5Hz are expected in this model, 3 related to 1st flapwise bending mode(1 for each blade) and 3 related to 1st edgewise bending mode; from the results only 5 frequencies with related damping ratios are found. For what concern the simulation performed with CompAero=0 the frequencies are correctly found but the damping ratios are not consistent.

For what concern the fixed frame, if aerodynamic is considered occurs that 1st flapwise collective mode is superimposed to the lower frequency 1st flapwise asymmetric mode. If aerodynamic is neglected the superimposition doesn’t occur but the damping ratios are not consistent compared with the ones related to azimuth angle distribution.

I hope these results are useful for the diagnostic of the problem and I hope the explanation was clear. In the file .xlsx are included all the results related to the simulations mentioned above with comments and recall to the models implied.

Thanks for your important support.

Simone Covre
ExcelReply5.xlsx (96.8 KB)

Dear Jason
I am using openfast 2.4 , when i run in matlab, it has an error as below
Error reported by S-function ‘FAST_SFunc’ in ‘Test01_SIG/FAST Nonlinear Wind Turbine/S-Function’:
FAST_InitializeAll:FAST_Init:FAST_ReadPrimaryFile:Invalid numerical input for file
“…\5MW_OC3Spar_DLL_WTurb_WavesIrr\5MW_OC3Spar_DLL_WTurb_WavesIrr.fst” occurred while trying to
read NLinTimes.

But it works well with CMD , i never change the fst file .
How can i solve it?

Best Regards!

Dear Rui,

Did you recompile the OpenFAST S-Function for v2.4 (this precompiled binary has not been provided)?

Regardless, as with any input file processor error, I would recommend enabling the Echo option to debug.

Best regards,

Dear Simone,

OK, thanks for the thorough analysis. Is this a proper summary of your results?:

There are two main concerns:

  • The structural damping of the ED+BD model shows azimuth dependence in the fixed frame of reference (after applying MBC3), particularly for the flapwise mode.
  • The damping obtained after azimuth averaging is not the same as the mean of the azimuthally varying damping.

General observations:

  • ED alone shows the expected behavior.
  • The problems in ED+PD are not influenced by numerical damping (the conclusions are the same with rhoinf = 0 or 1).
  • A precision of “ES17.9E3” is converged (going to a higher precision of “ES20.12E3” does not change the conclusions).
  • The problems are exacerbated when aerodynamics/damping is included.
  • The problems go away when structural damping in BD is eliminated (mu=0).

Because the problem appears to be tied to the structural damping, I would suggest eliminating aerodynamics/aerodynamic damping in further debugging, at least to start.

I would be curious if you see the same azimuth-dependence on damping when linearizing with zero rotational speed (again without aerodynamics/aerodynamic damping).

Best regards,

Dear Jason,

Your summary was perfect.

I performed the suggested simulations, in particular I neglected aerodynamics(CompAero=0) and I set RotSpeed=0 rpm, BlPitch(1,2,3)=0 deg, the incoming wind was a steady wind featured by 14 m/s speed. Since RotSpeed=0 rpm I performed 16 different simulation to cover the entire azimuth angles distribution. In BD module both structural and numerical damping were included.

The results are included in the attached .xlsx file and they have same layout of previous files. From the results it is clear that in ED+BD model the damping ratios keep on depending on the azimuth angle and the ones related to azimuth average show same inconsistencies found in previous reply results.

Thanks for your important help.

Simone Covre
ExcelReply6.xlsx (18.2 KB)