Hydrodynamic Drag loads formulation

Lorenzo.Cottura · October 30, 2020, 5:27pm

Dear Jason,

In this way I can create a 2D grid of the speed of the waves on the structure: MαNβVxi, etc. are the velocities of the waves on the moving structure, and therefore depend on the displacement of the platform.
If I use another substructure or a simulation with different wind, I will have different displacements of the system and therefore also the speeds of the waves impacting the structure should be different.
From what I understand, MαNβVxi, etc. it is integral with the structure, moving with it and therefore changing according to its displacement. This is why it is different from U-wind V-wind W-wind output from Turbsim, which instead are placed on a fixed grid, not dependent on the movement of the system.
My need is to have a 3D wave speed grid that is independent of the platform position: depending on the platform positions I will go to the 3D grid to see the speeds at the knots.
Do you think it is possible, at present, to create the 3D grid?
Thanks for your help,
Lorenzo.

Jason.Jonkman · October 30, 2020, 5:46pm

Dear Lorenzo,

Actually, MαNβVxi, etc. are the fluid velocities (not including structural or relative velocity) at the undisplaced position of the structure.

Best regards,

Lorenzo.Cottura · October 30, 2020, 6:07pm

Dear Jason,
Ok, I get it, so the 2D grid (y-z) is the one found in the coordinates (0, y, z). But if for example during the simulation the platform is at a coordinate x = 10m, I have the grid at x = 0m, how can I know the speeds that act on the structure at that instant (when the structure is at 10m)?
Using the 2D grid for any x coordinate value means assuming that, at a given instant, at a node n = (y_n, z_n), the velocity of the wave particle Vx (x, y_n, z_n) = Vx (0, y_n, z_n) whatever x is.
Thanks for your cooperation,
Lorenzo.

Jason.Jonkman · October 30, 2020, 6:56pm

Dear Lorenzo,

The wave kinematics in the strip-theory solution of HydroDyn are currently computed at the undisplaced position of the structure. This follows the long wave-length approximation where we assume the structural motions are small relative to the wave length.

I mentioned our ARPA-E ATLANTIS WEIS project in a previous post. Included in that project is an effort to add a grid of points where wave kinematics will be solved, and use that together with the displaced position of the substructure to calculate the wave kinematics at the displaced position via interpolation. This will be a more accurate approach for large motion of the structure relative to the wave length.

Best regards,

George.Chan · July 9, 2021, 4:17am

Dear Jason

I am George Chen, a PhD student major in wind energy in UNSW. Thank you very much for your consistent help！

I want to ask a conceptual question about the Hydrodyn Module:

Does this module take into account the hydrodynamic load above the mean sea level (MSL) due to the variation of wave profile?

I tried to simulate the hydrodynamic load on a OC3 monopile and extract its hydrodynamic drag and inertia load from “MiNiFDxi”,“MiNiFIxi” after definition of the “node location” in MEMBER OUTPUT LIST. However, I found the output items seems to quickly approach zero around the MSL. Since the wave height Hs set in the Module is around 7m in my simulation, I expect there will be hydrodynamic load acting upon the monopile shaft above MSL.

The relevant settings in Hydrodyn is shown as follow. Many thanks for your patience in advance in reading it!

---------------------- ENVIRONMENTAL CONDITIONS --------------------------------
1027 WtrDens - Water density (kg/m^3)
20 WtrDpth - Water depth (meters)
0 MSL2SWL - Offset between still-water level and mean sea level (meters)
---------------------- WAVES ---------------------------------------------------
2 WaveMod - Incident wave kinematics model
0 WaveStMod - Model for stretching incident wave kinematics to instantaneous free surface
1000 WaveTMax - Analysis time for incident wave calculations (sec)
0.25 WaveDT - Time step for incident wave calculations (sec)
6.50 WaveHs - Significant wave height of incident waves (meters) [used only when WaveMod=1, 2, or 3]
12.50 WaveTp - Peak-spectral period of incident waves (sec) [used only when WaveMod=1 or 2]
3.3 WavePkShp [gamma] - Peak-shape parameter of incident wave spectrum (-) or DEFAULT (string)
0.1257 WvLowCOff - Low cut-off frequency or lower frequency limit of the wave spectrum beyond which the wave spectrum
3.1416 WvLowCOff - High cut-off frequency or upper frequency limit of the wave spectrum beyond which the wave
0 WaveDir - Incident wave propagation heading direction (degrees) [unused when WaveMod=0 or 6]
0 WaveDirMod - Directional spreading function {0: none, 1: COS2S} (-) [only used when WaveMod=2,3, or 4]
1 WaveDirSpread - Wave direction spreading coefficient ( > 0 ) (-)
1 WaveNDir - Number of wave directions (-)
0 WaveDirRange - Range of wave directions (full range: WaveDir +/- 1/2*WaveDirRange) (degrees)
1234567 WaveSeed(1) - First random seed of incident waves [-2147483648 to 2147483647] (-) [unused when WaveMod=0, 5, or 6]
7654321 WaveSeed(2) - Second random seed of incident waves [-2147483648 to 2147483647] (-) [unused when WaveMod=0, 5, or 6]
False WaveNDAmp - Flag for normally distributed amplitudes (flag) [only used when WaveMod=2, 3, or 4]
“” WvKinFile - Root name of externally generated wave data file(s) (quoted string) [used only when WaveMod=5 or 6]
1 NWaveElev - Number of points where the incident wave elevations can be computed (-) [maximum of 9 output locations]
0 WaveElevxi - List of xi-coordinates for points where the incident wave elevations can be output (meters)
0 WaveElevyi - List of yi-coordinates for points where the incident wave elevations can be output (meters)

---------------------- MEMBER JOINTS -------------------------------------------
2 NJoints - Number of joints (-) [must be exactly 0 or at least 2]
JointID Jointxi Jointyi Jointzi JointAxID JointOvrlp [JointOvrlp= 0: do nothing at joint, 1: eliminate overlaps by calculating super member]
(-) (m) (m) (m) (-) (switch)
1 0.00000 0.00000 -20.00010 1 0
2 0.00000 0.00000 10.00000 1 0

-------------------- MEMBERS -------------------------------------------------
1 NMembers - Number of members (-)
MemberID MJointID1 MJointID2 MPropSetID1 MPropSetID2 MDivSize MCoefMod PropPot
(-) (-) (-) (-) (-) (m) (switch) (flag)
1 1 2 1 1 0.5000 1 FALSE

---------------------- MEMBER OUTPUT LIST --------------------------------------
1 NMOutputs - Number of member outputs (-) [must be < 10]
MemberID NOutLoc NodeLocs
(-) (-) (-)
1 9 0.133 0.233 0.333 0.433 0.533 0.666 0.677 0.688 0.699

---------------------- OUTPUT CHANNELS -----------------------------------------
“M1N1FDxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N2FDxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N3FDxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N4FDxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N5FDxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N6FDxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N7FDxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N8FDxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N9FDxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2

“M1N1FIxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N2FIxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N3FIxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N4FIxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N5FIxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N6FIxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N7FIxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N8FIxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
“M1N9FIxi” - Longitudinal, lateral, and vertical drag forces at output member 2 and output node location 2
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

As shown in the figures, 9 outputs of drag and intertia uniformly distributed load (UDL) seems to vanish right after reaching the MSL.

I want to ask whether is there anything wrong in my settings, or the HydroDyn Module ignores the hydrodynamic loads above MSL?

Thank you for your precious time, and I am eager to wait for your response!

Yours Sincerely
George

Jason.Jonkman · July 9, 2021, 1:42pm

Hi George Chen,

In the current version of HydroDyn, wave stretching theory for extending wave kinematics above the still water level (SWL) is not included. Once implemented, wave stretching would be enabled by the user through the WaveStMod switch that exists after the WaveMod switch.* With WaveStMod = 0, wave kinematics are always computed up to the SWL regardless of whether the wave elevation is at a crest or trough.

This topic has been discussed on our forum before, e.g., see: model Validation - #23 by Jason.Jonkman.

Footnote
*NREL did have a project to add wave stretching into OpenFAST, but the project ran out of funding before a successful completion, and the solution was not numerically stable for structural flexible substructures; and we have not yet been funded to restart and finish the project.

Best regards,

George.Chan · July 9, 2021, 3:09pm

Dear Jason,

Thank you very much for your instructive reply!

Yours Faithfully
George

TianYuan.Wang · May 9, 2023, 2:59pm

Dear Jason.Jonkman

For IEA-15-240-RWT-UmaineSemi (Ptfmsway is particularly large - #3 by Cuizhi.Zhu), official tutorial calculates quadratic drag using AddBQuad ( IEA-15-240-RWT/IEA-15-240-RWT-UMaineSemi_HydroDyn.dat at master · IEAWindTask37/IEA-15-240-RWT · GitHub). Could you please provide some reference on how to calculate AddBQuad (0.5CdRhoA), especially for the area A in roll, pitch and yaw directions.
In addition, does the AddBQuad have a reference point? I found AddBQuad directly multiplies with u|u| at PRP in HydroDyn.f90, but 0.5CdRhoA seems to have no pointing. Does the AddBQuad also apply to platform COG? Or how to calculate the AddBQuad at COG?

q(indxStart:indxEnd)   = reshape((/real(u%WAMITMesh%TranslationDisp(:,iBody),ReKi),rotdisp(:)/),(/6/))
qdot(indxStart:indxEnd)   = reshape((/u%WAMITMesh%TranslationVel(:,iBody),u%WAMITMesh%RotationVel(:,iBody)/),(/6/))
qdotsq(indxStart:indxEnd)   = abs(qdot(indxStart:indxEnd))*qdot(indxStart:indxEnd)
qdotdot(indxStart:indxEnd)   = reshape((/u%WAMITMesh%TranslationAcc(:,iBody),u%WAMITMesh%RotationAcc(:,iBody)/),(/6/))
m%F_PtfmAdd = p%AddF0(:,1) - matmul(p%AddCLin(:,:,1), q) - matmul(p%AddBLin(:,:,1), qdot) - matmul(p%AddBQuad(:,:,1), qdotsq)

Best regards,

Tianyuan

Jason.Jonkman · May 9, 2023, 4:51pm

Dear @TianYuan.Wang,

More often than not, the additional platform stiffness and damping are calibrated to known solutions rather than “calculated”; if you want to use equations involving 0.5Cdrho*Area I would expect that it would be easier to use the strip theory (Morison equation) formulation directly.

That said, if you know that draft force is applied at a given point that is offset from the platform reference point, you could calculate the equivalent AddBQuad matrix by a suitable transformation of the kinematics–that is express the motion at the given point based on rigid-body kinematics of the platform reference point (forces transfer directly with moments added due to the moment arm).

Regarding transferring matrices expressed relative to the platform reference point and another point (such as CoG), please see the following forum topic: Time domain analysis by AQWA.

Best regards,

TianYuan.Wang · May 10, 2023, 5:16am

Dear Jason.Jonkman
I would like to use strip theory, but the official tutorial calculates quadratic drag via AddBQuad and they give no information on the drag coefficient Cd ( A request for OpenFOAM case files of VolturnUS-S platform · Issue #136 · IEAWindTask37/IEA-15-240-RWT (github.com)).
Theoretically, we can get quadratic drag and velocity at COG based on the the quadratic drag, displacement and velocity at PRP. By solving the equation system of F_{quaddragCOG} = -0.5CdrhoAv_Cog|v_Cog| at each time step, we can ontain AddBQuad at COG. Unfortunately, the resulting AddBQuad at COG obtained at each time step is different.

Q11 = zeros(length(time),1);
Q33 = Q11;
Q44 = Q11;
Q66 = Q11;
Q15 = Q11;
for j = 1:length(time)
    Q33(j) = -COGQuadDrag(j,3)/(abs(COGVel(j,3))*COGVel(j,3));
    Q66(j) = -COGQuadDrag(j,6)/(abs(COGVel(j,6))*COGVel(j,6));
    AbsVV  =  [abs(COGVel(j,1))*COGVel(j,1)  abs(COGVel(j,5))*COGVel(j,5) 0;
               abs(COGVel(j,2))*COGVel(j,2) -abs(COGVel(j,4))*COGVel(j,4) 0;
               0 -abs(COGVel(j,2))*COGVel(j,2) abs(COGVel(j,4))*COGVel(j,4)];
    F      = -COGQuadDrag(j,[1 2 4])';
    x      =  linsolve(AbsVV,F);
    Q11(j) =  x(1);
    Q15(j) =  x(2);
    Q44(j) =  x(3);
end

For example, the Q11 results at the first ten time step are
1981278.42417695
986418.515674916
5418852.95329902
3454787.60711759
3184752.23508601
3094935.46173591
3061931.89612453
3054019.92542479
3060407.78292411
3074873.62856235

Q33 results at the first ten time step are
8728844.43833093
10662998.1085941
11317353.3910709
13239035.1559764
17160801.0458012
27435741.5191107
80450679.6279351
-1667776888.12429
-11280066.9452040

So I’m confused about the formula for AddBQuad. Why Q15 is non-zero in official tutorials? Could you give me some advice
In addition, could you please give me some guidance on how to calculate quadratic drag coefficients in Table 4-5 of report NREL/TP-5000-60601?

PS: I assume that AddBQuad is non-zero only on the diagonal, but the results are equally unsatisfactory

for j = 1:length(time)
    Q11(j) = -COGQuadDrag(j,1)/(abs(COGVel(j,1))*COGVel(j,1));
    Q22(j) = -COGQuadDrag(j,2)/(abs(COGVel(j,2))*COGVel(j,2));
    Q33(j) = -COGQuadDrag(j,3)/(abs(COGVel(j,3))*COGVel(j,3));
    Q44(j) = -COGQuadDrag(j,4)/(abs(COGVel(j,4))*COGVel(j,4));
    Q55(j) = -COGQuadDrag(j,5)/(abs(COGVel(j,5))*COGVel(j,5));
    Q66(j) = -COGQuadDrag(j,6)/(abs(COGVel(j,6))*COGVel(j,6));
end

Jason.Jonkman · May 10, 2023, 4:04pm

Dear @TianYuan.Wang,

I’m not really sure what you are referring to by “tutorial”.

The AddBQuad values documented for the OC4-DeepCwind semisubmersible in Table 4-5 of NREL/TP-5000-60601 were calculated by choosing coefficients that result in free-decay behavior that best fit the free-decay behavior calculated by a model with a more complicated distributed viscous drag formulation. I believe the same process was used to derive AddBQuad for the IEA Wind 15-MW RWT.

But again, if you already know the viscous drag coefficients that can be used in a strip-theory (Morison’s equation) solution distributed across the substructure, I would use that rather than use AddBQuad.

Best regards,

TianYuan.Wang · May 11, 2023, 12:44am

Dear Jason.Jonkman
I see, thanks very much

Best regards,

tianyuan

Topic		Replies	Views
help for Hydrodyn. Wave Energy	2	353	December 21, 2021
Viscous damping on the floating spar buoy General	3	634	October 27, 2021
Hydrodynamic implementation Wave Energy	63	24256	October 12, 2022
axial drag Computer-Aided Engineering Software Tools	1	983	June 2, 2020
internal axial forces and bending moments acting the support Structural Analysis	11	5593	May 29, 2017

Hydrodynamic Drag loads formulation

Related topics