Tilt Cosine Law

[mkrause-strath]

Hello,

I believe the way rotor effective velocities are calculated at the moment for tilt is wrong. I'm happy to be corrected though!

From FLORIS:

# Compute the rotor effective velocity adjusting for tilt
    relative_tilt = tilt_angles - ref_tilt

    rotor_effective_velocities = (
        rotor_effective_velocities
        * cosd(relative_tilt) ** (cosine_loss_exponent_tilt / 3.0)
    )

From this it can be seen that the cosine law is applied relative to the reference tilt at which the cp/ct curves are defined. I see three possible issues with this:

• Even though cp/ct curves are defined at a reference tilt the actual wind is coming from 0 degrees of course. So I think the actual correct way to calculate effective velocities would be to calculate the increase in turbine power going to 0 degrees and then the loss to the actual tilt. For example for a turbine with 5 degrees reference tilt at 15 degree tilt:
  • cos(15deg - 5deg) ^ (1.88/3.0) = 0.9905
  • 1 / (cos(0deg) ^ (1.88/3.0)) * (cos(15deg) ^ (1.88/3.0)) = 1 / 0.9976 * 0.9785 = 0.9899
• Any tilt angles below reference tilt are not supported at all in the current iteration (if modelling active ballasted floaters for instance). For instance an actual tilt of 0 degrees has a lower power than the reference tilt since the cosine law is applied - e.g. for 5 degrees reference tilt it is cos(0deg - 5deg).
• Finally, why is the cosine loss exponent 1.88 divided by 3.0. I don't understand why it isn't just the exponent?

I have read #1055 but I believe that is asking a slightly different question.

Regards,

Max

[misi9170]

Hi @mkrause-strath ,

Thanks for raising this! I believe that you're right, it would be more correct to use the ratio of the tilts, rather than the difference. I don't quite follow your math (or perhaps there are some typos?) but I believe what you are suggesting is that, letting $\theta_\text{ref}$ denote the reference tilt, and $\theta$ denote the "actual" tilt, the expression

$v_\text{rotor effective} = \cos(\theta - \theta_\text{ref})^{1.88 / 3.0}$

should instead be replaced by

$v_\text{rotor effective} = \left(\frac{cos(\theta)}{\cos(\theta_\text{ref})}\right)^{1.88 / 3.0}$.

Is that right?

When $\theta = \theta_\text{ref}$, these evaluate to the same thing; but when $\theta \neq \theta_\text{ref}$, they differ. As you point out, using the current method, any value of $\theta \neq \theta_\text{ref}$ will cause a decrease in the rotor effective wind speed, which leads to a reduction in power---even if the rotor has tilted "forward".

I made a short script to look at the differences between them:

import numpy as np
import matplotlib.pyplot as plt

tilt_ref = 5
tilt_test = 15
tilt_test_2 = -5
tilt_test_3 = 5

C_T = 1.88
base = 3.0

def cosd(x):
    return np.cos(np.deg2rad(x))

# Test 1: 15 degrees of tilt
# Current approach
correction_current = (cosd(tilt_test - tilt_ref))**(C_T / base)
print("Current approach correction: {0:.4f}".format(correction_current))

# Proposed approach
correction_proposed = (cosd(tilt_test)/cosd(tilt_ref))**(C_T / base)
print("Proposed approach correction: {0:.4f}".format(correction_proposed))

# Test 2: -5 degrees of tilt
# Current approach
correction_current_2 = (cosd(tilt_test_2 - tilt_ref))**(C_T / base)
print("Current approach correction (negative tilt): {0:.4f}".format(correction_current_2))

# Proposed approach
correction_proposed_2 = (cosd(tilt_test_2)/cosd(tilt_ref))**(C_T / base)
print("Proposed approach correction (negative tilt): {0:.4f}".format(correction_proposed_2))

# Test 3: 5 degrees of tilt (reference)
# Current approach
correction_current_3 = (cosd(tilt_test_3 - tilt_ref))**(C_T / base)
print("Current approach correction (reference tilt): {0:.4f}".format(correction_current_3))
# Proposed approach
correction_proposed_3 = (cosd(tilt_test_3)/cosd(tilt_ref))**(C_T / base)
print("Proposed approach correction (reference tilt): {0:.4f}".format(correction_proposed_3))


# Visualize differences
tilts = np.linspace(-20, 20, 100)
corrections_current = (cosd(tilts - tilt_ref))**(C_T / base)
corrections_proposed = (cosd(tilts)/cosd(tilt_ref))**(C_T / base)

fig, ax = plt.subplots(1,1)
ax.plot(tilts, corrections_current, label="Current approach", color="k", linestyle="--")
ax.plot(tilts, corrections_proposed, label="Proposed approach", color="C0", linestyle="-")
ax.plot([tilt_ref, tilt_ref], [0, 1.2], color="k", linestyle=":", label="Reference angle")
ax.grid()
ax.set_xlabel("Tilt angle [deg]")
ax.set_ylabel("Velocity correction factor [-]")
ax.legend()

ax.set_ylim([0.9, 1.1])

plt.show()

which produces

and

Current approach correction: 0.9905
Proposed approach correction: 0.9808
Current approach correction (negative tilt): 0.9905
Proposed approach correction (negative tilt): 1.0000
Current approach correction (reference tilt): 1.0000
Proposed approach correction (reference tilt): 1.0000

There are some differences between the value you reported in your question, but again that may be a typo.

I think you're right about this, and I'd be happy to support a bug fix if you'd like to open a pull request (alternatively, I can begin the PR if you'd prefer).

Regarding the use of (cosine_loss_exponent_tilt / 3.0) rather than simply cosine_loss_exponent_tilt: this is because of the cubic relationship between velocity and power. That is, if the velocity is changed by a factor of $a$, the power is changed by a factor of $a^3$. Since the cosine_loss_exponent_tilt is understood to be (from the actuator disk model) the amount that the power will change due to a tilt misalignment, the "effective velocity" change will be the cubed root of that to ensure that the correct change is applied to the power when using the effective velocity to query the power table.

Cheers,
Misha

[mkrause-strath]

P.S. this may be a stupid question but in three dimensions is it correct to just multiply the yaw and tilt if both are present? I.e. at the moment if there is both tilt and yaw and the exponent is 1.88 for both:

$v_{rotor effective} = (cos(\theta_{yaw}) * cos(\theta_{tilt})) ^ {1.88/3.0}$

I don't seem to be capable of the 3D thinking and trigonometry right now to come up with the answer.

Rephrasing it. If $e_x$ is the wind direction and $v$ is the hub orientation is:

$cos(a)=cos(\theta_{yaw}) * cos(\theta_{tilt})$ in the below picture?

(see https://en.wikipedia.org/wiki/Direction_cosine#/media/File:Direction_cosine_vector.svg for clearer picture)

[misi9170]

Hi @mkrause-strath ,

Great! You'll need to create a fork of FLORIS, commit the changes to a branch on your fork, and then submit a pull request (PR) from you fork back to the main nrel/FLORIS repo. See here and the associated link for details.

I'm not sure I totally follow your question; my gut reaction is that yes, a product of cosines is appropriate because the tilt misalignment and yaw misalignment are always perpendicular to one another. But I may be misunderstanding!

Best,
Misha

P.S. I'm heading out of office until mid-July, but I'll be happy to review your PR when I'm back

[mkrause-strath]

Hi,

I will create a fork then and go through the steps. I'll try to think about the product of yaw and tilt a bit more but I think you're right.

Hope you enjoy your holiday!

Max

[mkrause-strath]

Pull Request Here: #1126

[misi9170]

Thanks @mkrause-strath ! I'll try to get around to reviewing that later this week.