BeamDyn Coordinate Systems in OpenFAST

I have a beamdyn v1.01 input file with several key points (in “GEOMETRY PARAMETER” section) defining the blade structure from root to tip. I understand I may select any number of these key points by choosing them in the outputs section. My first key point is at z=0. I select this key point 1, as well as the rest, to be written in the *.out file. My question arises when post processing blade span loads. I have loads at the root, as well as at each blade span. Blade root loads and kp1 seem like they should be the same, but they are not. Even Fz & Mz loads are different. I’d like some help understanding these coordinate systems please. I see per the manual, that root loads are “expressed in r”, while the beamdyn spans are local, but the subtleties are lost on me.
Thanks!
Ryan

Dear Ryan,

Which version of BeamDyn in OpenFAST are you using? There have been many improvements made to BeamDyn over the past couple of years, including fixes to the sectional load outputs. When enabling BeamDyn, I would recommend using the newest master branch of OpenFAST.

Best regards,

Good morning, Jason -
We are using a development version from 20th Dec 2018. “OpenFAST-v2.0.0-25-g50c88c2f-dirty”

Thanks,
Ryan

Dear Ryan,

OK, can you clarify what exactly the problem is? Can you share your BeamDyn outputs and input file settings?

Best regards,

Jason-

My confusion is discerning the difference between the root loads and the blade span 1 loads. I’m attaching a short *.out time history, illustrating the progression of the first 4 span edge moments, along with the input and summary file for just those 4 nodes.

On one hand, I think blade node 1 and root loads are both at 0% span and so would expect the loads to be identical (less any twist in xy frame)

On the other hand, the moments are so different that loads appear to be resolved at different locations in space.

I’m trying to learn how to use blade root and blade node 1:n to map out blade loads accurately.

I hope this makes sense.
Thanks for your help on this puzzle.
Ryan
beamdyn_channels.xlsx (269 KB)

Dear Ryan,

What type of quadrature are you using? (1: Gauss; 2: Trapezoidal)

What do you have OutNd set to in your BeamDyn primary input file?

Best regards,

Jason -

I have quadrature set to 2: trapezoidal.
I also have NNodeOuts=42, and then get 1,2,3…42 for Outnds. My outfile records all of these channels reasonably; again the only apparent discontinuity is at the Root and node 1.

Aside, the individual blade file has 42 stations as well; matrices for 0:1 that correspond with key points.

---------------------- SIMULATION CONTROL --------------------------------------
False Echo - Echo input data to “.ech” (flag)
False QuasiStaticInit - Use quasistatic pre-conditioning with centripetal accelerations in initialization (flag) [dynamic solve only]
0.0 rhoinf - Numerical Damping Parameter for Generalized-alpha integrator
2 quadrature - 1: Gauss; 2: Trapezoidal
“DEFAULT” refine - Refinement factor for quadrature 2. DEFAULT = 1
“DEFAULT” n_fact - Factorization frequency: The Jacobian is computed every n_fact steps in N-R iteration. DEFAULT = 5
0.0005 DTBeam - Time step size
“DEFAULT” load_retries - Tolerance for stopping criterion
10000 NRMax - Max number of iterations in Newton-Ralphson algorithm
“DEFAULT” stop_tol - Tolerance for stopping criterion
“DEFAULT” tngt_stf_fd - finite difference for tangent stiffness flag
“DEFAULT” tngt_stf_comp - compare tangent stiffness using finite difference flag
“DEFAULT” tngt_stf_pert - perturbation size for finite differenced tangent stiffness
“DEFAULT” tngt_stf_difftol - tolerance for difference between analytical and fd tangent stiffness
True RotStates - Orient states in the rotating frame during linearization? (flag) [used only when line

(IP material)

---------------------- MESH PARAMETER ------------------------------------------
5 order_elem - Order of interpolation (basis) function (-)
---------------------- MATERIAL PARAMETER --------------------------------------
BeamDyn_B1.dat BldFile - Name of file containing properties for blade
---------------------- PITCH ACTUATOR PARAMETERS -------------------------------
False UsePitchAct - Whether a pitch actuator should be used (flag)
200 PitchJ - Pitch actuator inertia (kg-m^2) [used only when UsePitchAct is true]
2.0E+7 PitchK - Pitch actuator stiffness (kg-m^2/s^2) [used only when UsePitchAct is true]
5.0E+5 PitchC - Pitch actuator damping (kg-m^2/s) [used only when UsePitchAct is true]
---------------------- OUTPUTS -------------------------------------------------
True SumPrint - Print summary data to “.sum” (flag)
“ES10.3E2” OutFmt - Format used for text tabular output, excluding the time channel. Resulting field should be 10 characters. (quoted string)
42 NNodeOuts - Number of nodes to output to file [0 - 9] (-)
1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41,42 OutNd - Nodes whose values will be output (-)

Dear Ryan,

Have you changed the BeamDyn source code in some way? In NREL’s version of BeamDyn, NNodeOuts can be at most 9, but you’ve set it to 42. (?)

I would expect at the root, where the blade is cantilevered to the hub, that the only difference between the local (l) and root (r) coordinate systems would be the structural twist at the root (initial_twist). Thus, I would expect that Fzr = Fzl, Mzr = Mzl, Fxr and Fyr would be related to Fxl and Fyl through SIN(initial_twst) and COS(initial_twist), and Mxr and Myr would be related to Mxl and Myl through SIN(initial_twst) and COS(initial_twist).

Best regards,

Hi Jason -

We have modified the Beamdyn source code to accept 42 nodes. Your expectations about root loads and span 1 loads (at root) being the same, save twist angle, are helpful. I was expecting the same, but finding otherwise. I will run some tests with the original 9 node beamdyn code and hopefully see that blade root loads match span1 (kp1 at z=0m) loads appropriately. I will share my results.
Thank you,
Ryan

Hi Jason -

After a variety of experiments, I thought it perhaps easiest to use the 5MW NREL beamdyn model. I’m attaching a zip file with the an out file, a xls with some graphs, and beamdyn input parameters (with some adjustments, migrating from F8 to OF).
-My goal is to be certain of the span address & orientation of each local blade node.
-I understand from the manuals that node one is at the root, and as such should be equivalent values (after twist rotations)
-Results from the 5MW appear to suggest that B1Root loads are at a slightly different location than B1N1 loads. Fx,Fy,Mx,My loads are very similar, but differ as I would expect for N1,N2, N3 etc…down the blade. Mz and Fz loads show greater disparity between B1Root and B1N, suggesting I’m comparing incorrectly.

These are only clues that I gather in attempt to understand if B1N1 is at the same z=0 span, or conversely, to understand at what blade span, z, these B1N(beta) loads are prescribed. Perhaps another way to visualize the question: With 9 blade nodes and the root, do you have 10 points along blade span, from root to tip, describing loads? Or just 9, with root and N1 being redundant?

Thanks for you patience and help on this,
Ryan
NRELOffshr5MW_BeamDyn_CS_test.zip (2.21 MB)

Dear Ryan,

The 1st node in BeamDyn is indeed at the blade root, so, I would generally expect equivalent results between the root and 1st node (after twist rotations). My understanding is that the root and nodal reaction loads are calculated in a different way in BeamDyn, but that the results should be equivalent (after twist rotations). In looking at your results, I do see that the axial force matches between the root and 1st node (they don’t match in your plot, because you were accidentally comparing N1Fzl with RootMzr (instead of RootFzr)). I do see some difference between the N1Mzl and RootMzr, which I wouldn’t expect, but it is also a bit odd the way the blade-pitch angles look in this simulation, with pitch angle ramping up from 0 to nearly 80 degrees in the first second. I’m not sure what your ElastoDyn, AeroDyn, and ServoDyn settings are that would trigger this pitch ramp, but perhaps the difference in calculation approach between the root and nodal reaction loads shows up because of this extreme pitch ramp? Does this moment compare better when the pitch angle time series is more “reasonable”?

FYI: I would recommend enabling the quasi-static preconditioning at initialization (QuasiStaticInit=True) to reduce the impact of start-up transients on the BeamDyn solution.

FYI: There is yet-to-be merged pull request aimed at improving the nodal reaction load output calculations in BeamDyn–see: github.com/OpenFAST/openfast/pull/163.

Best regards,

Hi Jason

Your guidance is helpful. While I’ve found differences in the values between root and node 1, I expect that they are simply negligible numerical differences.

I ran an additional test of the 5MW offshore blade. I kept wind steady, and shut off the pitch controller. It is convenient to note that the all 6 nodes output have the same twist angle of 13.308deg. Vector sum channels of the RootMxy and Node 1Mxy match identically. This suggests to me that here is some numerical noise in my rotation of the the root loads, but certainly indicate the coincidence of these virtual guages.

I’m attaching a zip file including these results. Yes, root and Node 1 time histories overlay well, except Mz. Maybe that is because beamdyn is a more accurate model for torsion? On the 2nd tab of the attached spreadsheet, I show graphs of the 6 load channels down the blade, for 6 nodes plus root, at several instances in time. Ideally, I would see root and N1 loads be quite similar. I do find it confusing that the differences between root and N1 appear quite large when comparing the differences between N1 and N2 (only 0.2m out the blade). Certainly, the percentage difference is quite low, but I find the discontinuities along the blade spans discomforting. This confusion aside, I understand that Node 1 is at z=0, and will use these node channels to deliver accurate blade loads. I’ll keep the pending beamdyn pull-request in mind for evaluation as well.

In your opinion, do the root channels deliver any unique value, particularly over node 1? I would tend to choose node1 loads for hub design, etc, for accuracy and consistency.

Thank you for your help and patience as we worked through this.
Ryan

Aside: in the OpenFAST manual (release 1.0, Jan18, 2019) is missing the graphics for figures 6.7, 6.8, and 6.9.
5MW_SteadyNoPitch.zip (2.66 MB)

Hi Ryan,

I’m not sure, as I would not expect the loads to differ between the root and the first node. I’ve placed this potential bug in BeamDyn (as well as the documentation bug) as an issue on the OpenFAST github site: github.com/OpenFAST/openfast/issues/248. I suggest discussing this potential problem further there.

Best regards,

Hi Jason -
Thank you for opening the issue in the github site. I will move there and follow along/ participate as I’m able.
Ryan