MATLAB Need help with some matlab code - engineering project

[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 phiWe want this to recursively narrow down the value of a and aPrime until the difference between recursion cycles is negligible. 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. Note: The problem is somewhere in the section commented off between the big *PROBLEM HERE ^* and *PROBLEM HERE v* comments (there is two possible chunks). What I previously referred to as a and aPrime, are in the code as aBMT and aPrimeBMT.

        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

 

[uart]

Actually that's not recursion it's just a simple "iterative" method. (the term "recusion" has a specific meaning in computer science and is something other than simple "looping" or iteration).

It would help if you described the equations you're trying to solve in a bit more detail as well as telling us exactly what iterative method (fixed point method, Newtons method etc) that you are trying to implement. (BTW it looks like a fixed point iterative method, but you should be telling us this info as it a fundamental part of what you're trying to do)

 

[DyslexicHobo]

Yeah I guess I see how it's not exactly recursion. 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.

 

[uart]

Yeah I guess I see how it's not exactly recursion. 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.

Yeah that's called the "fixed point iteration" (FPI) method and it doesn't neccessarily converge unless the function you're trying to solve has certain properties (I think it needs local gradient less than unity). Your function is pretty messy and complicated so it's hard to tell if it's suitable for FPI. Given however that someone has told you to use the algorithm and said that it should work then it's mostly likely that you've made a programming error in implementing the function. Since we don't know exactly what function you're trying to implement it's obviuosly pretty difficult to tell you where your attempted implementation deviates from the (presently only know to yourself) desired function.

One thing that does stand out as a likely error however is the use of degrees in the trig functions. In real mathematics angles are always in radians.

I'd be pretty certain that you should be using phiRad in all the following expressions.

            aBMT = (1 / (1 + (4*sin(phi)^2) / (Cl * sigma * cos(phi))));
            aPrimeBMT = (1 / ((4*cos(phi)) / (Cl * sigma))-1)

 

[DyslexicHobo]

Wow uart, your post was extremely helpful! I can't believe I overlooked a detail like that. I was originally using phi in degrees because it was required to find Cl and Cd using the functions we found.