Control for scaled-down 5MW NREL

Dear All,

I am a Polish student (Lodz University of Technology) currently working on my Bachelor Thesis where I analyze a concept of multi-rotor wind turbines for deep off-shore application. The idea behind it is to replace the 5MW NREL wind turbine with a set of 7 smaller rotors (which in total give the same swept area as the baseline turbine), compare their energy production and design appropriate frame.

I attempt to use FAST to evaluate performance of an individual rotor of MRWT. I already managed to scale down the baseline turbine using formulas presented by NREL and now I tackle the problem of turbine control. I wonder whether I can apply control system implemented for 5MW NREL system, after making certain adjustments.

I studied “Advanced Control Design for Wind Turbines Part I: Control Design, Implementation, and Initial Tests”, I changed values like Ki, Kp, VS_RtPwr, VS_MaxTq in DISCON.f90 file but it does not work.

I am quite new in FAST, I have been working with it for no longer than 3 months. Could you tell me if my approach is even correct? Is it possible to adjust the control for 5MW NREL WT for scalled-down version of it (power seven times smaller). Is changing DISCON.f90 file enough? It is only one part of my Bachelor Thesis, that is why I do not want to develop the control system from the scratch.

I appreciate any help and guidance! If I did not put all necessary information concerning the problem, let me know. That is my first post here.

Regards,
Sebastian.Dziomdziora

Dear Sebastien,

Yes, I agree with your approach. You mention that it does not work, but I have hard time offering suggestions when you don’t provide further information on the problem.

After scaling the appropriate geometry, mass, and stiffness terms within the FAST input files, I agree that all you should need to change for the controller are the parameters in the DISCON controller.

Best regards,

Dear Jason,

thank you for your fast reply.

I scaled appropriate geometry, mass and stiffness terms as you scaled. For the purpose of the discussion let’s assume they are correct and I can move on to changing parameters in DISCON controller. The problem I am facing a problem computing values of Ki and Kp.

If I understood correctly, to get them I need sensitivity of aerodynamic power to the rotor collective blade-pitch angle ∂P/∂θ (I assumed natural frequency ω_φn and damping ratio ζ_φ to be 0.6 rad/s and 0.6-0.7 respectively).

In Table 7.1 (page 23, “Definition of a 5MW Reference Wind Turbine…”) you presented result for wind speed range from 11.4 to 25 m/s. Was the corresponding pitch angle found by simply running FAST simulation for different angles and selecting the one which yielded mechanical power the closest to rated?

You described the whole linearization procedure to obtain values of sensitivity in the ''Gain Scheduling" post (Gain Scheduling). I followed your procedure, I cannot find anywhere the D-Trnsmt matrix you mentioned. I opened the *.lin file with MBC matlab script called GetMats.m, then run cce.m script as well.

I am completely lost and would greatly appreciate if you could explain me the solution to the problem step-by-step. How the linearization for this purpose should be performed. I even downloaded the FAST v6.01 with AeroDyn v12.58 modified to apply frozen wake during linearization, however newer input files seem not to be compatible with the app.

Thank you in advance!

Regards,
Sebastian Dziomdziora

Dear Sebastian,

My Aug 24, 2010 post in the forum topic you referenced already contains a detailed procedure. To answer your questions:

  • Yes, the pitch angle that yielded mechanical power was found for each wind speed.
  • “D-Trnsmt” is the name of the D matrix (input transmission matrix) in the *.lin file.
  • I’ve answered the question regarding the input file compatibility in the aforementioned forum topic.

I hope that helps.

Best regards,

Dear Jason,

I followed all your guidelines and managed to calculate all necessary parameters for DISCON controller (or at least I thought I did). But it still does not work as it is supposed to do! Let me draft the current situation and I hope you will be able to help me ones again.

I am scalling down 5MW NREL baseline to be 7 times smaller (called later on 1MRA). The assumption was to obtain the same total swept area for seven 1MRA as for 5MW NREL. As result radius of 1MRA is 23.8m, Rated Mechanical Power 0.757kW, Rated Rotor Speed 32.02 rpm. I kept TSR, Cp, Rated Wind Speed, Rated Tip Speed the same.

Based on those input parameters I adjusted geometry concerning blades (AeroDyn.dat and blade.dat files), masses and it centers (using baseline WIndPact equations). I conducted linearization according to your suggestions in the post (Gain Scheduling) to calculate Ki and Kp.

Below I list all parameters in DISCON files that I changed:
PC_KI = 0.0072
PC_KK = 0.114598
PC_KP = 0.0169
VS_MaxRat = 2200.0
VS_MaxTq = 6771.82
VS_Rgn2K = 0.33415
VS_RtPwr = 782150.0
VS_SlPc = 22.45

Sensitivity ∂P/∂θ = -3.655E+03 kW/rad for θ=0 and θk = 6.565 degree

Ki = I_drivetrainΩω^2/(N_gear*(∂P/∂θ)) ;
Kp = 2I_drivetrainΩζω/(N_gear*(∂P/∂θ)) ;

where:

  • I_drivetrain = I_rotor (taken from .fsm file) + I_generator (I_generator for 5MW divided by 7)
  • ω = 0.6 rad/s, ζ = 0.7, N_gear = 36.65, Ω = 32.02 rpm

VS_MaxTq - set 10% above rated torque (obtained by dividing rated generator torque of 5MW by 7)
VS_RtPwr - set to 782kW instead of theoretical 757kW because simulations in FAST yielded such a value

VS_Rgn2K = 0.5ρPIR_1MRA^5(Cp_max/λ^3)(1/N_gear^3), where: Cp_max = 0.482, λ=7.55
VS_SlPc = 100
(Ω2/Ω0(Q=0) - 1), where Ω2 is 99% of Rated generator speed (1173.7), Ω0(Q=0) - speed where torque is 0 (found for region 2 1/2 contained between Ω1 (94% of Rated generator speed) and Ω2).

Everything seemed okey, however with Servo Mode turned on, for Rated Wind Speed 11.4 m/s, the Rotor Speed instead of being held at level of Rated (32.02 rpm) it drops to 19.90 rpm. The same with power from 782 kW to around 430 kW. The same situation happened when step wind profile was applied for speeds higher than rated; pitch actuation never takes place because required level of torque is not reached.

I have run out of ideas what can be fixed… Is the controller that sensitive, that even small errors in parameters lead to that kind of behavior? Do you have any idea what may be the cause?

Thanks for your help in advance!

Best regards,
Sebastian Dziomdziora

Dear Sebastian,

I haven’t reviewed all of the details in your post, but your overall approach sounds correct. I’m sure there is just a few details you’ve missed.

I noticed in your treatment of I_drivetrain that you missed the effect of the gearbox ratio–the equation should be I_drivetrain = I_rotor + N_generator*N_gear^2.

You didn’t explicit say it, but there are a few other parameters in DISCON.f90 file that would need to be changed for a different size turbine, including PC_RefSpd, VS_CtInSp, VS_Rgn2Sp, and VS_RtGnSp.

Best regards,

Dear Dr. Jason
If I want to scale down NREL 5MW to one sixty-fourth, how do I define PC_KI, PC_KK, and PC_KP in my new pitch controller?
Does it have any relationship with scaling proportion or I must use the linearization according to your suggestion in the post (viewtopic.php?f=4&t=394) to calculate Ki and Kp?
Besides, I have to use FAST_v6.01_AeroDyn_v12.58_WithFrozenWakeAssumptionDuringLinearization.exe, is the thought correct?
Thanks for your help.

Best Wishes
Jason.Lai

Dear Jason,

I have not done this myself, but if you are applying consistent scaling laws to all FAST variables, there is likely a way to also scale the pitch control gains accordingly (look at the PI equations to see how the variables are used).

You could also rederive the gains using linearization as you propose. However, you don’t need to use FAST v6.01; linearization with frozen wake is enabled in FAST v8.16 and OpenFAST v1.0 as mentioned in the aforementioned forum topic.

Best regards,

Dear Dr. Jason
I refer the paper “FAST Modularization Framework for Wind Turbine Simulation: Full-System Linearization.”
For each wind speed/blade-pitch angle, one rotor (360˚) revolution is simulated, and an OP is set and the linearized full-system input-transmission matrix (D ) is computed for each azimuth angle of the rotor in 36 steps of 10˚.

To date, I have done the azimuth average ( 36 Linearization times ie every 10 degree azimuth ) in FAST v8.16 and I have done MBC postprocess.
I read the paper and search some information in this forum but I still don’t know how to find dP/dtheta in my average D matrix?

According to your aforementioned reply: *Run FAST and–from the .lin linearization output file–record the value of the D-Trnsmt matrix associated with the “RotPwr” output and the collective blade pitch control input. This is dP/dTheta.
Could you describe it in more detail?

I want to verify the result of dP/dtheta versus pitch angle.
Thanks for your help.

Best Regards
Jason. Lai

Dear Jason,

In FAST v8.16, I would include AeroDyn output “RtAeroPwr”, which is the aerodynamic power (aerodynamic applied torque times rotor speed). Setting LinInputs = 1 in the FAST primary input file will mean ElastoDyn extended input “BlPitchComC” will be one of the inputs in the linear model. Thus–after running the FAST v8 simulation–the element of the full-system D matrix in the FAST linearization (*.lin) file associated with the RtAeroPwr row and BlPitchComC column contains the derivative of RtAeroPwr with respect to BlPitchComC i.e. dP/dTheta.

I hope that helps.

Best regards,

Dear Dr. Jason
Thanks for your suggestion and I have verified the result of dP/dtheta versus pitch angle.
I can scale down the model in FAST v7.02d successfully.
However, I meet the other question about scale-down NREL 5MW with a monopile substructure in FAST v8.16.
I want to get the result of dP/dtheta versus pitch angle in the scale-down model and then try to calculate the PC_KI, PC_KK and PC_KP.
Could you give me some comments about how to debug it?

[code]********************************************************************************


FAST (v8.16.00a-bjj, 27-Jul-2016)

Copyright (C) 2016 National Renewable Energy Laboratory

Running FAST (v8.16.00a-bjj, 27-Jul-2016), compiled as a 32-bit application usi
ng double
precision
linked with NWTC Subroutine Library (v2.09.00, 23-Jul-2016)

Heading of the FAST input file:
FAST Certification Test #19: NREL 5.0 MW Baseline Wind Turbine with OC3 Monop
ile RF
Configuration, for use in offshore analysis

Running ElastoDyn (v1.04.00a-bjj, 26-Jul-2016).

Running AeroDyn (v15.03.00, 27-Jul-2016).

Running AirfoilInfo (v1.01.00a-bjj, 5-Apr-2016).

Running BEM (v1.02.00, 29-Jun-2016).
Warning: Turning off Unsteady Aerodynamics because C_nalpha is 0. BladeNode = 1
, Blade = 1
Warning: Turning off Unsteady Aerodynamics because C_nalpha is 0. BladeNode = 2
, Blade = 1
Warning: Turning off Unsteady Aerodynamics because C_nalpha is 0. BladeNode = 3
, Blade = 1
Warning: Turning off Unsteady Aerodynamics because C_nalpha is 0. BladeNode = 4
, Blade = 1
Warning: Turning off Unsteady Aerodynamics because C_nalpha is 0. BladeNode = 1
, Blade = 2
Warning: Turning off Unsteady Aerodynamics because C_nalpha is 0. BladeNode = 2
, Blade = 2
Warning: Turning off Unsteady Aerodynamics because C_nalpha is 0. BladeNode = 3
, Blade = 2
Warning: Turning off Unsteady Aerodynamics because C_nalpha is 0. BladeNode = 4
, Blade = 2
Warning: Turning off Unsteady Aerodynamics because C_nalpha is 0. BladeNode = 1
, Blade = 3
Warning: Turning off Unsteady Aerodynamics because C_nalpha is 0. BladeNode = 2
, Blade = 3
Warning: Turning off Unsteady Aerodynamics because C_nalpha is 0. BladeNode = 3
, Blade = 3
Warning: Turning off Unsteady Aerodynamics because C_nalpha is 0. BladeNode = 4
, Blade = 3

Running InflowWind (v3.03.00, 26-Jul-2016).
Opening InflowWind input file: .\5MW_Baseline/NRELOffshrBsline5MW_InflowWind_S
teady8mps.dat
Processed 2 records of uniform wind data from ‘.\5MW_Baseline/Wind/uniform_1m
s.wnd’

Running ServoDyn (v1.06.00a-bjj, 26-Jul-2016).

Running ServoDyn Interface for Bladed Controllers (using Intel Visual Fortran f
or Windows,
14-Oct-2015).

Running HydroDyn (v2.05.01, 27-Jul-2016).
Setting WaveTMax to 0.0 since WaveMod = 0
Calculating second order difference frequency wave kinematics.
Calculating second order sum frequency wave kinematics.

Running SubDyn (v1.03.00, 11-Apr-2016).
Calculating Internal Modal Eigenvectors
Timestep: 0 of 78.75 seconds.

FAST_Solution0:CalcOutputs_And_SolveForInputs:SolveOption2:SrvD_CalcOutput:Runni
ng with torque
and pitch control of the NREL offshore 5MW baseline wind turbine from DISCON.dll
as written by J.
Jonkman of NREL/NWTC for use in the IEA Annex XXIII OC3 studies.
CalcOutputs_And_SolveForInputs:SolveOption1:ED_SD_HD_BD_Orca_InputOutputSolve:SD
_CalcOutput:Small
angle assumption violated in SUBROUTINE SmllRotTrans() due to a large UL input a
ngles. The
solution may be inaccurate. Simulation continuing, but future warnings from Smll
RotTrans() will
be suppressed.
Additional debugging message from SUBROUTINE SmllRotTrans():

Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 5, Blade = 1
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 6, Blade = 1
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 7, Blade = 1
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 8, Blade = 1
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 9, Blade = 1
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 10, Blade = 1
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 11, Blade = 1
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 12, Blade = 1
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 5, Blade = 2
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 6, Blade = 2
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 7, Blade = 2
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 8, Blade = 2
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 9, Blade = 2
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 10, Blade = 2
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 11, Blade = 2
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 12, Blade = 2
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 13, Blade = 2
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 5, Blade = 3
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 6, Blade = 3
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 7, Blade = 3
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 8, Blade = 3
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 9, Blade = 3
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 10, Blade = 3
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 11, Blade = 3
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 12, Blade = 3
Warning: Turning off Unsteady Aerodynamics due to high angle-of-attack. BladeNod
e = 13, Blade = 3

FAST_Solution:FAST_AdvanceStates:AD_UpdateStates:BEMT_UpdateStates(node 18, blad
e
1):UA_UpdateStates:UA_UpdateDiscOtherState:ComputeKelvinChain:Mach number exceed
s 0.3. Theory is
invalid. This warning will not be repeated though the condition may persist.
FAST_Solution:CalcOutputs_And_SolveForInputs:SolveOption2:AD_CalcOutput:BEMT_Cal
cOutput(node 17,
blade 1):Compute_UA_AirfoilCoefs:UA_CalcOutput:Mach number exceeds 1.0. Equation
s cannot be
evaluated.
CalcOutputs_And_SolveForInputs:SolveOption1:ED_SD_HD_BD_Orca_InputOutputSolve:SD
_CalcOutput:
Angles in GetSmllRotAngs() are larger than 0.4 radians.
ED_SD_HD_BD_Orca_InputOutputSolve:HydroDyn_CalcOutput: Angles in GetSmllRotAngs(
) are larger than
0.4 radians.
ED_SD_HD_BD_Orca_InputOutputSolve:SD_CalcOutput: Angles in GetSmllRotAngs() are
larger than 0.4
radians.
ED_SD_HD_BD_Orca_InputOutputSolve:HydroDyn_CalcOutput: Angles in GetSmllRotAngs(
) are larger than
0.4 radians.

FAST encountered an error at simulation time 6.25000E-04 of 78.75 seconds.
Simulation error level: FATAL ERROR

Aborting FAST.

[/code]
Thanks for your help again.

Best Wishes
Wen-Jeng

Dear Wen-Jeng,

My guess is your model is going numerically unstable because you’ve scaled down the model (effectively increasing the natural frequencies) without adequately reducing the time step. Warnings regarding a “small angle approximation violation,” warnings regarding “Mach number exceeds 1.0,” a simulation crash, or very large deflections that occur in the time series are good signs of a numerical instability. Or use “Search…” in the upper right of the forum to search for “small angle assumption” to see how we’ve solved similar problems in the past.

Best regards,

Dear Dr. Jason
Thanks for your suggestion and I have resolved my problem.
Besides, if I can get a measurement of model 6dof motions, how do I replace 6dof motions calculated in fast v8.16 by those?
Could you describe how to do it in detail?

Best Wishes
Jason. Lai

Dear Jason,

You cannot directly specify 6-DOF motion of the platform in FAST without customization of the source code.

For the 3 platform translation modes, you can use FAST v7’s Seismic module to prescribe time histories of motion–please see my post dated Aug 29, 2012 in the forum topic found here for more information: Software Releases.

Related questions have been discussed on this forum e.g.:

Best regards,

Dear Dr. Jason
I want to customize the source code to specify the 6-DOF motion of the platform in FAST.
And then I will verify the 6-DOF motion between FAST outputs and my scaling-down experiment results.
They have to match very well.
Is the thought correct?

Best Wishes
Jason. Lai

Sounds reasonable to me.

Hi Jason

I am currently trying to do the same as Jason.Lai, in real-time running. I understand that it is not possible to feed the platform into the openfast without customizing the code. I am planning to customize the code so that it can take real-time measured platform motion. Could you please give me some guidance on which part of the code would need to be modified? Any relevant suggestions are appreciated.

Best Regards
Saishuai.Dai

Dear Saishuai,

There are two ways to do this in my mind, as mentioned in the links above.

The first (and likely easiest) method is to increase the platform mass/inertia in ElastoDyn (m), convert the prescribed motion (a) into a prescribed load (F), and calculate the load such that F = ma, which is how Seismic solves the problem (the platform mass/inertia is increased so that the prescribed load (F) dominates other loads applied to the platform, i.e., loads transmitted from the tower, such that the motion (a) will follow what is prescribed). In this case, you could implement the logic for the force (F) in OpenFAST module ExtPtfm

The second approach talks about changing the structural module (ElastoDyn) so that the platform motion in prescribed rather than calculated. This approach would require rearrangement of the equations of motion in ElastoDyn, which is likely a large undertaking.

Best regards,

Dear Dr. Jonkman

I am trying to develop a wind tunnel test method for floating wind turbines using real-time simulation with FAST. In the RT simulation, aerodynamic forces are from wind tunnel measurement in Labview program rather than Aerodyn. Then the platform motion caculated by FAST should output into the Labview program. Now I have a NI RT sysytem, and the aerodynamic forces can be obtained from the Labview program. My two question are 1) how to modify FAST V8 to replace Aerodyn with the external load input in real time 2) is there Labview interface for FAST V8? Could you please give me some suggestions to accelerate the modification?

In addition, could you please send the " Unofficial FAST Theory Manual" to 35050310@163.com

Thanks very much.

Best Regards
Gary.Wu

Dear Gary,

FAST v8 and OpenFAST do not currently support an interface to Labview. This capability was last supported in FAST v7, but is something we’d like to add to OpenFAST in the future (dependent on funding).

You’d have to modify the FAST (or OpenFAST) source code to replace AeroDyn with a force computed externally. The approach to do this will be dependent on which version of FAST/OpenFAST you are using. Perhaps you want to revert to FAST v7 to take advantage of the existing Labview interface? Or perhaps you want to work to reintroduce the Labview interface to the newer versions of FAST, in which case I would recommend starting with OpenFAST. Please clarify which path you intend to pursue.

I’ve also e-mailed you the "unofficial FAST theory manual.

Best regards,