Hydrodynamic implementation

Wind & Water Wave Energy
41

Dear Jason,
Starting from the same equations of Riccardo (corrected as you suggested in the last posts), I am investigating the differences in the transition from AeroDyn = 0 to AeroDyn = 2 with regards to platform movements.

Using the reported forces and moments output of FAST, I calculated positions in the 6 DOFs. When AeroDyn = 0, I get the same position trends in the 6 DOFs to the FAST outputs. When I switch to AeroDyn = 2, however, the only positions consistent with FAST are x, z, ry, while the other 3 (y, rx, rz) have a totally different trend, and with greater oscillation amplitudes.
I was wondering if there is any other contribution in moving from AeroDyn = 0 to AeroDyn = 2

Thanks for your attention,
Best regards,
Lorenzo

42

Dear Lorenzo,

I’m sorry, but I’m not really sure I understand your question. By “AeroDyn = 0 or 2”, are you referring to “CompAero = 0 or 2”? And I’m not sure what you are comparing FAST to (your own model)? And I’m not sure what you mean by “any other contribution”?

Best regards,

43

Dear Jason,
Sorry, I was referring to CompAero = 0 or 2.
Adding all the forces and moments from outputs of FAST, then dividing by the mass matrix, I get the accelerations from which the positions in the 6 DOFs. Finally, I compare the positions thus obtained with the FAST outputs.
When I have CompAero = 0 the graphs are superimposed in the 6 DOFs, while with CompAero = 2 I have that the graphs of y, rx and rz do not replicate the FAST outputs well.
Since in both cases I use the same equation (in fact, I also take into account, in the case of CompAero = 0, the initial transient of the RotTrq, etc.), I wondered why in the case of CompAero = 2 it does not obtain equal trends in all DOFs but only in x, z, ry.
This is why I asked you if there are any other components of the equation, which I have neglected, that can motivate these inconsistencies in the positions y, rx, rz.
I hope I was clearer.
Best regards,
Lorenzo.

44

Dear Lorenzo,

OK, I think I understand now. But I’m not sure I have anything to add beyond what I’ve already discussed with Riccardo.

Best regards,

45

Dear Jason,
Thanks for checking.
I have tried to investigate further what can cause the differences between my outputs and those of FAST, and I have noticed that the cause is the resultant moment Mx.
In fact, I compared, for each DOF, that the equation was verified:
M * a = Fhydro + Fturbine + Fmoor - Frestoring - Fnonlinear
thus obtaining from the product M * a the resulting Fx, Fy, Fz, Mx, My, Mz which I had to obtain from the sum of the contributions to have the same accelerations as FAST.
To do this I used:

  • M: variable mass matrix with the movement of the platform
  • a: acceleration outputs of FAST (‘PtfmTAxi’, ‘PtfmTAyi’, ‘PtfmTAzi’, ‘PtfmRAxi’, ‘PtfmRAyi’, ‘PtfmRAzi’)
  • Fhydro: outputs FAST (‘HydroFxi’, ‘HydroFyi’ ‘HydroFzi’, ‘HydroMxi’, ‘HydroMyi’, ‘HydroMzi’)
  • Fturbine: forces at the base of the turbine tower, obtained starting from the aerodynamic forces / moments on the rotor (‘RtAeroFxh’, ‘RtAeroFyh’, ‘RtAeroFzh’, ‘RtAeroMxh’, ‘RtAeroMyh’, ‘RtAeroMzh’), transported in the fixed reference of the nacelle and then transported to the base of the tower.
  • Fmoor: starting from the forces in the 3 directions for each line obtained by MAP ++, the moments were obtained
  • Frestoring: where the contribution on Mx is given by m * x_g * ‘PtfmYaw’ - m * z_g * ‘PtfmRoll’ where m = total mass of the system, (x_g, y_g, z_g) coordinates of the center of mass with respect to the origin of the inertial reference system
  • Fnonlinear: as regards the moments, I used the formula omega x I dot omega (with omega = [‘PtfmRVxi’, ‘PtfmRVyi’, ‘PtfmRVzi’] and I: inertia matrix 3x3)

As you can see from the figures shown, I obtain summation trends of forces-moments that accurately overlap with the corresponding trend M * a in all directions, except on Mx.

WhatsApp Image 2020-09-26 at 17.42.42.jpeg
WhatsApp Image 2020-09-26 at 17.42.42.jpeg1366×651 100 KB

WhatsApp Image 2020-09-26 at 17.42.41.jpeg
WhatsApp Image 2020-09-26 at 17.42.41.jpeg1366×651 110 KB

This makes me make some considerations:

  • given that I get 5 out of 6 results, it is unlikely that there is an error of transport of the forces-moments in the turbine, of the calculation of the moments of the moorings from the forces or of the way in which I considered the restoring moment or the variation of the mass matrix . If there was an error I should have problems in other directions as well.
  • for the above reasons I have asked you in the past if there were other contributions to consider, in particular some effect that acts above all on Mx (and therefore can cancel this difference which I now find).

I hope I have clearly explained the problem to you,
Best regards,
Lorenzo.

46

Dear Lorenzo,

There should not be anything unique in OpenFAST regarding Mx. I agree with your calculations, with just two comments:

  • As I mentioned before, the equation you are using for Frestoring is linear, whereas ElastoDyn considers some nonlinearity here based on Eq. (2) in the Wind Energy paper I referenced above.
  • You don’t say this, but presumably the inertia matrix, I, is also variable with platform movement.

Best regards,

47

Dear Jason,
Yes, the inertia matrix varies with the movement of the platform.
I have consulted the paper you suggested. I think I understand that Equation 2 replaces Euler’s angles. However we use Euler angles to move to the inertial frame of reference (0,0,0), so this term should already be considered by us.
Furthermore, I don’t understand how Euler angles affect the non-linearity of the restoring moment.
Where can I find information about it?

Best regards,
Lorenzo.

48

Dear Lorenzo

The restoring moment can be written as

F_restoring = -r_g x { 0; 0; -m*g }

where:

r_g = [ TransMat^T ]*{ x_g; y_g; z_g }

and TransMat^T is the transpose (inverse) of the orientation matrix documented in Eq. (2).

Best regards,

49

Dear Jason,
To better understand the reason for the divergence of the previously analyzed results (see quote), I chose to carry out an analysis only on the hydrodynamic part, thus setting CompAero = 0 and defining the initial rotor speed = 0. In this way it is possible to neglect the loads at the base of the turbine throughout the simulation.

The equation of dynamics then becomes:
M * a = Fhydro + Fmoor - Frestoring - Fnonlinear
where each component is defined as in the quote, with the exception of the restoring moment where the changes you suggest have been inserted.
As I did earlier, I have compared forces-moments with M * a in all directions:

WhatsApp Image 2020-10-29 at 13.53.13 (1).jpeg
WhatsApp Image 2020-10-29 at 13.53.13 (1).jpeg1600×767 168 KB

WhatsApp Image 2020-10-29 at 13.53.13.jpeg
WhatsApp Image 2020-10-29 at 13.53.13.jpeg1600×767 203 KB

Comparing the results, I recorded the following errors on the rms values:
Fx: 2,5% Fy: 2,4% Fz: 0,4% Mx: 15,1% My: 6,7% Mz: -99%
As you can see, Mx has a small error especially in the first part of the simulation. The resultant of Mz, on the other hand, is of the order of ten while the corresponding M * a has values of hundreds: I cannot understand the reason given that HydroMz, MoorMz and NonlinMz have values of the order of unity, while the restoring moment is null in rz direction.
Is there any advice you can give me?
Thanks for the attention,
Lorenzo.

50

Dear Lorenzo,

I’m not sure. Do you see yaw motion in this simulation? If so, I would expect nonzero moments to show up in Fhydro, Fmoor, and Fnonlinear.

Another thing I notice is that the yaw moments are very small relative to the roll, and especially the pitch, moments. This suggests that that numerical round in the calculations may have more of an influence on the yaw moment.

Best regards

51

Dear Jason,
Yaw is very small (order of 10^-8 degrees), this justifies the low values of Mz in Fhydro, Fmoor, and Fnonlinear.
To compare, roll order is about 10^-5 degrees, from which we can see the proportion between Mx and Mz: in fact, Mx is of the order of 10^4, while Mz should be 10^1 (not 10^2, as obtained with M*a).
I agree with you that it could therefore be a numerical error which becomes more evident in Mz due to the low values. The most important thing for me is not to have forgotten any contributions.
Best regards,
Lorenzo.

52

Dear Jason,
To be completely sure that I have done everything right, I will rewrite for the last time, in as much detail as possible, the terms that make up the equation of hydrodynamics:
M * a = Fhydro + Fmoor - Fnonlinear -Frestoring
where is it:

  • M: variable mass matrix with the movement of the platform. It is defined starting from the mass matrix in the center of mass Mcg, whose value has already been confirmed by you in the topic OC3-Hywind RAOs - #25 by Lorenzo.Cottura.
    Then M = M_swl = TransMat^T^-1M_cgTransMat^-1 where
    TransMat =
    [ [ 1 0 0 0 -z_cg y_cg ];
    [ 0 1 0 z_cg 0 -x_cg ];
    [ 0 0 1 -y_cg x_cg 0 ];
    [ 0 0 0 1 0 0 ];
    [ 0 0 0 0 1 0 ];
    [ 0 0 0 0 0 1 ] ]
    and
    TransMat^T =
    [ [ 1 0 0 0 0 0 ];
    [ 0 1 0 0 0 0 ];
    [ 0 0 1 0 0 0 ];
    [ 0 z_cg -y_cg 1 0 0 ];
    [ -z_cg 0 x_cg 0 1 0 ];
    [ y_cg -x_cg 0 0 0 0 1 ] ]
    In these formulas, the coordinates of the center of mass r_g = [x_cg, y_cg, z_cg] are not fixed but vary with the movements of the platform, through the formula:
    r_g = [AngleMat^T ]{ x_cg; y_cg; z_cg }
    where { x_cg; y_cg; z_cg } is the initial position of the center of mass and AngleMat^T is the transpose (inverse) of the Eulerian matrix (or of the orientation matrix documented in Eq. (2) in your 2009 Wind Energy paper: onlinelibrary.wiley.com/doi/abs/10.1002/we.347).
    In this way a mass matrix M = M_swl 6x6 is obtained, variable and arranged as in the figure.
    Schermata 2020-10-30 alle 10.35.51.png
    It is also recalled that x_g, y_g, z_g in the figure, as well as the moments of inertia, vary from instant to instant depending on the arrangement of the platform. From now on, references to M will be references to the variable mass (example M (1,5) = m    z_g)

  • a: acceleration outputs of FAST (‘PtfmTAxi’, ‘PtfmTAyi’, ‘PtfmTAzi’, ‘PtfmRAxi’, ‘PtfmRAyi’, ‘PtfmRAzi’)

  • Fhydro: outputs FAST (‘HydroFxi’, ‘HydroFyi’ ‘HydroFzi’, ‘HydroMxi’, ‘HydroMyi’, ‘HydroMzi’)

  • Fmoor: starting from the forces in the 3 directions for each line obtained by MAP ++ ( example Fx_moor = Fx[1] + Fx[2] + Fx[3] + Fx[4] + Fx[5] + Fx[6]), the moments were obtained.
    For each line the forces in the 3 directions are obtained (eg Fx [1], Fy [1], Fz [1]).
    Defined C = (x_c, y_c, z_c) the vector connecting the point (0,0,0) with the anchor fairlead, we obtain the C_rot arm to be multiplied by the forces to obtain the moments:
    C_rot = [AngleMat ^ T] * {x_c; y_c; z_c}
    in fact, the arm used to obtain the moments from the forces varies with the variation of the position of the platform.
    Finally [Mx, My, Mz] = C_rot x [Fx, Fy, Fz] (example [Mx [1], My [1], Mz [1]] = C_rot [1] x [Fx [1], Fy [1], Fz [1]])

  • Frestoring: where Fx, Fy, Mz are 0,
    Fz=mg
    Mx = M(2,6)     g * ‘PtfmYaw’ - M(1,5) * g * ‘PtfmRoll’
    My = M(3,4)* g * ‘PtfmYaw’ - M(1,5) * g * ‘PtfmPitch’

  • Fnonlinear: defined the inertia matrix I =
    [[M (4,4), M (4,5), M (4,6)];
    [M (5,4), M (5,5), M (5,6)];
    [M (6,4), M (6,5), M (6,6)]]
    (which therefore depends on the movements of the platform), and defined omega = [‘PtfmRVxi’, ‘PtfmRVyi’, ‘PtfmRVzi’],
    as regards the moments I used the formula omega x I dot omega,
    as regards the forces I used the formula omega x (omega x r_cm * mass) where r_cm * mass = [M (6,2), M (4,3), M (5,1)]

I hope this will fully clarify the formulas I use.
Best regards,
Lorenzo.

53

Dear Lorenzo,

From my brief review, your equations make sense to me. Just two comments:

  • Presumably you are using TransMat^-1 and TransMat^T^-1 as stated in your equation for M in place of the stated TransMat and TransMat^T.
  • I proposed a nonlinear way of calculating F_restoring in my post dated Sep 30, 2020 above, which is different that what you have stated.

Best regards,

54

Dear Jason,

  • Yes, starting from the formula M_cg = TransMat^TM_swlTransMat I supposed M_swl = TransMat^T^-1M_cgTransMat^-1. Is it correct?
  • I thought my change to the restoring moment was the one you suggested. I will check further

Thanks for your reply,
Lorenzo

55

Dear Lorenzo,

Yes, I agree with the equations:

M_cg = TransMat^TM_swlTransMat

and

M_swl = TransMat^T^-1M_cgTransMat^-1

Best regards,

56

Dear Jason,
Investigating hydrodynamics, I found a difference regarding the radiation forces. Specifically, I performed the same test with fast, setting the calculation of the radiation forces for the first time with the convolution integral, the second with the state-space model.
The outputs [RdtnFx, RdtnFy, RdtnFz, RdtnMx, RdtnMy, RdtnMz] of both tests I then purified of the component Ainf * acc. Plotting the results, however, I found a big difference between the two tests (fig.1 shows only Mx, but all DOFs have similar differences).
How could this difference be due?

fig1.PNG
fig1.PNG1320×550 13.8 KB

57

Dear Riccardo,

I’m not sure what you mean when you say the following. The radiation outputs from the HydroDyn module include both added mass and radiation damping.

The state-space solution will only match the convolution solution closely when the state-space model is set up to provide an excellent match (likely requiring many states). How are you setting up the state-space model?

Best regards,

58

Hello Jason,
I discovered the same problem as Riccardo Caradonna. I am analysing the floating platform Volturn US and the 15 MW reference wind turbine. I found that the radiation-damping forces change if I simulate the convolution integral or the state-space approximation in a Free-decay test with the rotor spinning (so taking the gyroscopic forces) but with CompAero=false, CompServo=false.
I am using the same hydrodynamic data as you can find on Github.
The problem is dramatic in surge and yaw, but in roll, pitch, and heave the results are similar.
What could be the problem?
thank you in advance,
Massimo Sirigu

yaw
yaw659×511 18 KB

surge
surge659×511 38 KB

59

Dear @Massimo.Sirigu,

Presumably you are running the OpenFAST model of the IEA Wind 15-MW reference wind turbine atop the UMaine semisubmersible (Volturn US) found here: IEA-15-240-RWT/OpenFAST/IEA-15-240-RWT-UMaineSemi at master · IEAWindTask37/IEA-15-240-RWT · GitHub; is that correct?

I’m not familiar with the state-space model provided in this repository, but it sounds like you want this state-space model to have a better fit to the convolution-based approach. That means to me that you’d have to rerun the SS_Fitting pre-processor (SS_Fitting | Wind Research | NREL) and target a better fit than the provided state-space model has.

Best regards,

60

Hi Jason,
Yes, it is exactly that.
So should I trust more the convolution integral than the SS?
because I am working on a model of OFWT in simscape multibody using WecSim (a wave energy converter library developed by NREL and Sandia, that I rearranged for floating platforms) and the convolution integral here gets the same results as the SS, (similar in yaw and surge) but not the same as the convolution integral in FAST. But another problem is that all the radiation damping forces in wecsim are way smaller than the convolution integral in FAST, so maybe I forgot to dimensionalise the radiation damping.
I have to investigate more, thank you for your support anyway.