[MUSIC] In general, every recursive implementation of a function has an equivalent iterative implementation. To add insult to injury, the Iterative version is almost always faster than a recursive one. Because MATLAB doesn't have to make multiple function calls with the associated overhead of creating the stack frames and copying input and output arguments back and forth. So why bother with recursive functions? Well, for many real-world problems, the recursive solution is much simpler conceptually, and consequently, it's easier to implement. Let's see a fun example, the Sierpinski triangle, which is named after the Polish mathematician who invented it over 100 years ago. We're going to go into considerable detail with the Sierpinski triangle because it's a great way to gain a deep understanding of recursive functions, but it's also a very interesting shape. Here's an illustration from Wikipedia. It has the overall shape of an equilateral triangle, subdivided progressively into smaller and smaller equilateral triangles, as you can see here. Or you can look at it as very small equilateral triangles, arranged to produce progressively larger and larger equilateral triangles. It's an example of a class of shapes called a fractal. Meaning that if you zoom in, you'll see the same pattern repeating over and over. And here's its recursive definition. The base case is a Sierpinski triangle of depth 0, which is one equilateral triangle. The recursive case is a Sierpinski triangle of depth n, which is three Sierpinski triangles, each of depth n minus 1, arranged as shown here with one at the bottom-left, one right next to it at the bottom-right, and one on top sitting on the other two. This simple recursive definition produces elaborate shapes, and we can write a recursive MATLAB function to implement it. In fact, we've already written one for you, let's take a look at it. Here it is in the editor window. Its file name is highlighted over here in the current folder window, and its name is sierpinski.m. And as you can see, it contains two functions, our main function named sierpinski right here, and a subfunction named sier. I should warn you that we're going to analyze Sierpinski in great depth. In fact, before we're done, you might say that we're analyzing it to death. But stay with us, there's a good reason. As we said, it's a fun example, but it's much more than that. It's not only fun, it's a nearly perfect example for teaching you the power of recursion. The power to do very complicated things with very simple algorithms. It's also nearly perfect for showing how recursion operates simultaneously on multiple levels. And finally, it's perfect for demonstrating how recursion can turn on you and chew up computer resources. I want you to understand recursion, not just become acquainted with it, really understand it. And we're going to get you there with Sierpinski. Before we analyze it, let's make sure it works. The main function requires an argument specifying the depth, let's try depth 5. That looks pretty good, and let's try the minimum depth of 0. Just a simple garden variety equilateral triangle with the side of length equal to 1, which is the base case of the definition. And the next level up, which we get by calling sierpinski with a depth of 1, Gives us three equilateral triangles, sides equal to one-half, arranged as we saw for the cursive case, one at the bottom-left, one at the top middle, and one at the bottom-right. Okay, it works, let's analyze it. First, the main function is named sierpinski, but it doesn't really do Sierpinski. It just sets the stage for its subfunction, sier, which is the real star of the show, and then calls that function to perform the Sierpinski magic. That magic involves plotting hundreds or thousands of equilateral triangles, each one just the right size and in just the right position to produce the Sierpinski pattern. Let's go over what the main function does. First, it takes one argument, depth, but it doesn't do anything with it except pass it on to sier down here. Next, it sets s, which will be the size of each side of the Sierpinski triangle, to be equal to 1. And sets c, which will be the position of the center of the Sierpinski triangle, to be the origin of the figure axis right over here, 0 for x, 0 for y. Then it executes the clear figure command, clf, whose job it is to prepare a nice, empty figure ready for plotting. And it uses the axis function to adjust the plotting area so that the entire Sierpinski triangle will be within the field of view. And It includes the equal attribute so that what it plots won't be distorted. It uses the title command to display a title like this one, Which it constructs with the sprint f command. And it uses hold on to allow multiple triangles to be plotted on the same figure. Then it calls sier, which does all the hard work. And finally, it turns hold off, in case the user wants to plot something else. And now let's see what sier does. First of all, it takes three input arguments. The first is d, which is the depth of the triangle. The second is s, which is the size of one side of the triangle. And c, which is the position of the center of the triangle. Then it starts to work, first, like any good recursive function, it decides whether to enter the base case or the recursive case. You'll remember that Rfact did that too by determining whether its input argument n was equal to 0. Well, sier does a similar thing, it determines whether its input d is equal to 0. If d is equal to 0, the base case is run. And that base case simply calls the plot function to draw just one equilateral triangle of size s centered at c. I would invite you to study these arguments to plot on your own, to convince yourself that these arguments do in fact produce an equilateral triangle of size s centered at c. if d is nonzero, then the recursive case is run. The job of the recursive case is to draw three Sierpinski triangles Arrange like this three ordinary triangles here, here and here. We might call them sub-Sierpinski Triangles because together they make up one full Sierpinski Triangle. The size of each of these three sub-Sierpinski Triangles is set to one-half. So that they will just fit inside the full Sierpinski Triangle side equal to 1, which is the value of s that was input to sier. So the first step in the recursive case is to cut s in half. Next, the height of the sub-triangle is calculated. And then three recursive calls are made. That will produce three Sierpinski Triangles at the positions indicated in these comments. Bottom left here, top middle here, bottom right here. To do that, each one must be shifted relative to the center of the full Sierpinski Triangle. Down into the left here, up to the middle here, and down into the right here. And I would invite you to stare at these three shifts here, here and here to convince yourself that they do just that. And by passing itself the correct sizes and positions, that's all sier needs to do to plot two three sub-Sierpinski Triangles correctly. But I've left out what might be the most important parameter in each one of these recursive calls. The first one. This first parameter doesn't have anything to do with the size of the positions of the sub-triangles. But it's very important because it determines how many triangles will be plotted within these sub-triangles by determining how deep that recursion will go. The depth is reduced by one for these calls as you can see, and it will be reduced again in the next level of calls and so on. And all the calls that receive a depth of zero will run the base case. And we'll plot one equal lateral triangle with no recursive call. In fact, notice that we never call plot in a recursive case, each triangle is drawn in the base case. Okay, it took me a few minutes to say what sier does. Still though, I don't know about you, but it's amazing to me that these 11 lines of code can produce something so complicated and there I say beautiful. That guy fair enough. Sierpinski was no Michelangelo or Wong Mang. He wasn't an artist at all, but you gotta admit he came up with some pretty impressive drawings. Let's see another one. This picture looks like a big triangle with smaller triangles inside it, and smaller triangles inside them. But in fact, it comprises only the smallest triangles and their placement in the figure forms the larger triangles. If you want to see the order in which the triangles are drawn, we can put a break point on this line in our base case where we call it plot function. Before I do that though, I'm going to insert one line temporarily right here before we call the sier It will draw a red outline of the entire Sierpinski Triangle. So you can see where in that big triangle, each little tiny triangle was plotted. Okay, let's put the break point in. And run it. And now let's click Continue. And there is the first little triangle at the bottom left, I'll click a second time and a third time. And now we've seen Sier draw three little triangles. Of course, each one was drawn in this base case, which is executed when the depth D is 0. And these three depth zeros Sierpinski Triangles taken together form one depth one Sierpinski Triangle. The next three clicks will give us three more depth zero triangles. For a second depth one Sierpinski Triangle. And three more will complete. A third depth one Sierpinski Triangle. Together these three depth one Sierpinski Triangles form a depth two triangle, and so on. At the lowest depth which is zero, the order is bottom left, top, middle, bottom right. And at the next depth which is one, the order is bottom left, top middle, bottom right. And this depth to triangle is the bottom left part of a depth three triangle, which is the bottom left part of the whole depth four triangle. And this depth two triangle comprises three depth one triangles each of which comprises its own set one, two, three of depth zero triangles. It was drawn because sear was called once in the main function with its first argument equal to four. And then once inside itself in this first recursive call with its first argument equal to three. And then once again in this first recursive call, with its argument equal to two. And then three times here, here and here with it's argument for equal to one. And finally, three times here, three times here, and three times here with its first argument equal to zero. So Sier has already operated at five depth levels from four down to zero to draw these triangles. And it shows you that in a recursive function, the order of things is not the simple line by line order that we are accustomed to for non-recursive functions. Well it is, but with recursive function a single line may be operating simultaneously at many different levels in the recursion. So you have to think and think hard, which is exactly why I like recursion. To encourage you to think hard about recursion, I'm going to let you watch it draw an entire depth for triangle, one tiny triangle at a time from beginning to end. And instead of pounding on the continue button over and over, I'm going to clear the break, Quit debugging. And just slow things down by inserting a pause of two-thirds of a second, Right here before the plot statement. And here we go, It has a nice one, two, three rhythm like a lullaby, in fact I'm getting kind of sleepy watching this [MUSIC] I don't know whether watching see or drawing triangles one by one helps you understand recursion or not. I think it helps, and I hope it does. But if it doesn't, at least it might help you if you have trouble falling asleep. So back in the code over here, in order to draw all these triangles, the main function sierpinski calls sier just one time, passing in the value 4 as a depth, after which sier makes three depth 3 recursion calls. Here, here, and here, and here's an interesting question, how many total times a sier call before the second depth 3 call? The immediate reaction might be to the first being the original depth 4 call here and the second being the first depth 3 call here. But that ignores the fact that this first depth 3 call draws an entire depth 3 sierpinski triangle at the bottom left. And in order to do it, it makes three depth 2 calls and each of those three depth 2 calls makes three depth 1 calls. And each of those nine depth 1 calls, makes three depth 0 calls. And those 27 depth 0 calls finally, plot the tiny triangles. If you count them out that's 41 calls of sier before we get to the second depth 3 call here. And then this second depth 3 call, draws another entire depth 3 sierpinskit riangle at the top middle which requires 40 more calls. And then finally, the third depth 3 call, draws the third depth 3 sierpinski triangle. For a grand total of 121 calls to sier to draw the full depth 4 sierpinski triangle. That's a lot of calls of sier in one execution considering the code has only four calls, one, two, three, four, and no loops. To see the depth 3 sierpinski triangles appear in the order that they're drawn, I'll get rid of the pause, And save the file. And put a conditional break at the beginning of sier, we're going to make it stop when d==3. Okay, let's input 4 to sierpinski one final time. Okay, we're waiting at the conditional break, and we can see in the workspace over here that d is equal to 3. So we're ready to draw a level three triangle. When I hit Continue, It appears, then I hit Continue again, and the second one appears. And when I hit it again, there's the last one, and there's our level four Ssierpinski triangle, comprising three, one, two, three, level three triangles. Okay, let's remove the break. I'll just click twice on it, and let's run sierpinski with a larger depth say, 7. Whew, lots of detail there. And we can look at that detail by coming up here, hitting this little plus sign and then clicking, clicking, clicking, clicking, and we see all those little tiny triangles. I clickHome, we go back to where we were before. If we go even deeper, Let's do 8. MATLAB slows down quite a bit, it takes a lot of CPU time. Let's see exactly how many seconds. You can look up the function CPU time, but it just returns the total CPU time that MATLAB has used since you started it up. Over here in the workspace, you can see that t naught equals 213.38, meaning that from the moment I started this current MATLAB online session until just before I call sierpinski here, it had used 213.38 seconds of CPU time. And to get the time for the command itself, we just call it before and after then subtract. Okay, it's time to wipe all this stuff off my command window. And I'll clear the workspace too. This relatively long 7.54 seconds of CPU time was required because sier was called thousands of times. To figure out how many times, we note the for a call at a depth of d, The number of calls to sier is equal to 1, plus 3, plus 3 squared, plus 3 cubed, plus 3 to the fourth power and so forth, all the way to 3 to the eighth power if d is equal to 8. That can be written formally this way and we can calculate this sum for d equals to 8 in one statement. And if you're good at summation formulas, you can work out this formula for the sum. And we can check it with MATLAB for D equal to 8, and we get the same number, so the formula checks out. Either way that's 9,841 calls of sier, and each call produces one Sierpinski triangle, including the biggest one, all the intermediate ones, and all the tiniest ones. So in this figure, you can see 9,841 triangles, even though only a subset of them were actually plotted, the tiniest ones. And how many of the tiniest ones are there? Well, since only the calls at the deepest level actually call plot, the number's 3 to the 8th power, 6,561 tiny triangles. Which is approximately two-thirds of the 9,841 calls to sier, you remember from our study of the stack with rfact, the for each of these 9,841 calls, a frame is added to the stack. And when the called function returns, the frames popped off the stack, you may also remember that with rfact the stack reached a height equal to the number of calls. For seer, however, the situation is much different. The case we're looking at with 9,841 calls of seer resulting from an original call with a depth of eight, the stack reaches a maximum height of just nine. Whenever it reaches that maximum height, the latest instance of the function returns and its frames popped off the stack before the next call. The maximum stack height for seers is always equal to one more than the value its first argument D. And the maximum height for rfact, is always equal to one more than the value of its first argument. Well, its only argument, which is n. The functions are similar in that regard, the reason that the stack got so tall with rfact was that the inputs we gave it were so large, like 69,898, which caused rfact to issue 69,898 recursive calls. When you give rfact an argument of only 8, the maximum stack height is only 9, and the total number of recursive calls is only 8. When you give sier an argument of only 8, the maximum stack height is also only nine, but the total number of recursive calls is 9,840. So what's the difference between rfact and sier that causes so many more recursive calls in sier than in rfact? Well the difference is that our fact executes only one recursive call in the recursive section of its code, whereas sier executes three,one,two and three. When there's just one call the number of recursive calls will increase by 1 when the level of recursion increases by 1. If two recursive calls are executed, the number of recursive calls will double when the level of recursion increases by 1. If there are three as there are for sier, the number will triple, and so on. For functions like rfact with just one recursive call and the code, the number of calls increases linearly with the depth of recursion. For functions like sier or with two or more calls, the number of calls increases exponentially, which means that it will increase not by a fixed amount, but by a fixed factor. By the way, I don't know about you, but I'm getting very tired of reporters and commentators, trying to sound smart by using words they've heard actual smart people use. Like exponential, when they're in fact just guessing what the words mean and guessing wrong, here's a good example. Forecasters fear that the destruction of just one large oil field in Saudi Arabia, could cause the world's oil prices to increase exponentially. It'd be so nice if they would just try sticking to words they understand like, oil prices can increase a whole lot. And if that doesn't sound smart enough, they could say increase drastically or dramatically. Increasing exponentially mean something specific, namely increasing by a constant factor, or at least a relatively consistent factor. The factor can be 3 or 2 or 1.01, it just has to be greater than 1, and you have to have some sequence of three or more values. One jump in a value can never be exponential, no matter how big it is. Exponential growth is very common in recursive operations, and it always has the same upward curving pattern, we can see it in a plot of CPU time versus depth for sier, let's get some CPU times to plot Here I'm making one throwaway call to Sierpinski, before I start timing things. I'm doing that because MATLAB often does some setup on the first call and a sequence of calls of the same function. We don't want that setup time to count in one call, and not in the others. After the setup is done, the for loop calls sierpinski with progressively larger and larger depths, storing the times in the vector T. I'm going to hit return. Well, I'm not going to make you wait for this, let me fast forward time for you a little bit. Okay, now let's plot T, And I'm going to turn the grid on, This is just the sort of curve you'd expect when times increase exponentially. It starts out increasing slowly, but at the end it has really taken off. From the shape of this plot, you can see that exponential growth can slow things way down, and I'm going to prove that to you by doing something that we recommend that you don't try at home, especially if you're using an installed version of MATLAB. Usually we encourage you to do everything that you see me do, but not this time. If you do what I'm about to do, you're going to be sorry you did, for one thing you're going to hang up MATLAB. If you do what I'm about to do with an installed version, you'll probably have to reboot your computer, it sounds like clickbait, right? You're not going to believe what happens when Mike Fitzpatrick does this, but in this case, I'm serious, you might even lose some work so don't do it yourself, just watch me do it. As we pointed out in our introductory course, we're professionals, and by the way, I do all my own stunts, and here's what I'm going to do. I'm going to set my computer on fire, okay, ready? Well just kidding, I'm not actually going to set any computers on fire, but I am going to increase the temperature of the computer that my session of MATLAB online is running on wherever it is. I can assure you that, because it's about to get very busy, I hope it's fan is working. Because I'm going to run sierpinski with a depth of 12, last warning, don't do it yourself. Here goes. Well here again, I'm not going to make you wait. I'm going to skip you forward in time a few minutes. And we're back. As you can see, increasing the depth by three increased the time a lot, a whole lot. A depth of 12 required 507.78 seconds. That's about eight and a half minutes. And that's a long time. I know, I had to sit here and wait for this thing to finish. And it hadn't finished yet. Because we zoom out, we can see that we're still waiting for the triangles to be plotted. Why does this thing take so long? Well, first off since there are three recursive calls in the code, we would expect to see an increase in the CPU time by a factor of three for each increase in the depth of one. Well, we've increased the depth by three. So, okay, this is what I was warning you about. We've crashed MATLAB, or as this message says, there was an issue with MATLAB. Well, the issue is that MATLAB could not plot the triangles and it recommends that we contact technical support. Well, that's not so easy for you but it is very easy for me because I have a MathWorks guru on call. So I've done that for you. I paused my recording session and contacted him. And he said that the problem is that MATLAB's graphics system, while it can do extremely complicated plots, is not designed to handle this many plot calls in one application. How many plot calls is too many? Well, he didn't say. But applying our formula of 3 to the 12th power gives 531,441 triangle plots. So we know that's too many. The hang up occurred while the figure was being updated after sierpinski had done all its work, generated all its plotting instructions and returned. The display system then tried to process all those plotting instructions, and that's when MATLAB crashed. Well, that's an explanation, not a solution. So is there a solution? Yes there is, and we'll show it to you in a few minutes. But in the meantime, let's get back to our analysis of this exponential growth. As I was saying when the lights went out, since there are three recursive calls in the code, we would expect to see an increase in the CPU time by a factor of three for each increase in the depth by one. And we've increased the depth by three beyond the last point in our plot. So we would expect to see the time increased by 3 times 3 times 3, which is a factor of 27 relative to the 19 seconds shown in the plot for the depth of 9. That's very close to what happened here. The elapsed CPU time increased from about 19 seconds to about 508 seconds, a factor of 26.7. Time to clean up again. And let's close the figure too. And I've cleared the editor window too. So we have a clean slate. It's time to fulfill our promise to show you the solution to the crashing problem. The solution is provided in this m file which I've added over here. It's called sierpinskiE.m with an uppercase E at the end which stands for efficient. We'll double-click it. And before I explain it, I'll run it. Not bad, instead of 507.78 seconds, it took 1.56 seconds. And as an added bonus, it actually plots and just to confirm that there is an actual full sierpinski triangle here, let's zoom in on the thing. As you zoom in, it just looks like you're seeing a bunch of copies of the same thing that you saw before we started zooming. But eventually, if we get in close enough, and by the way, I'm zooming with the scroll wheel on my mouse, you can finally distinguish the individual tiny triangles. And you can also believe that this function did draw 500,000 of them in a second and a half. So what's the difference? Well for openers, sier does no plotting, none at all. Only the main function now calls plot. Up here where it's always plotted the outline and down here where all half million triangles are plotted in just one call to plot. The approach this time, which is based on a version that Acos created, is that sier just calculates the points that need to be fed to plot instead of plotting them and returns them as an output argument pts, pts stands for points. And it's a two by an array of points in which the x values are in its first row, and the y values are on its second row. When main gets its hands on those x and y values, it passes them directly to plot, which plots the entire sierpinski triangle. Acos' approach takes advantage of the fact that plotting is far more efficient both in execution time and memory space when everything is plotted at once. The same figure takes far longer and far more memory when many individual features are plotted one by one in individual calls to plot. But the downside of this approach is that you don't get to watch the triangles appear one by one in rhythm with Brahms' lullaby. But the need for that actually turns out to be surprisingly rare. Let's look at sier's base case. Let me close the figure so we can see all of it. In the base case, instead of plotting a triangle, the x and y values are calculated for the triangle here and here. And then appended to a set of points that were passed into the fourth argument to sier when it was called. Those points were calculated in a previous activation of sier except for the case in which the main function was the caller. The main function doesn't do any calculation of points. So it just passes in an empty array for that fourth argument. If you examine sierpinski and sierpinskiE closely, you'll see that the four values for x and the four values for y plotted in the base case of sierpinski are exactly the same as the x and y values appended to the points array here. But here in sierpinski, a value of nan is added after the last x value. And after the last y value. These not a number of values serve as a signal to plot not to continue the line from the end of one triangle to the start of the next. You can look that feature up in the help page for plot, it's way down at the end. And now let's look at sier's recursive case. The first call gets the points that were input up here as its fourth argument. And then it does its work, which is to append all the points for its Sierpinski triangle of depth d minus 1 to the points that it received, and return the result. That result is assigned back to pts. And the second recursive call gets as input the set of points that were just returned by the first recursive call. It appends the points for its own depth d minus 1 Sierpinski triangle to those, and its result is also assigned back to points. Finally, the third call does the same thing. And as a result, pts now contains three depth d minus 1 Sierpinski triangles. Which together make up one depth d Sierpinski triangle. Nothing's been plotted, not one point in one triangle. But all the points for all the tiny triangles that make up a depth d triangle are right there in the output argument pts. The speed up with his plot everything at once approach is enormous. But that doesn't affect the shape of the growth curve. The number of recursive calls is still going to increase exponentially with increasing depth, and so is the execution time. Let's plot the time for depths from 1 to 14. I'm going to speed you up a little bit here. There, here again, it starts out increasing very slowly, but then begins to accelerate, and really takes off at the end, just like any good exponential growth. So both the plots that we've seen have had the classic exponential form, but there's a big difference in their time values. The first Sierpinski function hit 19 seconds at a depth of 9. And this sierpinskiE version is much less than 1 at a depth of 9, and comes in at less than 14 seconds for a depth of 14. Well, I'm almost done with this sierpinskiE function, but there is one last detail to deal with. And that's bad input, which is to say negative or fractional deaths, either of which will send this function into infinite recursion. We could handle those cases in the same way that we did with Rfact, by checking the inputs and quitting if they're bad, but there is a simpler way. Just change the base case from d equal to 0 to d less than or equal to 0. I'll let you play with that. All right, I'm going to clean up, And congratulations, you have made it through our entire analysis Sierpinski triangles. And I hope you can now appreciate how simple recursive code can yield complicated results, and how it can operate on more than one level at a time. And how it can gobble up computer time with exponential growth. If you now appreciate these three ideas, then you understand recursion, really understand it, and our job is done. And you've even learned a little something about how to plot quickly. Looking back to our first recursive example, the factorial, you'll remember that after we finished a thorough examination of the recursive calculation with Rfact, we showed you how the calculation can be done more simply with an iterative version. If you're hoping to see us do that with Sierpinski, you're going to be disappointed. Because Sierpinski is an example of a calculation for which the recursive version is far simpler than any iterative version. In fact, that's another reason we chose it as an example. We're not going to show you an iterative version. In fact, I don't even know an iterative version. It would be very painful to try to figure out all the correct coordinates with all those thousands of tiny triangles with loops. It can be done, I can assure you of that, because it can be proven that any recursive function can be done iteratively. But writing the iterative version is not worth the effort. And that alone should convince you that learning recursion is worth the effort. [MUSIC]
