Earlier this week, you learned that we are going to use a nonlinear search to find the value of input current that causes future cell terminal voltage to encounter some limiting value. And we're going to use that with the enhanced self-correcting cell model. And you'll learn that we're going to use a bisection algorithm to perform that nonlinear search. In the previous lesson, you learned what the steps are for the bisection algorithm, and then how to implement that in octave. In this lesson, you're going to learn how to apply all of this to an enhanced self-correcting cell model, to compute the power limits. And we begin by recognizing that we still need to specify a function that will be used by the bisection algorithm when searching for the maximum current for discharge, and the minimum current for charge based on a future terminal voltage encountering some limiting value. So we're going to use the enhanced self-correcting cell model, which has a state equation that is linear, if the input current is constant. And we're assuming that the current is held constant over the future time window. So that is, the next state is equal to a matrix A, times the present state pLus a matrix B, multiplying the present input. So in this equation, A and B are both constant if the input current is constant, but they're also functions of the level of input current for the enhanced self-correcting cell model. So we actually need to define octave functions that compute the A matrix and the B matrix for any desired level of input current. So I share with you some code that does this using implicit inline functions that are quite convenient because I don't need to create separate M Files function files for them. So the first line of code defines the function A to be a function of the input current ik. And it says that the output of this function is the diagonal matrix, where the first element is 1 and the second element is the exponential of -1 divided by RC. And the third element is the familiar exponential term from this hysteresis equation. This organization that we're using for the A matrix is assuming that the state vector of our cell model has the state of the charge state on top. It has a diffusion resistor current state in the middle. And it has the dynamic history state on the bottom. The next line of code defines an inline function B, which has input ik as well. So function B computes a matrix that is three rows tall and two columns wide. It assumes that the input is actually not a single scale or value, but a two component vector. The top component of u is the actual cell current. And the bottom component of u is the sign of the cell current. So this allows us to pre-compute the nonlinear thing, which is the sign of the input current. And then multiply the current and the sign of current by a matrix to compute its impact on how this state evolves, over time. When we do this, the first column of B multiplies the actual current, and the second column of B multiples the sine of the input current. So the top left component of the B function, the matrix produce by the B function, is the update for the state of charge equation. And the middle left component scales the input current to update the diffusion resister current. Finally, the bottom right component scale in the sign of the input current, to update the history's estate. One objective of our bisection search is to ensure that the future voltage is between specified limits. So I need to compute a future voltage. To compute the future voltage, I also must be able to compute a future state. And I repeat here an equation I shared with you earlier this week that computes the future state based on present state and a cost in input current over that time interval. And we desire to compute this equation as efficiently as possible, since it's going to be executed many, many times in our bisection loop. So first at every point in time we are bisecting multiple times. And we're bisecting for multiple cells, and we're bisecting over the entire drive cycle profile. So this is executed a lot. We want it to be super efficient. So in this equation, the first term computes this A matrix, which is really the implicit function you learned about on the previous slide to the power of k delta T. And then multiplies that by the presence state vector. Now the matrix A raised to the power of k delta T in general is quite a complicated operation. I have to take a A matrix and multiply it by itself a whole bunch of times. But, in our case it's simplifies because our A matrix is diagonal. And whenever we take a diagonal matrix to a power, the same result can be computed by erasing every diagonal element to that power, and that makes computing the first time really quite straightforward, and quite fast. The summation term requires some more analysis. We're going to first factor the constant B matrix outside of the parenthesis but we still need a way to compute this summation of the A matrix raised to different powers. You might recognize that this is simply a finite geometric summation of the A matrix. And I expect that you've seen these results for scalars, raising a scalar to different powers and then summing up a result. I think you have seen that in a calculus class before, but you may never have seen this result when A is a matrix. So how do we compute a finite geometric summation when I'm using a matrix instead of a scalar? So the steps on the slide show you how to do it. And I'm actually not going to discuss them in detail. Hopefully there is enough information in here that you could decode all of the steps of the method if you chose to do so. I think that the most difficult step is the one where we moved from the right side of the first line of the equation into the second line of the equation. In order to validate this result, I recommend that you write out the summation terms from the right side of the first line of the equation. For some example, value of k delta T that's small enough, for which you can do that, like maybe 3, something like that. Then you compare that to the multiplication listed on the second line of the equation to convince yourself that this is actually working. The main result of this analysis is that we can write the summation in a closed form quite simply. And that's very nice. I don't need to have a for loop to compute this result. And that is the bottom right result here is the final answer. But I also need to mention that in this derivation here, we assume that the A matrix. I don't think it actually has to be diagonal but here it is diagonal, and we are assuming that the A matrix never has the element 1,0 .0 on it's diagonal, never can have unity values on the diagonal. If any diagonal a the A matrix is equal to 1, then the bottom result is going to give us a 0 divided by 0 result, which is indeterminate. But we can repeat this analysis for the special case of an element being equal to 1.0. And come up with a special case result that any time the letter A has a value of 1 in the diagonal we can simply use that special case result. And I'll show you what that is on the next slide. This slide implements the results from the equations in the previous slide in octave. So it very efficiently simulates cell voltage and cell state k delta T samples into the future. The inputs to the function include the constant level of input current and the present state and the functions that implement the matrices A and B for different levels of input current. It also includes the number of samples k delta T into the future for which I desire to compute the voltage and state. And finally, it requires the present temperature and the present cell model. In principle, if I had the temperature and the model structure, I could query that model for all of the individual cell parameter values that I might want using the get param ESC toolbox function that you learned about in the second course. But instead of spending time extracting those values from the model every single time this particular function is invoked, which is a lot of times, I extract those values outside of this function once only. And then I pass them to this function in order to make the code execute a little bit more quickly. The first line of actual code and the function invokes the A function. Remember the implicit A function that you learned about it in this lesson earlier. Using the desired value of input current. And that gives us a constant a matrix that applies for the next k delta T samples for this level of input current. We also invoke the B function to find a B matrix that applies for the next k delta T samples. And then it computes the vector dA, which is just the diagonal of the a matrix, which is important in what follows. We compute a variable ADT, which is the finite geometric summation of the A matrix using the equation on the bottom of the previous slide. But we need to be careful for some special cases. If the A matrix has a value of 1 on its diagonal then the previous result does not work because we encounter a 0 divided by 0 problem. So in that case, we repeat the previous analysis for the special case where the A matrix has a value of 1 on the diagonal, and find that the summation value for that special case is simply equal to the KDT variable. So we need to use a special case always for the state of charge element because the a matrix always has a 1 in that spot saying that the future state of charge is equal to 1 times the present state of charge plus something related to the input. We also need to invoke the special case for the dynamic history's estate whenever the input current is equal to zero. So if the input current is equal to 0 in the ADT matrix, then it has one form, and if the input is not equal to 0, then the ADT matrix has a second form. And that's what that conditional statement is doing here. After computing the ADT matrix, we compute the future state in a very straightforward way. The future state is equal to the diagonal elements of the A matrix raised to the KDT power, all multiplying the present state. Plus the ADT matrix times the B matrix multiplying the vector containing the input and the sin of the input current. The future voltage, once I know what the future state happens to be, is simply the future open circuit voltage evaluated at the future state of charge plus the voltage minus the diffusion resistor voltages and minus the omega voltage drop. So at this point you've learned how to simulate a cell KDT seconds into the future quite efficiently. We could use this inside of a bisection method to solve for when the future voltage is equal to some specified voltage limit. But we also want to make sure that the future state of charge is within some range as well. We don't want to violate some future state of charge limit. So we could use the bisection method to guarantee the future voltage. And then use a second bisection to guarantee to a future state of charge. But we could in fact do it together. And I'm going to share with you how we could do that together. First consider the discharge case. I am defining a function call bisectDischarge here. When the output is 0, then I have encountered either a minimum voltage or a maximum state, I'm sorry a minimum state of charge or both. If the output of the function is less than 0, then both the voltage and the state of charge are above their minimum values, so they're okay. If the output of the function is greater than 0, than either the future voltage or the future state of charge or both below their minimum values and that's not okay. So here is how it works. The first line simulates the cell voltage and cell states k delta T samples into the future. Then we compare the future voltage and the minimum voltage. The minimum voltage minus the future voltage must be negative to meet the design constraints. We also compare the future state of charge and the minimum state of charge. The minimum state of charge minus the future state of charge must also be negative to meet the design constraints. We want to know if either one these terms is positive, because that would violate the design constraints. So if either is positive, then we've violated maybe a voltage limit or state of charge limit, or both. So we need to take the maximum of the two results to see if that's happened. That guarantees that if either results is positive, then the upper the function is positive as well. But if both results are negative, then the output function is negative, and we know that we have not violated constraints. So what we have done is we have created a function that's always negative if the future state and voltage are within limits. And it's always positive if it is outside of limits. So this function is suitable for using with a bisection method to look for where the function is exactly equal to 0. Which is when I'm encountering either a voltage limit or a state of charge limit or both. We can create a similar function when we are looking for state of charge limits. In this case, we need to be careful because the sign is going to be different. So this function returns an answer that is positive when both future voltage and state of charge are inside of limits. And it returns a value of 0 when either the future voltage or state of charge is exactly equal to a limit. And it returns a negative value, if the future voltage or state of charge are outside of bounds. So the main point, again, is that the output has one sign when all quantities are inside of bounds and it has the opposite sign when any of the quantities is outside of the limits. So again, we could use this function inside of a bisection method to seek the value of current that pushes us right up to the limits but not any further. Here I share with you quickly how we can use this with the bisection method. And I will share with you the entire code base in the next lesson. But this should be enough to give you the general idea. First, we define an inline function, h, that takes input argument x. Then sent function invokes either the bisect discharge function using the current input equal to x and all of the parameter values taken from the work space. Then we can use this h function or the bisection algorithm to find a limiting value of input current between imin and imax for some tolerance. To summarize this lesson, you've learned that to find current limits, we will use the bisection algorithm with the enhanced self correcting cell model. You've learned how a bisection will be used inside of an overall approach to find the current limits and hence the power limits. To use the bisection method with the enhanced self correcting cell model, we must be able to simulate the cell KDT seconds into the future. And you've learned an efficient way to do this in simulation. You've also learn how to create a bisection cost function that will simultaneously look for limits of future state of charge and on future voltage. And you've seen how we are going to invoke this bisection function from the previous lesson, using the cost function from this lesson. So all of the pieces have now been defined. And in the next lesson, I will share with you how to put everything together to make an overall power limits estimator based on the bisection search method.