OpenFAST / openfast

Main repository for the NREL-supported OpenFAST whole-turbine and FAST.Farm wind farm simulation codes.
http://openfast.readthedocs.io
Apache License 2.0
671 stars 453 forks source link

SubDyn: Floating system issues #1956

Open RBergua opened 8 months ago

RBergua commented 8 months ago

Bug description I'm simulating a floating system in SubDyn and the behavior is not the expected one. The floater looks as follows: image

The problem seems to come from the elements that connect the vertical stem to each leg. I have tried several elements available in SubDyn (e.g., cables, rigid links, and slender beams) and all of them return different results and the behavior is not the expected one.

First, I simulated the platform in ElastoDyn only considering a lumped mass and inertia at the proper center of mass for the floater. I did the same exercise in SubDyn with excellent agreement. And this is also aligned with the output that I obtained analytically. For reference, these 3 approaches take the mass, center of mass, and moments of inertia from the SubDyn summary file for the "detailed" design.

The problem is that when I try to simulate with the "detailed" floater as rigid (i.e., SttcSolve = False, Nmodes = 0 in SubDyn), the output is all over the place. Below you can see a plot with the most interesting outputs: i. Reference: output from ElastoDyn (lumped mass and inertia), SubDyn (lumped mass and inertia), and analytical solution. ii. SubDyn detailed with pretensioned cables (the cables have an initial pretension T0). iii. SubDyn detailed with cables (the cables don't have initial pretension). iv. SubDyn detailed with slender beams (the cables are removed and slender beams with equivalent axial stiffness and mass per unit of length are used). v. SubDyn detailed with rigid link elements (the cables are removed and rigid links with the same mass per unit of length are used). vi. SubDyn detailed without tilted elements connecting the vertical stem and the legs.

image

As it can be observed, I can only recover the expected behavior if I get rid off all the elements connecting the vertical stem with the legs (24 elements in total).

Side note: in order to reduce possible differences coming from the glue code, I simply use a hydrostatic stiffness matrix (heave, roll, and pitch) in HydroDyn. There are no Morison elements. All the tests use the same hydrostatic approach.

Interestingly, if I look at the Guyan modes for these different SubDyn models, I can observe large differences. For this floating system, I would expect the Guyan modes being the six rigid body modes (0 Hz). However, this is not true for most SubDyn models with elements like rigid links, pretensioned cables... See below for reference: i. Reference: image

ii. SubDyn detailed with pretensioned cables (note that the last 3 Guyan modes are at 0.75 Hz, 0.75 Hz, and 0.93 Hz): image

iii. SubDyn detailed with cables (no initial pretension): image

iv. SubDyn detailed with slender beams: image

v. SubDyn detailed with rigid links (note that the last 3 Guyan modes are at 0.58 Hz, 0.58 Hz, and 0.70 Hz): image

vi. SubDyn detailed without tilted elements connecting the vertical stem and the legs. image

Another interesting fact is that despite using equivalent elements in terms of mass, the summary file from SubDyn returns different center of masses for each of these runs: i. Reference: CMz = -9.92 m. ii. SubDyn detailed with pretensioned cables: CMz = -9.61 m. iii. SubDyn detailed with cables (no initial pretension): CMz = -9.91 m. iv. SubDyn detailed with slender beams: CMz = -9.92 m. v. SubDyn detailed with rigid elements: CMz = -9.68 m. vi. SubDyn detailed without tilted elements: CMz = -9.82 m.

Specially surprising is how the SubDyn model with vs without initial cable pretenstion, experiences a significant jump in the center of mass reported by the summary file (-9.61 m vs -9.91 m).

OpenFAST Version I'm using OpenFAST dev version.

jjonkman commented 8 months ago

Dear @RBergua,

I would guess the issue is related to the Guyan modes not being pure rigid-body modes as expected for this floating case. Are you getting warnings from the eigensolver when running SubDyn for cases ii-vi? Presumably these are run in single precision; if so, does switching to double precision change the result?

Perhaps this is unrelated, but can you clarify what GuyanLoadCorrection you are using?

Best regards,

RBergua commented 8 months ago

Dear @jjonkman,

All the simulations have been performed with GuyanLoadCorrection = True.

I'm not receiving any warning from SubDyn/OpenFAST. See below for reference: image

I have also rerun the different models using an OpenFAST solver compiled in double precision. The results that I obtain in time domain are the same. Regarding the Guyan modes, they are still not showing the expected values. The eigenfrequencies obtained are very similar to the ones obtained when compiling in single precision.

Below you can see the Guyan modes for the detailed SubDyn model with pretensioned cables: image

Note that the Guyan modes when using cables without initial pretension (iii) or the slender beams (iv), are very similar to the ones obtained for the detailed model without tilted elements connecting the vertical stem and the legs (vi). However, only the output from the model vi is equivalent to the reference output (model i).

jjonkman commented 8 months ago

Dear @RBergua,

What do the ill-conditioned Guyan modes look like? (Can you share the JSON files?)

I assume you'd get the correct answer from SubDyn if all members were beams and all joints were cantilvered; is that correct? If so, which features (universal joints, cables, rigid links) result in ill-conditioned Guyan modes?

Best regards,

RBergua commented 8 months ago

Dear @jjonkman,

The wrong Guyan modes look like rotational rigid body modes, but they have some inherent elastic contributions. For example, I have included one dummy rigid rigid link between the top of the vertical stem and the tip of the leg 1. The rigid link connects to a duplicated node of the leg tip location and it's just used as a reference to know the motion that would be described by a rigid system. image

Attached you can also find the JSON file: USFLOWT_6legs_still_water.SD.CBmodes.json

Guyan mode 6 illustrates pretty well the issue.

Moreover, I have performed some additional tests and observed that the issue is related to the universal joints used at the base of the legs. That explains why the model vi (no tilted elements connecting the vertical stem and the legs and cantilever joints at the base of the legs) returned the expected values.

I may be able to overcome this issue by using spring elements instead of universal joints. For example, I can model the kinematic joints by means of the penalty method.

Interestingly, once the universal joints are removed, if I rerun the model with cables without initial pretension, I can reproduce the expected behavior ("Reference"). However, if I use cables with initial pretension, the platform pitch frequency, as well as the equilibrium position, are different. This is specially surprising as all these tests are performed for a rigid floater (SttcSolve = False and Nmodes = 0 in SubDyn). Important: using an initial pretension for the cables (despite the universal joints being removed), also results in wrong Guyan modes. image

Summarizing, there seems to be 2 issues:

  1. The use of universal joints (workaround: spring elements to model the kinematic joint).
  2. The use of cables with initial pretension.
jjonkman commented 8 months ago

Dear @RBergua,

Thanks for isolating these two issues further.

For the universal joints, I'm curious if this issue exists anytime a universal joint is used or if there is something unique about your complicated model? The simplest model I can think of would be a triangle made up of three beams, two cantilvered joints, and one universal joint. Does the same issue exist for that model? And is the issue isolated to universal joints, or also to pin and ball joints?

For the pretensioned cable elements, I'm a bit surprised by this because I'm aware of several SubDyn models that have been developed that make use of pretensioned cable elements. So, I'm curious why some SubDyn models work and others don't.

Best regards,

RBergua commented 8 months ago

Dear @jjonkman,

It's hard to tell if the simple model with one universal joint would return the expected results. Note that if the Guyan modes are wrong, then the output in time domain is likely wrong. However, as the initial model iii shows (SubDyn detailed with cables that have no initial pretension and universal joints), having the proper Guyan modes does not ensure that the behavior in time domain will be the expected one. The proper Guyan modes with wrong time domain behavior were also observed for the model iv (slender beams approach and universal joints).

Based on tests performed for the USFLOWT, if I define pin joints, the results are very close to the expected ones but not perfect. See below for reference: image image

Maybe the assumption that we make about SttcSolve = Flase and Nmodes = 0 resulting in a rigid system is somehow not true?

Interestingly, if I define all base legs as cantilever joints except one leg defined as kinematic joint (e.g., universal joint), the outputs are the expected ones. However, if I define two universal joints, I already observe the solver behaving different and having difficulties. image

This would indicate that the element itself is not to the problem. But rather the multiple uses of it and how this is assembled on the global matrices (?)

Regarding the use of pretensioned cable elements, I also observed problems during the TetraSpar project (https://wes.copernicus.org/preprints/wes-2023-103/). That numerical model uses slender beams rather than pretensioned cable elements. For reference, using the slender beams returns the same behavior as cable elements without initial pretension. I observed that providing different initial pretensions resulted in different final tensions for the system in equilibrium conditions. I'm not sure if there was another project using pretensioned cable elements and they checked the outputs against analytical values to ensure that the behavior was the expected one... image

Note that the axial stiffness for the keel lines in the OC6 phase IV project was very large (3.618E8 N/m). Including or not including an initial pretension has a very minor impact in the equilibrium position of the keel. But it has a significant impact in the loading (as shown above).

For reference, the TetraSpar with pretensioned cable elements (e.g., T0 = 4.25E6 N) returned Guyan modes close to 0 Hz, but the behavior in time domain was not the expected one. image

jjonkman commented 8 months ago

Thanks, Roger. Indeed, both issues sound more complicated.

I'm hoping that @ebranlard can comment, or perhaps @deslaughter, as he is currently reviewing SubDyn in preparation for introducing new functionality such as support for larger floater yaw motion.

It will probably be important to share a model where the issues is clear to aid in debugging.

Best regards

luwang00 commented 1 month ago

Hi @jjonkman, @RBergua and I have investigated this a little bit more. More specifically, we are looking into why the Guyan modes for a floating structure including pre-tensioned cables can deviate from the expected rigid-body modes. I think we have a plausible explanation.

The non-zero Guyan mode stiffness and natural frequency for the rotational degrees of freedom (Modes 3, 4, and 5) came from the unbalanced geometric stiffness of the cable pretension. These are the stiffness terms associated with u_x and u_y (in red boxes) in the cable element stiffness matrix below. Note that these terms are simply equal to T0/Le, where Le is the element length. They represent the stiffness normal to the cable due to the rotation of T0.

image

Consider a very simple example where we have a vertical beam. The top node is the TP node and the bottom node is a free node. Parallel to the beam and connecting the same two nodes is a pre-tensioned cable element. See sketch below. When we rotate about the TP Node, T0 is also rotated. The horizontal component of T0 on the lower node after rotation creates a stiffness about the TP node opposing the rotation. A quick experiment with SubDyn using this setup shows that the spurious rotational stiffness for the Guyan modes (roll and pitch) is indeed very close to T0·Le, where Le is length of the element (the moment arm) as expected. (The stiffness is exactly equal to T0·Le if the beam is replaced with a rigid link.)

image

The reason for this spurious moment is because we don't have a structure that is in equilibrium. In the absence of gravity and other forces, the cable will simply shorten until the preload is zero. Of course, in this example, we also have a beam element coinciding with the pre-tensioned cable. The axial stiffness of the beam will help keep the cable in tension. This also means the beam should have a pre-compression that is missing in the current SubDyn implementation. The pre-compression in the beam is equal and opposite the cable pretension. This additional force is shown in the sketch below in orange. Note that just as with the cable pretension, the beam preload will also lead to geometric stiffness associated with the displacement of the lower node. Since the force is equal and opposite, the contribution from the beam preload will cancel the contribution from the cable pretension, resulting in zero stiffness for the Guyan mode, which will return to purely rigid-body modes as expected.

image

In summary, the issue we are observing is due to only including the geometric stiffness coming from the cable pretension, but not the geometric stiffness associated with the preload in other elements, which inevitable arises for a structure in equilibrium with pre-tensioned cables. This leads to an imbalance in Guyan mode stiffness, which should be all zeros for a floating structure.

One possible solution is to solve for the preloads in the rest of the structure due to the presence of any pre-tensioned cable elements assuming equilibrium and include the geometric stiffness coming from these preloads in the overall stiffness matrix. However, this is not always possible with arbitrary structure and cable pretension input. It might be impossible for the structure to be in equilibrium. The simplest example would be the setup above with the beam removed, leaving only a pre-tensioned cable with a free end. Alternatively, if the beam is replaced with a rigid link, it is also not easy to solve for the preload on the rigid link and incorporate its effects.

A short term solution would be to disable pretension for cables or neglect geometric stiffness. Neither solution is ideal, but at least it should prevent grossly unphysical results.

jjonkman commented 1 month ago

Thanks, @luwang00, for digging into this and explaining the issue, which makes sense to me.

With the presence of rigid links, it sounds like this will be a hard problem to solve in general (or perhaps we need to prevent SubDyn models from having both rigid link and pretensioned cable elements). I don't believe the cable-alone issue is really a problem because we don't allow SubDyn to have internal rigid-body modes anyway.

I'm curious how much error there would be in the short-term solution if we neglect the geometric stiffness associated with the pretension, given that the cable stiffness itself is still present (the EA/L0*1 terms in the stiffness matrix)? Perhaps it would be worth @RBergua checking that in his model to see if the results are plausible.

Thanks,

RBergua commented 1 month ago

I have been performing several tests for the USFLOWT design considering only one leg. All tests used the SubDyn stand-alone version. I apply one time-serie of vertical force to study a free-decay over the leg. Like that we can look at the deflection at the tip of the leg (where the force is applied) and also the dynamic behavior after the leg is suddenly unloaded.

First, I considered fixed-bottom conditions. See the representation of the system below: image

The results for this fixed-bottom condition are as follows: image

From a dynamic point of view, the cable without pretension and the cable with pretension and geometric stiffness = OFF are very similar (as expected). The only noticeable difference is the mean value (the pretensioned cables result in a slightly different equilibrium position since the resultant force in the vertical direction is different than 0).

Comparing the approaches with vs without geometric stiffness, it can be observed that when the geometric stiffness is included, the system is slightly stiffer (as expected). This results in slightly smaller displacements at the tip of the leg and slightly higher frequency for the leg bending motion.

Now we can move to floating conditions. In this case, the stem (vertical beam) has been made very stiff to avoid any deflections and try to make comparable results to the fixed-bottom condition. See the representation of the system below: image

The results for this floating condition are as follows: image

The outputs for this floating condition are pretty well aligned with the ones observed for the fixed-bottom condition. Although not ideal, disabling the geometric stiffness for the cables allows to use the pretension in floating systems within OpenFAST.

As a reminder: the problem for a user is that if it uses pretension for the cables, the geometric stiffness of the cables will result in wrong Guyan modes (the 3 rotations are not rigid body modes). This have terrible implications for the motion of the system (wrong global equilibrium position and wrong system dynamics). The plot below compares the response with and without cable pretension for the USFLOWT system (gravity-only condition and still water): image

jjonkman commented 1 month ago

Thanks, @RBergua, for the great demonstration of the issue and for outlining its impact on your USFLOWT model. Going forward in the near time, is your plan to use the model with pretension and without geometric stiffness?

At some point, we'll need to correct SubDyn so that the geometric stiffness associated with pretension is consistently accounted for across the formulation. I'm just not sure how much effort that would be and what project could fund that fix yet.

Best regards,

RBergua commented 1 month ago

Yes, @jjonkman. Our plan is to use the model with pretension and without geometric stiffness. Till now, we were using pretension = 0 N to effectively get rid off the geometric stiffness and have the proper Guyan modes (i.e., proper global dynamics of the floating system). By means of an analytical post-processing, we were including the effect of the pretension in the SubDyn output forces (e.g., cables, legs and stem) just as an offset based on the T0 value.

By removing the geometric stiffness, now we are able to include the initial pretension in the cables and avoid the analytical correction in the post-processing that I was commenting. Moreover, our analytical correction for the pretension didn't take into account the slightly different orientations that the leg and cables experience during the simulation. So, now our results will be slightly more accurate than before and easier to obtain.