OpenFAST linearization issues

[RBergua]

Bug description Performing an eigenanalysis after running an OpenFAST linearization results in the wrong system eigenfrequencies for a rigid system.

To Reproduce I'm illustrating this with a very easy example: the system considered is a rigid floating wind turbine modeled as a lumped mass with inertia (ElastoDyn) and a stiffness matrix in HydroDyn representative of the hydrostatics and the mooring system.

For simplicity, the system in ElastoDyn only includes dummy blades and a dummy tower. Both are rigid (DOFs are disabled) and can be considered massless. The mass properties are given by a lumped mass and inertia defined by the platform at the (0,0,0). See below for reference:

In HydroDyn, only one dummy potential flow body is included. The potential flow files (e.g., *.hst, *.1., and *.3) only contain zeros and are used to enable the use of the stiffness matrix (AddCLin):

For reference, the HydroDyn model doesn't contain any member. No added mass is considered either and there is no damping. Since HydroDyn is used, SeaState (still water conditions) is also enabled within OpenFAST.

From this system, I expect 6 elastic modes: surge, sway, heave, roll, pitch, yaw.

However, when I run the linearization of the system and I perform the eigenanalysis, I only get 5 modes and one of them corresponds to rigid body motion (0 Hz). The eigenfrequencies found are as follows: 0.0000 Hz 0.0158 Hz 0.0164 Hz 0.0341 Hz 0.0395 Hz

The system is simply a diagonal mass matrix (6x6) and a diagonal stiffness matrix (6x6). I should get 6 elastic modes.

Expected behavior Given the very easy nature of the system, the eigenfrequencies can be computed analytically as f = (1/(2pi))*sqrt(K/M).

Analytical eigenfrequencies: Surge: 0.0158 Hz (observed in OpenFAST) Sway: 0.0164 Hz (observed in OpenFAST)
Heave: 0.0720 (missing/wrong in OpenFAST) Roll: 0.0354 Hz (missing/wrong in OpenFAST) Pitch: 0.0359 Hz (missing/wrong in OpenFAST) Yaw: 0.0376 Hz (missing/wrong in OpenFAST)

OpenFAST Version The same issue has been observed with OpenFAST v 4.0.2, OpenFAST dev, and OpenFAST with tight coupling.

For reference, the linearization is performed at t = 0 s in OpenFAST and the output format has been changed to "ES17.9E2" to increase the precision.

[jjonkman]

@RBergua -- Thanks for reporting this issue. Just a few comments/requests:

• Can you share the full set of OpenFAST input files to reproduce the problem?
• What is Gravity set to? Given that you are linearizing at t = 0, the model is likely not in static equilibrium given that PtfmMass * Gravity /= AddF0(3) when Gravity = 9.80665 m/s^2, which can impact the system stiffness. To get around this issue, you can set Gravity and AddF0(3) both to zero.
• Can you confirm that the stiffness matrix (K) from HydroDyn derived through the OpenFAST linearization analysis is correct (equaling AddCLin). You should be able to find this in HydroDyn's D matrix.
• Can you confirm that the mass matrix (M) from ElastoDyn derived through the OpenFAST lineariation analysis is correct? You should be to see this in ElastoDyn's D matrix, as discussed the following topic on our forum: https://forums.nrel.gov/t/openfast-2nd-order-linearization/2249.
• Can you confirm that the full-system state A matrix derived through the OpenFAST linearization analysis is correct, including -M^-1 * K in the lower-left quadrant?

Best regards,

[RBergua]

1. I have shared the OpenFAST input files with you via a Box folder.

2. By removing the gravity acceleration and the external force representative of the buoyancy (AddF0(3)), I get one elastic mode (heave) instead of the rigid body motion of 0 Hz that I had before. However, still one elastic mode is missing (I only get 5 out of 6 modes in the system) and 2 eigenfrequencies are not the expected ones. See below for reference: 0.0158 Hz 0.0164 Hz 0.0341 Hz 0.0395 Hz 0.0720 Hz

3. The HdroDyn's D matrix equals the stiffness matrix defined in AddCLin, with the caveat that the signs are wrong. I get this same opposite sign phenomena when linearizing HydroDyn as standalone (and it has been like that since the beginning of the linearization capability in OpenFAST).

4. The ED output file from the linearization is available here: VolturnUS_linearization_without_hydrostatic_OpenFAST_v5.1.ED.txt

I requested the platform accelerations as output in ElastoDyn. Accordingly, I have the platform forces and moments as inputs (columns 37-42) and the accelerations as outputs (rows 1-6). To get the mass matrix, I should compute the inverse of this 6x6 submatrix (see below). However, this submatrix is singular (the determinant is 0).

Looking at the outputs from eigenanalysis of the OpenFAST linearization, I only get the proper translations (surge, sway, heave) that correspond to the upper-left 3x3 submatrix. In the lower-right 3x3 submatrix I have the problems with one coefficient in the diagonal being 0. I guess this is the mode that I don't even find (that's why I get 5 modes instead of 6). And somehow, the other 2 modes related to rotations are also off...

Am I missing something?

[RBergua]

For completeness... If I compute the expected lower-left quadrant analytically (M^-1*K), I have:

And what I have from the OpenFAST linearization (A matrix) is:

As can be observed, the 3x3 submatrix shown in red has issues. This is aligned with the problems observed in the D matrix (see previous message) that relates to the mass matrix.

[jjonkman]

Thanks, @RBergua,

In a side conversation with @luwang00, it sounds like OpenFAST linearization is producing the correct result for this case (when Gravity = AddF0(3) = 0) when using OpenFAST v4. So, this linearization issue may be an issue unique to v5 (the tight coupling branch). Probably something for @deslaughter to look into.

Regarding HydroDyn's D matrix, the D matrix should show -AddCLin because HydroDyn will compute the hydrodynamic force on the right-hand side of M * qdd = F, that is, F = -AddCLin * q in this case.

Best regards,

[RBergua]

After talking to @luwang00, I realized that the issue was indeed with the platform pitch having an initial angle of 6 deg.

For reference, OpenFAST v4.0.2, OpenFAST dev, and OpenFAST tight coupling return the same results: 0.0158 Hz 0.0164 Hz 0.0354 Hz 0.0359 Hz 0.0376 Hz 0.0720 Hz

And quite interesting fact: Linux systems handle properly the eigenanalysis when the system is not in perfect equilibrium condition (e.g., you can have initial platform pitch = 6 deg or gravity enabled with the vertical force enabled and not having perfectly equal magnitude) and you will still get 6 modes that make physical sense. However, Windows systems will not handle this properly (at it has been illustrated in this thread).

This has been also checked by @luwang00. So, Linux and Windows systems will return different eigenfrequencies if there is any imbalance in the system (e.g., heave or pitch DOFs in the cases studied).

[jjonkman]

Regarding Linux versus Windows, is the state matrix "A" changing or is the eigenanalysis changing for the same "A" matrix?

[luwang00]

Hi @jjonkman, the "A" matrix is different on Linux. For instance, on Linux, the "A" matrix in the .lin file has all expected entries even if gravity and AddF0(3) are not zeroed. This leads to the correct list of eigenfrequencies. With 6-deg pitch, the "A" matrix on Linux becomes slightly different compared to the one without pitch offset, as is to be expected, which leads to a slightly different but similar set of 6 eigenfrequencies.

[xiaolan0195]

Dear Dr.Jonkman,

I would like to ask about the issue of linearlization. I want to linearize based on this formula（As shown in the above figure）. The selected state variables x are [rotor speed, generator speed, torque angle θ, electromagnetic torque Tg, pitch angle], the control variable u is [electromagnetic torque Tg, pitch angle], and the output variable y is [rotor speed, power]. May I ask if this can be linearized using OpenFAST to obtain the corresponding ABCD state space equation matrix?

Best regards,

[jjonkman]

Dear @xiaolan0195,

I'm not sure I understand what your last two equations are, but otherwise OpenFAST can account for these equations in its linearization, with the generator and drivetrain DOFs enabled in ElastoDyn. Note that in ElastoDyn, the rotor speed, omega_r = qd_g + thetad and the generator speed, omega_g = qd_g * N_g, where qd_g is the generator speed state expressed on the low-speed shaft side of the generator and the rest of the variables match yours.

Best regards,

[xiaolan0195]

尊敬的 ，

我不确定我是否理解您的最后两个方程是什么，但除此之外，OpenFAST 可以在其线性化中解释这些方程，在 ElastoDyn 中启用发电机和传动系统自由度。请注意，在 ElastoDyn 中，转子速度 omega_r = qd_g + thetad 和发电机速度 omega_g = qd_g * N_g，其中 qd_g 是发电机低速轴侧表示的发电机速度状态，其余变量与您的变量匹配。

此致敬意

Thanks！Tg is the electromagnetic torque of the wind turbine,it is the first-order system of electromagnetic torque. The last equations is the first order pitch servo system .It is the first-order system of pitch.I tried running the linearized fst file and obtained the following values. The image shows the description of the relevant variables within Xop.The 20th row 'ED First time derivative of Variable speed generator DOF (internal DOF index = DOF_GeAz), rad/s' is Generator speed(omega_g)；The 21th row 'ED First time derivative of Drivetrain rotational-flexibility DOF (internal DOF index = DOF_DrTr), rad/s' is Drivetrain torsion rate(thetad); So the rotor speed can through omega_g get,omega_r=omega_g/ N_g+thetad? . Is my understanding correct?

[jjonkman]

Dear @xiaolan0195,

The first-order systems for torque and pitch servo-dynamics will not be included in the linearized model generated by OpenFAST because they are not inherent in the OpenFAST model, but rather, can only be included if you implement them yourself in a user-defined controller.

As I said in my prior post, the generator DOF (DOF_GeAz) in ElastoDyn (ED) expresses the rotation of the low-speed shaft at the entrance to the gearbox, so, omega_g = qd_g * N_g. And omega_r = gd_g + thetad.

Best regards,

[deslaughter]

@luwang00 @RBergua I'm looking into the difference in linearization between Windows and Linux, could you confirm that OpenFAST was compiled in double precision on both systems?

[RBergua]

@deslaughter The runs in Windows with OpenFAST v4.0.2 and v5 were performed with a compile in single precision:

[deslaughter]

I think using single precision for linearization is a large part of the issue. I was able to get a non-singular mass matrix when using double precision, which is the default when compiling with CMake. This may explain the differences between Windows and Linux. With 0 gravity and 0 initial platform pitch, I got the following frequencies in Hertz

0.01579178, 0.01642382, 0.03644564, 0.0369004 , 0.03755791, 0.0719742

Do these look more like what you expected?

[luwang00]

@deslaughter The Linux version used is double precision I believe.

[RBergua]

@deslaughter The expected frequencies in Hertz are: 0.0158, 0.0164, 0.0354 0.0359, 0.0376, 0.0720.

Note that the roll inertia is 4.7674E10 kgm^2 and the roll stiffness is 2.3647E9 Nm/rad. Both magnitudes are placed at the (0,0,0). So, the expected roll frequency is: f = (1/(2*pi))*sqrt(2.3647E9/4.7674E10) = 0.0354 Hz

Similarly, the pitch frequency is 0.0359 Hz.

Somehow, roll and pitch frequencies seem to be off in your outputs?

[deslaughter]

@RBergua, thanks for providing the expected frequencies. I still had a non-zeroAddF0(3) which was giving incorrect frequencies. Now when using the tight coupling branch compiled in double precision on Windows, I get the expected frequencies. Does that resolve the outstanding issue?

[jjonkman]

@deslaughter -- So, is the recommendation that we need to use double precision when performing a linearization analysis, and perhaps also when using tight coupling? I thought we had tried to ensure all variables that require double precision use double precision even when compiling in single precision. Perhaps there are still variables that require double that are not yet declared that way when compiling in single precision?

[deslaughter]

@jjonkman I think the recommendation is to use double precision when performing a linearization analysis. All of the linearization glue code in tight coupling is in double precision so the round off is likely in ElastoDyn/HydroDyn. I can try to track down exactly where it's coming from.