 Hi, in this lecture we will learn to implement a slight variation of Gaussian elimination method called Thomas algorithm which is applicable if the coefficient matrix A is a tridiagonal matrix. We will first see how to do Thomas algorithm and then we will also learn to implement Thomas algorithm in this lecture. Let us consider a tridiagonal system. A tridiagonal system has nonzero elements only at its diagonal part and then the lower diagonal part and the upper diagonal part. Other than these all other entries are 0. This is shown in this system. Here I am only showing the left hand side of the system Ax equal to b. So, this entire thing should be equal to the b vector on the right hand side. Here you can see that the diagonal elements are nonzero denoted by beta i's and similarly the lower diagonal elements alphas are nonzero shown here and the upper diagonal part. Gamma's are also nonzero entries in the coefficient matrix. Remember few of them may be 0, but generally we will assume them as nonzero elements other than that everybody else will be 0. So, how the structure of the matrix will look like therefore we have diagonal element which is taken to be nonzero then the lower diagonal elements and then the upper diagonal elements. All other elements are 0's here and this is how a tridiagonal system will look like. So, the first element of the matrix is denoted by b 1 and the first row second element is gamma 1 and all others are 0's and then in the next line the first element is alpha let us call it as alpha 2 then beta 2 gamma 2 and then all other elements are 0 like that it goes and finally at the last row you will have beta n and previous to that you will have alpha n that is what you can see here and this is how the structure of the tridiagonal system will be. Tridiagonal systems often occur in certain applications especially when you try to build a numerical method for certain PDEs we apply finite difference method or even finite element methods the PDE can be approximated by a linear system and that linear system in some cases will turn out to be a tridiagonal system and therefore solving tridiagonal system is a very important problem Thomas algorithm is often used in such cases Thomas algorithm is a slight variation of the Gaussian elimination method the underlying idea is exactly the same but we will take the tridiagonal structure and adopt the Gaussian elimination method for this. Let us try to understand the Thomas algorithm the first equation is given like this beta x 1 plus gamma x 2 equal to b 1 we have to rewrite the first equation like this where e 1 is given by gamma 1 by beta 1 and f 1 is given by b 1 by beta 1 now once you reduce the equation 1 then let us go to the equation 2 in equation 2 we have only 3 terms all other terms are 0 right we have alpha 2 x 1 plus beta 2 x 2 plus gamma 2 x 3 equal to b 2 now what you do is you have to do an elimination process just like what we do with the Gaussian elimination method to eliminate the variable x 1 here once you do that then the resulting equation which involves only x 2 and x 3 on the left hand side will be finally return in this form with an appropriate definition of e 2 and f 2 you can check this how we are getting it just do the elimination process similar to the step 1 in the Gaussian elimination method then you will get this equation you have to replace this equation with the second equation of our original system. Now you continue like this in general the jth equation is written as x j plus e j into x j plus 1 equal to f j and once you write this then the next equation that is j plus first equation is again put in this form with e j plus 1 given like this and f j plus 1 given like this. So, this is the general form of the Thomas algorithm once you finish the elimination from 2 to n minus 1 remember the first equation is done directly there is no elimination process for that we are only rewriting the first equation from the second equation onwards we have the elimination process just explained in the case of second equation you have to carry over a similar idea for all the other equations starting from the second equation up to loss, but one equation that is up to n minus 1 equation nth equation again can be done in a rather direct way and you can get the nth equation like this note that we have to do an elimination process in the nth equation also, but once you do that the nth equation is just involving x n and therefore, we can get x n directly from this equation right once you get x n from the nth equation you can substitute that into the n minus 1 equation and get x n minus 1 once you know x n minus 1 and x n you can substitute them into the previous equation and get x n minus 2 and so on. So, in other words after doing this elimination process you can get the solution of the tridiagonal system by doing the backward substitution process right. So, this is the idea of Thomas algorithm let us try to see how to code this method let us first see what are all the inputs for our code well the dimension of the system n is the first input I am just taking it as 4 and then instead of defining the entire matrix we can just define the lower diagonal elements upper diagonal elements and the diagonal elements recall that we have used the notation alpha beta and gamma for lower diagonal and upper diagonal elements right the same notation I am using here also recall that the equation will look like beta 1 gamma 1 0 0 then alpha 2 beta 2 gamma 2 0 and then 0 alpha 3 beta 3 gamma 3 and then 0 0 alpha 4 beta 4 this is the tridiagonal system in the case of n is equal to 4 right just keep this structure in mind and see how I have defined alphas, betas and gammas I am defining alpha as a list or array where I have put 0 at the first position because there is no alpha in the equation 1 right that is why I put 0 here and then I am taking alpha 2 as 1 alpha 3 as 2 and alpha 4 as 3. Similarly beta is defined like this and you can see that beta has 4 components it appears in all the equations right and similarly gamma has only 3 components it starts from the first equation and goes up to the third equation lost equation will not have gamma that is why I have put 0 on the lost component of gamma and then I am taking b as an input and I have just given it as 1 2 3 4. Now let us go for the elimination process if you recall we have to obtain 2 parameters one is e and another one is f in the method e is sitting here and f is coming from the right hand side elimination process right. Let us see how these things have to be incorporated in the code well we will use this formula for e 1 and f 1 and then for e 2 to and e 3 we will use this formula again for the lost equation we have to plug in this expression. Let us see how to do that we will first initialize e and f as a 4 dimensional vector and we will initialize x as again a 4 dimensional vector you can also write x is equal to 0 comma 0 comma 0 comma 0, but just to show that we can also initialize these vectors in this form I have just written like this these two are equivalent in python and now we have to define e 1 right that is nothing but gamma 1 by beta 1 that is what we have seen in the method that is what is written in this line remember e 1 mathematically means e naught in python because python starts its index from 0 and goes up to n minus 1 I am always emphasizing this because this is very important and often it is a confusing terminology therefore I am just emphasizing it again and again similarly f 1 is defined as b 1 divided by beta 1 and that is what given here now once we have the first equation second equation onwards till n minus first equation we can go in a loop that is what we are theoretically writing here from j equal to 2 to n minus 1 we can define the parameters e and f in a loop let us see how to do that for that we will use the for loop that runs from 2 to n minus 1 right theoretically j runs from 2 to n minus 1 and in python equivalently it has to run from 1 to n minus 2 right that is what is happening here remember range of 0 comma 10 will create a sequence starting from 0 1 2 up to 9 only it will not include 10 that is what we have seen in one of our python session that we have to remember this confusing terminology that when you say range from some number 2 some number n it will go from this number till one number less than this right similarly here it starts from 1 and goes up to n minus 2 only it will not include n minus 1 in the sequence that the range generates therefore in order to create a sequence like this we have to use this command and we are taking j varying from 1 to n minus 2 mathematically that is equivalent to saying 2 to n minus 1 that is what precisely we want in our method and I am just going to define e j as if you recall what is the expression for e j plus in gamma j plus 1 divided by beta j plus 1 minus alpha j plus 1 into e j and a clever observation is that the denominator in e j plus 1 is the same as the denominator in f j plus 1 also therefore what you can do is instead of defining e j and f j with the denominator calculated separately each time you can just calculate the denominator once store it in the variable d and then use it here as well as use it in the denominator of f in that way what you are doing is you are saving some operations right you can also put this expression explicitly here as well as here but that will unnecessarily involve one additional arithmetic operation just to avoid that I have computed once and then plugged in those values in the expression of e and f now you can see that we have found e naught f naught and then all the other is right from j equal to 1 to n minus 2 we are left out with the last equation n minus 1 that is the last equation as far as the python is concerned theoretically nth equation is the last equation right and from the last equation we can directly get the last component of x that is x n which is in python's notation x n minus 1 and if you recall that can be directly obtained from this last equation here you can see that x n can be written as this divided by this quantity right so that is what I am doing here x n minus 1 which is nothing but our x n is equal to the right hand side divided by what was sitting in the left hand side is now brought to the right hand side once you have found x n minus 1 then getting the other expressions are quite easy you found x n minus 1 which is mathematically x n that you will put in the x n minus first equation to get x n minus 1 that is what I am doing here now I am creating a loop j ranging from n minus 2 that is the first element of the range goes till 0 so it has to go till 0 n minus 3 and so on up to 0 it has to go that is why I put in range minus 1 and the increment is minus 1 here it means it will keep on reducing the number till it reaches 0 so that is what is meant by this command and now for each j what we are doing is x j is equal to f j minus e j into x j plus 1 so that is how the equation is going because if you recall we have written x j plus e j x j plus 1 equal to f j that was the equation and you know x j plus 1 because you know x n from there you can get x n minus 1 and so on so this loop will go and do the back substitution process and finally it prints the vector x here I have just put this command in order to get only the x values rounded to 5 decimal digits for that I have used this command now whatever the python prints from this command will be restricted to 5 digit rounding let us try to run this program and see the output you can see that the program is executed and this command had printed the solution of the tridiagonal system as this you can see that it is restricted to only 5 decimal digits in each component of x because of this command and that is the Thomas algorithm implementation before ending this class let me also discuss an important program which I wanted to discuss in our previous python session but because of the lack of time I could not discuss it there but I would like to discuss it here you can see that in one of our theory class we have taken this function and we have calculated the value of this function at the point x is equal to 1 lakh using 6 digit rounding now how to implement this calculation as a python code that is the question let us try to understand how to implement this calculation on python well to define the function you can also use this command which will create a function f and assigns this expression symbolically and you can evaluate the function value by just plugging in f of x here right so that can give you in one go the value of the function f which is equal to x into square root of x plus 1 minus square root of x remember square root is not a python command it is inbuilt in one of the libraries or modules whatever is say and that module is numpy therefore I am calling numpy and then executing the square root here and we want to find the value of this function at the point x equal to 1 lakh so I have also assigned the variable x equal to 1 lakh now this line will get the value of f at 1 lakh with of course 52 digits in mantisa because generally python computes values in double precision but what we want is to compute the value of f at x equal to 1 lakh using 6 digit rounding how to do that well you can do it using this command the command is float open bracket format open bracket whatever value you want to give and then this point 6g suppose if you want to use 7 digit rounding then you have to use 0.7g or any n digit rounding you have to use dot ng that n value you have to plug in here that is the format what we are doing is first we are finding f l of x right if you recall in one of our previous lectures we have explained how to do arithmetic operation using n digit rounding you first have to round the number x and then you add that number with 1 to get this x plus 1 expression now what I am doing is I am rounding x right to 6 decimal digits storing that value in the variable x 1 and then I am adding 1 to the value x 1 then I am again finding the 6 digit rounding of x 1 plus 1 that is what I am doing here and that value is stored in x 2 and then I am taking square root of x 2 that is I am computing this term here after computing square root of x 2 again I am taking 6 digit rounding of that and storing it in x 3 right once you get x 3 you keep it on one side and go to do this second term remember the argument of the second term is f l of x that is stored in x 1 therefore I am taking square root of x 1 and then doing a 6 digit rounding of that and storing it in x 4 now x 3 precisely has the 6 digit rounding of this term and x 4 has 6 digit rounding of the second term now what I am doing is I am subtracting these two terms and again taking the 6 digit rounding of that term right and that is stored in x 5 so finally you have the complete 6 digit rounding calculation of this term stored in the variable x 5 and then remember already f l of x is stored in x 1 so finally you are doing x 1 into x 5 and then doing the 6 digit rounding of that term and that is stored in x 6 therefore the final answer that we want that is we want the value of this function at the point x equal to 1 lakh with 6 digit rounding calculation is finally stored in x 6 and that is what I am printing in this line and I am also printing the value of the function f directly with 52 digit rounding because it is a double precision calculation and let us see the output of this program the output of the function with 6 digit rounding is 100 whereas the more accurate answer is 158.113 if you go back to our class where we discussed this problem there also we have shown the output of this calculation as 100 which we have now calculated using a python code with this we will end this class thank you for your attention.