Control for scaled-down 5MW NREL

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 (http://forums.nrel.gov/t/gain-scheduling/286/1) 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: http://forums.nrel.gov/t/software-releases/522/1.

Related questions have been discussed on this forum e.g.:
http://forums.nrel.gov/t/imposing-estimated-motion-to-bottom-of-the-tower/1361/2
http://forums.nrel.gov/t/rigid-wind-turbine-with-harmonic-motions/685/6

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,

Dear Dr. Jason

Thank you very much for your time and quick reply.

I am using FAST V8. After learning the code, I plan to modify FAST_Prog.f90 of FAST V8 according to the Labview interface of FAST V7. In FAST_Prog.f90, I will replace FrcT0Trb and MomX0Trb by the measured tower-base loads and output the OtherState%AllOuts( PtfmTDxi、 PtfmTDyi、 PtfmTDzi、PtfmRDxi、 PtfmRDyi、 PtfmRDzi) as platform displacement into the Labview program. Do you think this modification is correct and feasible?

Best regards,
Gary.Wu

Dear Gary,

Just a few comments:

  • If you are going to upgrade a version newer than FAST v7 to reintroduce the Labview interface, I would start with OpenFAST (github.com/OpenFAST/openfast), where all active FAST development is taking place. If you are already familiar with FAST v8, then it should not be hard to upgrade to OpenFAST, which is based on, and contains improvements relative to, FAST v8.
  • FrcT0Trb and MomX0Trb are local variables within the ElastoDyn module, and so, are not available to the the FAST glue code, including FAST_Prog.f90. These variables also contain the entirety of applied and inertial forces of the rotor + drivetrain + nacelle + tower, which it doesn’t sound like what you want.
  • It sounds like what you want is to replace the aerodynamic applied loads computed by AeroDyn with those derived from the wind tunnel. If these loads will be expressed as a single lumped load (3 forces, 3 moments), you should set the input to the ElastoDyn module (HubPtLoad%Force, HubPtLoad%Moment on the HubPtLoad input mesh).
  • I’m not sure what you mean by “OtherState%AllOuts”; I don’t see that such a variable exists in the source code. What do you plan to physically model in the wind tunnel, e.g., is tower flexibility considered? If so, I would expect that you’d want to output the motions at the tower-top rather than the motions at the tower-base/platform anyway.

Best regards,

Dear Dr. Jason

I will also try to start with OpenFAST to develop the wind tunnel test method. Now I have understood the meaning of FrcT0Trb and MomX0Trb clearly based on your guidance. In the wind tunnel test, tower flexibility will not be considered. In this case, I think there exists a kind of fixed relationship between the motions at the tower-top and tower-base. I will ask you for help again when encountering new problems.

Thank you very much for the great guidance.

Best Regards
Gary.Wu