Engineering project - simple code - need guidance

[DyslexicHobo]

This is some Matlab code that my group has been working on. None of us really code all that much, and have only had one year of Comp Sci. It doesn't help that this is all our our first times using Matlab, either. :/

The problem is calculating a value for a and aPrime. We're trying to use recursion to calculate the correct value for them. The way it's supposed to work:

wRelative : defined in terms of a and aPrime
phi : defined in terms of a and wRelative
alpha : defined in terms of phi

a is calculated based on an equation we were given in terms of phi
aPrime is calculated based on an equation we were given in terms of phi


What I'm trying to do is calculate a value for 'a' (and aPrime). I use 'a' to calculate phi, and then use phi to calculate 'a'. I want it to constantly re-calculate 'a' and phi so that it converges to a specific value. However, instead of the values converging, they just alternate between 2 values (values with a fairly large difference). Can someone try and pinpoint our error? Thanks.

What I'm trying to do is calculate a value for 'a' (and aPrime). I use 'a' to calculate phi, and then use phi to calculate 'a'. I want it to constantly re-calculate 'a' and phi so that it converges to a specific value.

*If anyone can tell me why my code is posting improperly, that'd be great. For some reason it's not retaining the proper spacing in the "attach code" window.

while(aBMT < 1)% (abs(aDif) > .1) && (abs(aPrimeDif) > .1)) Finds values for a and a' such that they are within .1% of the last recursion
%*PROBLEM HERE v*
wRelative = sqrt(windSpeed^2 * (1-aBMT)^2 + (w^2*r^2) * (1 + aPrimeBMT)^2);
phiRad = asin((windSpeed * (1-aBMT))/wRelative);
phi = phiRad * (180/pi);
alpha = phi - beta;
%*PROBLEM HERE ^*
%FINDING CL
if(alpha>0 && alpha <= 8)
Cl = -0.0001*alpha^5 + 0.0016*alpha^4 - 0.0093*alpha^3 + 0.0219*alpha^2 + 0.0928*alpha + 0.0006;
elseif(alpha > 8 && alpha <= 27)
Cl = -.00001*alpha^2 + 0.0542*alpha - 0.5037;
elseif(alpha > 27 && alpha <= 90)
Cl = -.00000009*alpha^4 + .00003*alpha^3 - 0.0036*alpha^2 + 0.1761*alpha - 1.8521;
else
Cl=0;
end

%FINDING CD
if(alpha >= 0 && alpha <= 8)
Cd = -.00001*alpha^3 + 0.0003*alpha^2 - 0.0003*alpha + 0.0134;
elseif(alpha > 8 && alpha <= 90)
Cd = -.0000000004*alpha^5 + .0000002*alpha^4 - .00003*alpha^3 + 0.0018*alpha^2 - 0.0196*alpha + 0.1616;
else
Cd = 0;
end

sigma = numBlades * chord / (2 * pi * r);
%*PROBLEM HERE v*
aBMT = (1 / (1 + (4*sin(phi)^2) / (Cl * sigma * cos(phi))));

aPrimeBMT = (1 / ((4*cos(phi)) / (Cl * sigma))-1)
%*PROBLEM HERE ^*
end

 

[CycloWizard]

I think there's a fundamental flaw in the logic as you described it. If you have a functional relationship between x and y, then reversing it to solve for x in terms of y should always return the same value. Maybe you meant to say you compute phi from aPrime? Anyway, I'm taking a look at the code and will report back.

edit: I got the loop to end. My values for a and aPrime kept varying, which isn't what I gathered from your post. Of course, I just made up values for the parameters:

aBMT=0;
aPrimeBMT=1;
windSpeed=10;
w=5;
r=5;
beta=1;
numBlades=3;
chord=1;

without a little more background, it's hard to say what you're really trying to do. I could give you some pointers if you're trying to find an optimal value for a or something like that.

 

[DyslexicHobo]

I posted for help on another forum, too (PhysicsForums). A poster there said that what I'm trying to do is called fixed point iteration.

I didn't go into too much detail on what exactly the code is supposed to be doing, because it's pretty specific to a project that I'm working on. It's a wind turbine simulation for my engineering class, and I don't really have the best computer science background (AP comp sci in high school). Basically, we're designing a wind turbine with the following parameters: chord length (cross-sectional width of the blade... NOT the thickness), beta (angle of attack), radius of the blades, and number of blades. Using these 4 parameters, we need to calculate a value called a and a prime, which helps give the amount of total energy able to be harvested from the wind. I won't go into too much detail, but it concerns how the wind behaves as it nears the blade. 'a' refers to the linear momentum of the wind, and aPrime refers to the rotational inertia of the wind.

But as far as optimization goes, a and aPrime *should* converge upon a single value. I found one problem was that we were using phi in degrees instead of radians for the calculations, which could have been a big part of the problem. I'm going to try to fix this tomorrow and I'll post the results.

Thanks for the help!

 

[esun]

Fixed point iteration should be very straightforward. This is typical in solving transcendental equations, e.g., 5 - x = exp(x). If you typed the equations in correctly (and quite frankly, I can't tell since I don't know what equations you're using), iterating shouldn't be a problem. The degrees/radians thing definitely sounds like it could be the problem, though.