Help with a Script using nested for and fzero function
显示 更早的评论
Hi guys,
I need to implement a script that will find the zeros of a function and after that plot the value of alpha(x axis) vs beta(y axis).
the function is given by:
(cos(phi)*sin(beta-acos(cos_phi))-a)+(a-cos(phi)*sin(alpha-acos(cos_phi))*exp^(-(beta-alpha)/tang(acos(cos(phi))))
The program that i wrote till now is:
for a=0:0.2:1
for cos_phi=0:0.2:1
for alpha=0:0.1:pi
beta_calc= @(beta,alpha,cos_phi,a) (cos(phi)*sin(beta-acos(cos_phi))-a)+(a-cos(phi)*sin(alpha-acos(cos_phi))*exp^(-(beta-alpha)/tang(acos(cos(phi)))))
teste=fzero(@(beta) beta_calc(beta), [pi/2 2*pi])
plot(alpha,beta)
end
end
end
But obviously is not running and im stuccoed :(
Can someone help me?
I am new on matlab so this is a kind of high level for me...
Thanks in advance!
采纳的回答
Walter Roberson
2017-10-6
Code attached.
I had to assume that tang() means tan()
The search range for fzero has to be constrained so that beta >= alpha: otherwise the sign of the term in the exp() can change, leading to the problem that the endpoints of the fzero search are not of opposite sign.
You will notice that all of the graphs are exactly the same. If you work through the formula you posted, you will find that it is of the form
(something - a) + (a - something_else)
which is equivalent to (something - a + a - something_else), and so the a's cancel out.
16 个评论
Jucimar Carpe
2017-10-6
Walter Roberson
2017-10-6
I would suggest that you made a mistake in the formula. Perhaps you could post an image of the formula?
Jucimar Carpe
2017-10-6
Jucimar Carpe
2017-10-6
Walter Roberson
2017-10-6
Your expression uses acos(cos(phi)) which might not be the same as phi, which is a mis-optimization on your part. You need to work with phi. You can take cos(phi) at the time you assign phi, the way I show. This is a single execution of cos() for each phi, and that is more efficient (and more correct) than the multiple acos(cos_phi) that you have.
I plotted on different figures because it is impossible to tell the graphs apart for different a values. Look at what you have:
((cos_phi*sin(beta-acos(cos_phi))-a)+(a-cos_phi*sin(alpha-acos(cos_phi))*exp(-(beta-alpha)/tan(acos(cos_phi)))))
so that is
((cos_phi*sin(beta-SOMETHING)-a)+(a-cos_phi*sin(alpha-phi)*exp(-(beta-alpha)/tan(phi))))
which is
((cos_phi*SOMETHING-a)+(a-cos_phi*SOMETHINGELSE*exp(-SOMETHINGOTHER/ANDANOTHER)))
which is
((SOMETHING-a)+(a-SOMETHINGELSE))
and that rearranges to
SOMETHING - a + a - SOMETHINGELSE
which is SOMETHING - SOMETHINGELSE + (a - a)
and so a has no effect on the calculation.
This is why you need to recheck the calculation. I recommended that you post an image of the formula so that we can look over it to see where you might have made a mistake in your (). For example perhaps (cos_phi*sin(beta-acos(cos_phi))-a) should instead be (cos_phi*sin(beta-acos(cos_phi)-a)) which moves the -a part into the sin() call.
Jucimar Carpe
2017-10-6
Jucimar Carpe
2017-10-6
Walter Roberson
2017-10-6
As expected, you had a bracket in the wrong position.
I have attached a revised version that adheres to the given graph as closely as feasible, including the tick labels and grid, and a good approximation of the colors.
I did not, however, replicate the groups of plots. I was able to determine that the rule for any given value of a was that the plot starts at alpha = asind(a) and runs until 90 + asind(a) (so, the two places where sin(alpha)=a), but I did not manage to figure out what the rule was for the height, so I was also not able to text() appropriate labels into place.
Parts of the graph are mostly empty. There are places where the fzero() is not able to find a solution because the endpoints are the same sign.
This happens in particular for a = 0.8, for which the function starts slightly negative at beta = pi/2, rises to slightly above 0 and falls again, ending notably below zero at beta = 2*pi. Even if one were to adjust the zero-crossing algorithm to be able to find a solution, one would have to ask which of the two solutions would be appropriate.
Jucimar Carpe
2017-10-7
Walter Roberson
2017-10-7
Your sample graph shows plots of a range of alpha_1 values from 0 to 180. I use that entire range as alpha values. However, for alpha values greater than 90, fzero is not able to find solutions for some combinations, and the other solutions are essentially on top of each other.
Jucimar Carpe
2017-10-7
Jucimar Carpe
2017-10-7
Walter Roberson
2017-10-7
The correction to cos_phi_vals and the change to 0.0001 was the only real change in your version.
The num_* is there for efficiency: there is no point in re-calculating all those things every iteration of the loops, especially with nested loops.
Using rad2degree is good for readability but it does have more overhead than *180/pi.
Using .*exp instead of *exp does not matter because you are working with scalars at that point.
Using (tan(phi)) instead of tan(phi) makes the code more difficult to read: although it is not wrong, when the reader sees that you have emphasized the grouping like that, they have to stop and think about all the various circumstances under which the extra () level could make a difference, all the reasons why you might feel it is important to emphasize the grouping, only to realize that there are no reasons...
In the example plot, it is not clear to me whether, for any given a value, only asin(a) to 90+asin(a) are to be iterated over, or if instead the whole 0 to 180 range is to be iterated over but then the result is to be mapped (linearly?) into the asin(a) to 90+asin(a) range for presentation purposes.
Jucimar Carpe
2017-10-7
Walter Roberson
2017-10-8
Okay, we are getting much closer. See attached.
I replaced the direct calls to fzero with calls to a routine intended to search for multiple zeros. The routine divides the given range up into the given number of subdivisions, and looks for all transitions in sign under the assumption that the user specified enough subdivisions that all of the transitions will be present. It finds extracts the transition intervals and does fzero() on each of the intervals to get the more exact location. The routine returns all of the detected transitions.
With the aid of this routine, I was able to find that each combination has between one and three zeros. Some of the zeros are not crossings, being points where the curve just touches zero from below and then goes negative again; I found that I was able to detect those by adding in the special points that are simple fractions of pi (e.g., 1/8 * pi). Or at least I was able to detect them well enough that I was able to find at least one zero for each combination.
I decided to plot all three zeros. In the graph that is produced, the largest of all of the zeros is associated with the ^ marker, including if that is the only zero. In cases where there are more than one zero, the smallest of the zeros is associated with the v marker. In cases where there are three zeros, the middle zero is associated with the + marker.
You will notice that a number of points are drawn below 90, including most of the lower and middle markers. In the plot you see from running the code, none of the upper markers are below 90. That is, however, only because I tweaked the code. In the situation where cos(phi) is sufficiently close to 1, the markers are below 90; in particular for cos(phi) = 1 exactly, the markers are at asind(a), which is never more than 90.
For alpha above roughly 130, there are lower beta greater than 90 for most a values. For a = 0.8 there are also middle beta greater than 90 -- so if the cut-off for drawing is to be [90 360] then you will still get some values plotted.
You might notice that the solid lines (upper beta) appear to vanish at sufficiently high alpha. Those lines are not absent: they are merely all the same positions, so the line drawn last is the one whose color is visible.
Jucimar Carpe
2017-10-8
更多回答（0 个）
类别
在 帮助中心 和 File Exchange 中查找有关 Get Started with MATLAB 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
