Dear Dr. Jason Jonkman
I am trying to reproduce the RAO for OC3-Hywind model given in Fig 3. of the article entitled, “Investigation of Response Amplitude operators for floating wind turbine”, by using FAST V8. I have considered all the condition same as given in the article e.g. only 6 rigid body degree of freedom are enabled. The wave model is considered 3 (white noise spectrum), wave significant height is 2 m. I use the strip theory only (disabled WAMIT) to calculate hydrodynamic forces. The drag coefficient is also Cd =0, as given in the article. I also specified the additional linear damping as given in the article and also specified the hydro-static stiffness in the input. I used fft matlab (PSD estimates using FFT) command to obtain the auto and cross spectral densities for wave elevation and the output response, and then took the ratio of the two to calculate RAO. However, there is difference in the RAO i obtained (both in frequency and amplitude). I did the simulation for 8000 sec and excluded the initial 2000 sec as given in the article.
I studied all the previous discussion regarding the RAO on this forum but still could not figure it out. I would highly appreciate if you can help me in finding where i am making a mistake. I have attached the RAO i obtained from the simulation. I also attached the snapshot of Hydrodynamic input file.
i calculated PSD taking the absolute value squared of the DFT (calculated by fft() in matlab) and then use “sgolayfilt()” for smoothing purpose. I did not know about the sampling frequency so i use 2*pi rad as sampling frequency. I will attach anything else if required. I will highly appreciate your valuable comments.
Hydrodynamic input snapshot is given here.
Dear Zahid,
I don’t really see anything wrong in your HydroDyn input file. One thing I would probably change is to set WaveNDAmp = False to ensure that the input wave spectrum is smooth. Presumably you also have the six platform DOFs enabled your ElastoDyn module and some mooring module enabled?
Regarding your results, the first thing I see is that your natural frequencies are not correct. In rad/s, I would expect to see RAO peaks in surge at 0.05, heave at 0.20, and pitch at 0.21. When you look at the time series, do you see these frequencies correctly? If so, I would expect a problem in your MATLAB script, but I’m not sure exactly what that might be.
Best regards,
Dear Dr. Jason Jonkman
I highly thankful for your help and time. I followed your suggestions and found that the frequencies corresponding to peaks matches to the natural frequencies of the structure and the difference was result of some error in the Matlab script. After correcting that i re-plotted the RAOs for surge, heave and pitch. I also set WaveNDAmp = False. Also i set MAP for mooring lines and enable only six DOFs of platform in ElastoDyn module. I have attached the plots i obtained. The frequencies corresponding to peaks matches to that of the given in the article. However, there if significant difference in the magnitude of the RAOs. I also attached the snapshot of the procedure i used to obtain PSD, in case if there is any error in it. i would appreciate your valuable comments.
I also tried to calculate PSD by using the same methodology as given in the previously attached image, but i remove the “square” from last term “abs(xdft)”. In that case the surge and pitch were approximately equal in magnitude to that of given in the reference article, but there is significant difference in the magnitude of heave RAO. I also attached those plots. I dont know much about PSD, so i will be thankful for your comments on the PSD procedure i am using, and any reference material will also be highly appreciated. If you want to take a look of the time histories of wave elevation and surge, heave and pitch, i will also attached them.
Dear Zahid,
I’ve attached a simple MATLAB script that I developed that calculates both 2-sided and 1-sided PSDs with frequencies in both rad/s and Hz directly based on fft. Does this use of this script solve your problem?
Best regards,
Jason_PSD.m.txt (2.73 KB)
Dear Dr. Jason Jonkman
I am highly thankful for your help and the script you have provided. I tried using the matlab script however, i still found significant difference in the magnitude and could not figure out the reason for it. I was thinking that there might be some mistake i am making in the hydrodynamic input. Therefore, i attached screen shot of the hydrodynamic coefficients (as this part was not visible in the last screen shot of hydrodynamic input) and will appreciate if you can take a look at it. I have also attached the time history plot of wave elevation, surge heave and pitch. I will highly appreciate any suggestion/advice. Heartfelt thanks for the help and time.
Dear Zahid,
I would not expect your strip-theory-only solution to match the potential-flow-only solution employed in the ISOPE paper because (1) you have set the axial added mass and dynamic pressure coefficients to zero and (2) you have not set the transverse added mass and dynamic pressure coefficients to match those from the potential-flow solution for surge/pitch. If you want to better match the ISOPE paper, I would (1) use the potential-flow-only solution rather than the strip-theory-only solution or (2) set your added mass and dynamic pressure coefficients in the strip-theory-only solution to better match the potential-flow-theory-only solution.
Best regards,
Dear Dr. Jason Jonkman
I am highly thankful for your response. I will try to find suitable coefficients for strip theory only solution and will also use potential-flow-only solution for comparison. I have another small query and will appreciate your response. I wondered if the we can get similar RAO if we use the JONSWAP/Pierson-Moskowitz spectrum (WaveMod = 2). Will it cause significant difference?
Dear Zahid,
The RAO is derived from the time series by dividing the cross-spectral density of the input-output by the auto-spectral density of the input in the frequency domain. This would imply that the input is irrelevant. However, I would expect numerical problems wherever the energy in the input is small i.e. at low and high frequencies when using the JONSWAP/Pierson-Moskowitz spectrum.
Best regards,
Dear Dr. Jason Jonkman
I am thankful for your help and sorry to bother you once again. Until now i was using only strip-theory in the HydroDyn, however now to better match the RAO i have to use Potential flow model only. I followed the HydroDyn user’s guide and theory manual to understand properly the methodology of using Potential flow model only. I set the “PotMod = 1” to enable the WAMIT. Following the Section 6.8 (Floating platform) of HydroDyn manual which says that “for a potential-flow-only model do not create any strip theory joints or members in the input file”. I made all the hydrodynamic coefficients zero and removed all the member joints and members etc. when i run the analysis it gives me an error and aborted. The error was:
“FAST_InitializeAll:HydroDyn_Init: Both meshes must be committed before they can be mapped.
HydroDyn_Init: Both meshes must be committed before they can be mapped.”
I searched on this forum to find any help regarding this error and fortunately found your kind response to one of the user, given as
(Link: Foundation Stiffness and damping in FAST v8 - #11 by Patrick.Quasten)
“Dear Patrick,
I was able to reproduce your error when I set up a similar model myself. It appears that there is a bug in HydroDyn preventing it from functioning properly when a model has no strip-theory nodes (NJoints = NMembers = 0). We’ll fix this in a future release of HydroDyn. In the meantime, I suggest that you set up your HydroDyn module with two joints and one member, but with all member properties and coefficients set to zero.”
I followed your suggestion and made NJoints = 2 and NMembers = 1, with all the properties equal to 0. But the error i got was
“A member can not have zero length.”
I also tried it by assigning small length e.g. 0.5 and run it but then the error i obtain was “MDivSize must be greater than zero”.
I will be highly thankful for your time and your help.
Dear Zahid,
To avoid an error, simply ensure that the two joints are not coincident e.g. set the first one at (0,0,0) and the second one anywhere but (0,0,0) e.g. set the second joint arbitrarily at (0,0,-1).
Best regards,
Dear Dr. Jason Jonkman
I am thankful for your help and by following the suggestions i finally obtain good results. I want to understand the white noise spectrum (WaveMod = 3) and have few points in mind. I would appreciate if you give some valuable comments.
- For WaveMode = 2 (JONSWAP/P-M), we need to define both Hs and Tp. But for WaveMode = 3 (White-Noise Spectrum) it is given that only Hs is required. I was studying your Phd work to understand the theoretical background of White-Noise Spectrum. In the report entitled, “Dynamics Modeling and Loads Analysis of an Offshore Floating Wind Turbine” it is given in Eq. 2-16 (page 23) that two-sided spectrum is related to 1/2 (one-sided spectrum), which is JONSWAP/P-M spectrum. So, it made me confused why the Tp is not required for WaveMod = 3, when is it related to one-sided spectrum.
Dear Zahid,
The wave peak-spectral period (Tp) defines the peak of the wave spectrum. There is no peak in a white-noise wave spectrum, so, Tp is not necessary. The equation for the one-sided power-spectral density of white-noise waves (with frequencies in rad/s) is:
WaveS1Sdd = WaveHsWaveHs/(8( MIN( Pi/WaveDT, WvHiCOff ) – WvLowCOff ) )
where the MIN() is used in case WaveHighCtOff is specified above the highest frequency resolved by the time step.
Best regards,
Dear Jason,
Can the Matlab script you uploaded be used to calculate both the auto- and the cross- power densities?
If so, would you mind to give an example on how to use it for calculating RAOs?
To save your time, I attach at the end of this post an example script I developed today together with the FAST signal file.
To my surprise, the final calculated RAO is a 2049*5000 matrix with all elements being the value of 0.0000 + 0.0000i.
Did I make mistakes in the example script?
Thanks a lot.
Best regards,
Yingyi
DriverFastPlotRAO.m.txt (806 Bytes)
Test25.zip (1.08 MB)
Dear Yingyi
I only took a brief look at your script; it looks like pxx and pxy are different sizes; once you fix that and make them the same size, you would need to use a element-by-element divide (./) rather than a matrix divide (/).
I don’t have the time to update your script for your. However, I’ve attached a MATLAB script for calculating RAOs that my NREL colleague Amy Robertson developed. This script calculates the auto- and cross-spectral densities directly using FFT (similar to my script above). Perhaps you and others will find it useful.
Best regards,
RAO.m.txt (2.44 KB)
Dear Jason,
Thanks. That really helps.
Two further questions regarding this issue:
(1) In the RAO.m script, it seems that the input signals are from multiple files, indicated by the function explanations. One of the input parameter “nens = number of ensembles to average (i.e., 4 or 6)” also confirms this point. But I don’t see where we can input the multiple signals since the inputs “time”,“U”, “Y” are all 1-D arrays. Should the multiple signals be combined into one signal before using the RAO.m script? By the way, what does “ensemble” mean?
(2) In Amy’s publication (nrel.gov/docs/fy14osti/61154.pdf), Table 3 defines the load cases run in OC4 Phase II. For the LC 3.7 (RAO estimation with wind), the wave condition is Banded white noise with PSD =1 m2/Hz for 0.05-0.25 Hz. My question is, where can we set the three parameters (PSD =1 m2/Hz for 0.05-0.25 Hz) in the HydroDyn input file?
Thanks a lot.
Best Regards,
Yingyi
Dear Yingyi Liu,
Here are my answers to your questions:
(1) Correct, you should concatenate multiple simulation outputs into a single 1-D array. You can think of “ensemble” as the mean of the outputs.
(2) Set WaveMod = 3, WvLowCOff = 0.052pi rad/s, WvHiCOff= 0.252pi rad/s, and WaveHs = SQRT( 18(0.25 - 0.05) ) m.
Best regards,
Dear Jason,
Your answers are very helpful, now the script works.
I run a simulation of LC 3.7 and find three more questions:
-
Can WaveSeed(1) and WaveSeed(2) be set as any arbitrary numbers? I’m using 4 ensembles for one calculation, and I set the corresponding WaveSeeds as
123456789 and 1011121314 (in Ensemble 1)
378281615 and 1736193745 (in Ensemble 2)
461917361 and 1793479262 (in Ensemble 3)
635438273 and 1308917163 (in Ensemble 4)
Are there any mistakes in the above WaveSeeds I set?
-
In the lower-frequency region (0.05-0.06Hz), the attached RAO of Pitch (post-process of 4600s simulation output data, no data discarded) is a little bit noisy, is it a common issue or does it mean that my calculation is not so good? Is it necessary to discard the first 1000s or more time-series data because of the transient effect (having not reached a stable status)?
- I find the FAST calculation stops inconceivably at the ending of ChkptTime which is specified in the FAST primary input file (see the screenshot below). Do you know why it happens?
Dear Yingyi Liu,
Here are my answers to your questions:
-
Yes, the WaveSeeds can be set arbitrarily.
-
Yes, white-noise-generated RAOs can be a bit noisy. Averaging over several time series will help. Yes, you should eliminate start-up transients from post-processing. But 1000 s sounds like a long time. Section 6.8.2 of the HydroDyn User’s Guide and Theory manual provides some guidance for choosing proper initial conditions to minimize start-up transients for floating wind systems: wind.nrel.gov/nwtc/docs/HydroDyn_Manual.pdf.
-
It appears that the stack size is not large enough for generating checkpoint files with the long simulation with wind and waves you are running. Do you need to produce checkpoint files? If not, I suggest that you set ChkptTime > TMax.
Best regards,
