 Good morning to all of you, I welcome all of you to the last day of workshop on competition fluid dynamics. In the yesterday's lecture I was in the middle of a topic which is solution of new stoke equation on a staggered grid. I had started with an introduction and completed finite volume method on a staggered grid where I had shown you that let me start the lecture on whiteboard and because I have as I said that in my last lecture that this is a this is probably the most difficult and challenging part of teaching CFD. And I had discussed the issue of pressure velocity decoupling and the role of staggering of grid on resolving this problem of pressure velocity decoupling in the lecture slides. But let me do the same thing slowly in whiteboard so that each one of you are able to understand and appreciate this issue of pressure velocity decoupling and how staggered grid is a remedy to that. Before that I would like to make some general statement saying that right now you may be thinking that we are discussing staggered grid so when staggered grid is used and when collocated grid is used. Let me tell you that this staggered grid was the first remedy to the problem of pressure velocity decoupling. However, when people as I said that the way if you go to the development of CFD initially people were not able to solve the new stoke equation due to this pressure velocity decoupling problem. Then the first remedy which came up was a staggered grid. Now when this staggered grid was good as far as standard coordinate system or problems are concerned like if some CFD problem is in standard cylindrical, spherical or Cartesian coordinate system then this formulation of staggered grid works well. However, when you move into more challenging industrial problem in those situations the geometry becomes much more complicated and you have to go for complex geometry formulations using finite volume method. Now when you go for complex geometry formulation it is very difficult to or it becomes almost impossible to use this staggered grid approach. So with that problem people came up with a second remedy to the pressure velocity decoupling which is called as a collocated grid method. Due to limitation of time I am not taking that topic here. However, we teach the collocated grid method also in our CFD course here, but I would just like to mention there the grid point same grid point is in collocated grid same grid point is used for pressure, velocity, temperature. So that is for all the variables which is sitting at the centeroids of the control volume. However, to resolve the pressure velocity decoupling problem instead of linear interpolation to calculate the mass flux at the face center of the control volume there is another type of interpolation which is done which is called as momentum interpolation. So let me start with explaining what is pressure velocity decoupling in a whiteboard. So when you discretize the continuity equation you get terms like if you take a control volume let us say this is the east face center, this is the west face center. So when you want to calculate the total mass transfer across the control volume let us say in the x direction it becomes rho u e delta surface area in the east face is delta y minus rho u w into delta y this is the total mass transfer in the x direction. Similarly, you can do in y direction I am showing you here only one of the direction. So what you get is rho delta y u e minus u w now if you do linear interpolation where this is the node p this is the node capital E and this is the node capital W then you get u capital E plus u capital P divided by 2 minus u capital P plus u capital W divided by 2. Now when you do the subtraction this u p cancels down and you end up with an expression in mass conservation equation which is this. Now what you realize is that to calculate to do the mass balance your equation is such similarly when you do in y direction you get rho delta x v capital N minus v capital S divided by 2. So what you realize is that to let us say this is the node N and this is the node S. So to do mass balance the nodes which are involved are this the u velocity of the east and west neighbor and v velocity of the north and south neighbor. Similarly when you can when you apply the momentum conservation so this is for mass conservation similarly if you apply a momentum conservation then at the same point so this is the u velocity grade point. Now when you want to calculate total pressure force you get p small e here delta y p small w into delta y and the net pressure force in the x direction is p small w minus p small e into delta y. Now when you do the linear interpolation it becomes p capital W plus p capital P divided by 2 minus p capital E plus p capital P divided by 2 into delta y and this p p again cancels down. So the pressure force which is acting for when you sorry this will be capital W and this will be capital E. So the pressure force which acts in the x direction comes from not node p but neighboring nodes capital E and capital W. Similarly when you calculate the pressure force in the y direction when you apply the law of conservation of y momentum the pressure of the north and south comes into the expression. So the point is to calculate the velocity at a particular grid point the pressure of neighboring grid points are used. So this is what is called as pressure velocity decoupling because to calculate pressure at a particular sorry velocity at this node the pressure of the same node is not being used. So there is a breakage in communication between the pressure and velocity of the same grid points which results in what we call as pressure velocity decoupling. And why this is happening this is happening because of this interpolations. I will write down this is the first interpolation here I am showing this is the second one this is the third one and this is the fourth one. So first and second are the interpolations for the normal velocity at the phase center. Third and fourth are interpolations for the pressure at the phase center. So the here the issue is this decoupling is occurring due to this interpolation. So if you want to resolve this problem of pressure velocity decoupling we have to somehow do some other interpolation which is done in a collocated grid method called as momentum interpolation for normal velocity. Otherwise the best thing is that we avoid interpolations and that is what is being done in a staggered grid. How we avoid interpolations I will show you in the next page is the way to avoid interpolation is I say that this is the grid point for pressure this is the grid point for u velocity this is the grid point for v velocity. Now this if this is p subscript p then there are two vertical phase centers. So what I am doing is here is note that if I want to avoid this interpolation which is of normal velocity at the control volume of pressure I am prescribing proposing the u velocity on the vertical phase center and v velocity on the horizontal phase center. Now the convention which we are following the running indices of this corresponds to the running indices of the cell center at which it is lying on the east phase. So this will become up similarly this become vp this become uw and this become vs. So when you apply mass conservation now you can see it becomes rho into delta bar up minus uw if you want to calculate net out flux of mass plus rho into delta x vp minus vw. So in this you see that we are not doing any interpolation for the normal velocity. So we have eliminated that now when you go to the note that this is the control volume for mass conservation. Now control volume for x momentum conservation will be a control volume where now this should be sitting at the centroid of the control volume. So now if you want this square to be sitting at the centroid of the control volume it is staggered in the x direction by delta x by 2. Now this becomes up from here you can see this is pp then this becomes pe. Now v velocity sitting at the top this is vp this will be vs this will be vse this will be ve. E means i plus 1 s means j minus 1 se means s is j minus 1 e will be i plus 1. So now if you see this control volume for u velocity in which you apply law of conservation of x momentum now you see that what will be the net pressure force in the x direction p e divided by delta y this is for x momentum. Similarly if you generate control volume for y momentum then it will be let us say this one. So for y momentum control volume you get pressure now you can see for this control volume pp is here p capital N will be here and the net pressure force for y momentum will be p capital P minus p capital N into delta x. So now you see that this values at the phase centers we are avoiding any interpolations which we are doing in previous slides. Now with this it is very clear step by step how the interpolations for pressures at the phase center in case of momentum equation and in linear interpolations of normal velocity in the mass conservation is avoided. Now let me go to a case when I said that this is good in standard coordinate system but this is very difficult in a complex geometry problem. In a standard coordinate system the control volume we have let us say this is the control volume in Cartesian. Let us say this is the control volume in cylindrical. In a standard coordinate system you can easily find the locate this because there is a coordinate along which you have let us say in this case green square are along horizontal direction. In this case the green or the square here not green there is no color here this square here are along the angular directions. Now this red inverted this inverted triangles which corresponds to the v velocity here in Cartesian they are in y direction whereas in this case they are in radial directions. This is a r phi case but when you have a complex geometry what is meant by complex geometry I will discuss. In that in a complex geometry what happens is the control volumes are not aligned along the any standard coordinate directions you may have coordinates like may have control volumes like this. Let us suppose this is P this is capital E this is capital W. Now in this you can see that this is x coordinate direction this is y coordinate direction this is phi coordinate direction this is on a radial coordinate direction. Now here when you join let us say this points what you realize is that they are not in any standard coordinate direction and in this case if I want to define velocities here which velocity I should define that is the question because the norm actually we want we do when we want to define grid point it is for normal velocity. So the normal velocity on this face is on this direction if you go to this face it is in some other direction when you go here the normal velocity is let us say here when you go here normal velocity is in this direction. So face to face the normal velocity direction is varying in a standard coordinate system note that there are only two directions of the normal velocity horizontal and vertical in carton and angular and radial in cylindrical. So when you want to use this staggered grid formulation for a complex geometry problem then you cannot use this staggered grid. So then people came up with a formulation which is called as a collocated grid method where instead of using linear momentum sorry linear interpolation for normal velocity in the mass conservation they used an interpolation which is called as momentum interpolation. I had written one chapter on finite volume method in a book computational fluid flow and heat transfer edited by professors Murli there from IIT Kanpur and professor Sundar Rajan from IIT Madras and in that chapter that is one chapter which is added in the revised edition and in that chapter this collocated grid method for complex geometry is discussed we teach this topic complex geometry formulation using collocated grid method to our students here however due to limitation of time I am unable to cover that topic here. So I had mentioned that if you use a fully explicit method then to get u velocity you can take x momentum equation and when you take x momentum equation and use fully explicit method where you take not only u velocities and v velocities of previous time level you take pressures of the previous time level also but whenever you take pressures of the previous time level you can solve u momentum equation explicitly such that if there are 25 yellow circle 25 interior points you will get 25 equation and in each equation there will be only one unknown. So you can solve for u velocity and this u velocity will be we are using a staggered grid formulation. So in a staggered grid formulation the u velocity when you go for mass conservation is the normal velocity for the pressure control volume. Then if you use bi-momentum equation and obtain the v velocity explicitly then you get that v velocity x is the normal velocity for pressure control volume. So let us go back to this window and I will show you an animation. So using the pressures of the previous time level velocities of the previous time level by solving x momentum equation let us suppose you have obtained the u velocity v velocity of the next time level but they are not obeying mass balance. As they are not obeying mass balance this velocities are called as predicted velocity. Note that here the methodology which is being discussed it is like predictor-corrector method and it is iterative image. Why it becomes iterative? So I had mentioned there is a issue that we do not have explicit equation for pressure. We do not have any conservation law where we have del p by del t. There is no explicit term for temporal variation of pressure and we want to capture temporal variation of pressure. Let us say to create a movie of how pressure changes with respect to time. So when you go to solve mass conservation equation continuity equation in fact we do not have any pressure term and we want to use del u by del x plus del v by del y for an equation as an equation for pressure. So mathematically it is converted into an equation for pressure but here I am showing you, here I am discussing the same formulation philosophically and later on so I am starting with some physical discussion and then I will go to the mathematical expression. So using the predicted velocities which acts less, acts as normal velocities on the face center of the pressure control volume we end up with this velocity which are predicted velocity they are not obeying the mass balance. So we need to tune the pressure. So let us suppose you have a knob, so right now more mass is going out as less mass is coming in. So to do a mass balance you want to increase the inflow. If you want to increase the inflow and if you have a control knob you would, what you would like to do? Increase the pressure, what is the effect of increasing or decreasing the pressure as far as inflow is concerned? If you are increasing the pressure of this control volume it will resist the inflow. If you decrease the pressure of this control volume it will create a suction. So more mass will come in. So in this case what do you want? You want more mass to come in. So let us suppose you decrease the pressure by minus 0.5. So one atmosphere is the pressure of the previous time level you have a pressure correction which is minus 0.5. So p plus p prime 1 minus 0.5 which is 0.5. This completes your first iteration then you again predict. So with this change in the pressure you again predict the mass flux and let us suppose with then updated pressure of the cell center you got updated or let us say corrected mass fluxes, mass flow rates not mass flux sorry. This is the second iteration. Now what is happening? More mass is coming in. Now we should create suction. You should decrease the pressure or you should increase the pressure. Now you should increase the pressure of this control volume. So let us say the tuning which you have done is 1.5. So you are increasing the pressure by an amount 1.5 atmosphere. Then the total pressure changes from 0.5 plus 1.5 which is 2 atmosphere. Then what will happen? Third picture will come which will be the corrected value. So you start with the prediction then tune the pressure you get a corrected value of mass fluxes. From that mass flow rates you calculate the mass imbalance. From the mass imbalance you tune the pressure you have to see whether the mass imbalance is positive or negative. So this way it goes on and on. And this is the philosophy with which we convert the continuity equation as an equation for pressure or pressure correction. I will show you symbolically about the equations. So the first equation is in a functional form is prediction of mass fluxes as a function of velocities and pressures in the previous time level. And what is this mass source? Mass source is net mass going out. So in the first picture what is the net mass going out? 4 minus 2 kg per second which is 2 kg per second. In the second iteration or second picture what is the net mass out flux? 2 kg minus 6 kg minus 4 kg. So mass source is net mass out flow. Note that when the mass out flow is positive when the mass source is positive pressure correction is negative. What is the pressure correction in the first figure? I have shown through that yellow arrow minus 0.5. What is the mass source in the first picture? 2 kg per second. So when it is plus 2 kg per second pressure correction is minus 0.5. Now let us go to the second picture. Here mass imbalance is mass source is 2 minus 6 minus 4. So when mass source is negative pressure correction is plus 1.5 atmosphere. So note the relationship and later on I will show you mathematical expression which follows this relationship. So meet equation which I am showing you. First I had shown you the first for mass flux prediction. Here I am showing in a symbolic form. Later I will show you detailed expression. Second to calculate the mass source. So the prediction of mass flux which I had just written as in the first picture 1 kg, 2 kg all these are calculated by one equation. The mass imbalance there is an equation. Right now I have discussed philosophically that let us say you have a knob by which you can tune the pressure and by certain amount which we call as pressure correction this tuning is done by an equation. And that equation philosophically you can appreciate and understand that it should be function of mass imbalance. So later on I will write out symbolic here. Later on I will show you detailed expression. Once you tune the pressure let us say on the top phase in the first picture it is 2 kg per second. In the bottom phase in the bottom picture at the top phase it is 1 kg per second. So there is a mass flux correction minus 1. And on the left surface between the two picture the mass flux correction is plus 2 kg per second. So this mass flux correction which is plus 2 or minus 1 shown in the fourth step in this slide. I will show you that this mass flux correction I will show you equation as a function of pressure correction. Because this correction occurs when you tune the pressure. Then you as I said that you start with a prediction then do a correction then get an updated value of. What is updated value? Like on the left phase on the first picture it is 1 kg and 3 kg per second in the second picture is an updated value. And the correction is plus 2 kg per second. So you get an updated value of the mass fluxes as in the second picture. And with this you get a new mass source that is step 2 new pressure correction new mass flux correction as a function of pressure correction. So this is this way the picture changes always. So it goes in an iterative mode. This is right now shown philosophically. Let us go to the slide and discuss the detailed expressions. Right now here the equations are symbolic in nature and I have discussed with figures and animation the philosophy. Now let us discuss the equations. Here as I am not giving you the lecture notes so whatever I have discussed that is written in detail. I will not discuss this because it is the same thing which I had discussed just now. Now let us go to the expressions. Here let me go I would like to mention that the continuity equation should always be at a time level and plus 1. So whether you are following an explicit method let me tell you in neo stokes equation there is no fully explicit method. Why? Because it does not work. In a neo stokes solver either you have what is called as semi explicit method or semi implicit method method. What you may know more is that semi implicit pressure linked equations which is called as simple algorithm, simpler algorithm, simple algorithm. As a semi explicit method is more towards explicit but there is some reason one term which is made implicit so that is why it is called as semi explicit. Semi implicit method is more close towards implicit. The more will be clear as this lecture proceeds but in both the methods note that the pressure is taken implicitly because we want to make sure that the velocities which we have predicted for the next time level it should obey the mass balance. So if we say that we have predicted from x momentum the u velocities and v velocities from the y momentum for the next time step then it should obey this mass imbalance for the next time. Now in an explicit method as I said why it is semi explicit because pressure is taken implicitly. Note that advection and diffusion are taken explicit in a semi implicit manner even this advection and diffusion are taken at a new time level. So the difference between the original proposition in semi explicit and semi implicit is in a semi explicit advection as well as diffusion terms are explicit. In a semi implicit these terms are implicit that is they are calculated from new time level values. So we come up with an original proposition and follow up predictor-character approach which I had discussed philosophically in animation through an animation philosophy of pressure correction. So this is the original proposition but this has a starting trouble. What is the starting trouble? Which term is causing the starting problem? Whether this term will given any starting problem? This is the rate at which x momentum is lost by the fluid in the control volume using the velocities of the previous time level which we already know. So this is there is no problem in calculating this. There is no problem in calculating total viscous force acting in the x direction using the velocity from the previous time level. Now when you go to this term do you know the pressure of the next time level? No. So this is the term which causes the starting trouble. So how we can overcome this starting trouble? What we know pressure? Pressure what we know is pressures of the previous time level which is this but if you use this you get a velocity which you have to check whether it is obeying the mass balance or not. So the pressure of the time it is not. So this velocity is called as predicted velocity. Now here I have shown you through arrows. The difference of the velocity actual velocity of the or updated velocities of the next iteration un plus 1 minus u star actually let me tell you that this un plus although I have written as un plus 1 but actually this becomes un plus 1 only if you get a picture where mass balance is obeyed because it is not that you do one calculation for prediction and one calculation for correction and you are done. You are not done in first iteration. You have to keep doing it. So actually this un plus 1 is an upgraded value of u star unless you get a u star value in a particular iteration at which the mass balance is obeyed then that u star becomes the un plus 1 which is the velocity u velocity for the next time. However the difference of the u velocity between two consecutive iterations is what is called as velocity correction. In my previous slide as I said that the mass flux is changed from 1 kg per second to 3 kg per second on the left wall. So that 1 kg per second is let us say u star. 3 kg per second is not un plus 1 it is an updated value of u star but symbolically here we are writing as un plus 1 and the difference between the two which is 2 kg per second is u prime. I can show you this is u star divided by rho into delta y because it is a mass flow rate. This is updated value of u star this is not un plus 1. So this is basically 3 divided by rho into delta y. This is u star this is an updated value and what is the correction that difference between the two. So this is the velocity correction but if you want to convert into mass flux correction you have to multiply by density and surface area which is delta y and that mass flux correction was 2 kg per second on that left wall. Here again pn plus 1 is an updated value of pressure. I will go back and show you this you started with this pressure of the old time level and you got an updated value of pressure by this pressure correction. Now this updated value of pressure becomes the correct pressure for the next time step that is pn plus 1 if and only if in the next picture it is obeying the mass value. Otherwise this is an updated value and the difference between this and this is pressure correction which is shown here. So this note that I am showing you which equation right now it is velocity correction is a function of pressure correction but you just need to multiply by density and area to convert into mass flux correction as a function of pressure correction. Now let me go back and show where is symbolically represent mass flux correction as a function of pressure correction where is that symbolically represented in which step, third step mass flux correction as a function of pressure correction. Here I had shown symbolically here I am showing you mathematically it is just that you have to multiply by density and surface area but this is the mass flux correction in the x direction for vertical faces. Then this is the slide for horizontal faces same thing is done the only difference is now we are doing this for y momentum equation note that for x momentum equation the control volume which is used is different from what is used in the y momentum equation. So in this slide we calculate the mass flux correction in the vertical direction. So that way we have completed the equation for mass flux correction as a function of pressure correction. So if you tune the pressure it causes this must change in the mass flow rate at the faces by tuning the pressure the consequence of the tuning the pressure. If you create suction mass flux correction will be positive if you increase the pressure mass flux correction will be negative. So as there is no specific question for calculation of pressure there are actually two class of methods for calculation of pressure one which is called as pressure correction method which is most commonly used simple method simpler simple r pressure correction method. There is another method where instead of using pressure correction pressure is used as a variable called as projection method. So how what is the expression between the velocity predictions and the correction these are the expressions. But note that there is a statement which I had written at the bottom which I had emphasized earlier in earlier slides. The pressure and velocity correction is continued iteratively till a divergence free velocity field is obtained then only the right hand side of the above equation becomes the velocity and the pressure for the next time step otherwise it is just an updated value. So this does not happen in one iteration with that we had we had done with the mathematical expression for mass flux correction as a function of pressure correction. We need other expression also for that so we have already completed the discretization of x momentum y momentum. Now we are left with continuity because now what is the objective here for from x momentum and y momentum we have been successful in getting an expression for mass flux correction as a function of pressure correction. Now we want to get an expression for pressure correction which is pressure tuning as a function of mass imbalance correct which is let me go back to the slide what I am saying is that from x and y momentum equation we have been successful in getting detailed expression not symbolic expression shown here for this. Now we want to get an expression for pressure correction as a function of mass imbalance. So this we have got from momentum equation this we will obtain from continuity equation this equation is very simple mass imbalance is just doing a balance. So you know this expression there is nothing difficult about it this I had already shown this is you predict u velocity using the u velocity of the neck you predict u velocity for the next iteration using the velocities and pressures for the of the previous time step when it is a first iteration. So and this mass flux correction I had already shown you the expression. So there are in fact two main expression mass flux correction as a function of pressure correction which comes from momentum equation and the second pressure correction as a function of mass imbalance which comes from continuity equation mass conservation equation which I will show you now okay. So this is the discrete form of mass conservation equation note that the normal velocities are sitting at the places of the control volume in which we apply the mass conservation and we need not do any interpolations. If u p n plus 1 in previous slide for x momentum I had shown you this expression what I am showing you here is basically this is p p prime minus p e prime into delta y. So u p n plus 1 is equals to u p star plus p p prime minus p e prime into delta y. Similarly you can do for v p n plus 1 that is what I am showing you okay and this is the predicted velocity and this whole term is called as the velocity correction from which you can calculate the mass flux correction. Now the objective is that we obtain expression for u p n plus 1 v p n plus 1 which I am showing you here. Similar expression can be obtained from u w for u w n plus 1 and v s n plus 1 and we substitute into this expression and then we get this expression. Now how is the nature of this equation let me ask you what are the terms in this equation what are the terms in this equation predicted velocities and pressure correction. Now these predicted velocities are which velocities they are normal velocity they contribute the mass flux mass flow rates. So we will get pressure correction as a function of mass imbalance from this equation. However what is the nature of the equation do you find any similarity of this linear algebraic equation with heat conduction equation with volumetric heat generation okay. Now let us suppose let us look into this term do not you think that this term has analogy with k dT by dx on east phase center do not you think this is analogy with not plus k minus k dT by dx what is the discrete form of minus k dT by dx it is minus k p capital E minus t capital P divided by delta x as it is minus k when you make it plus k it becomes t at capital P minus t at capital E divided by delta x. So instead of temperature here we are having pressure correction and instead of conductivity we are having delta t actually this delta t is analogous to minus of k okay. In fact this equation although right now it is an algebraic form but it can be converted into a differential form also and it is a Poisson equation it is like a steady state conduction equation note that it is analogous to steady state conduction equation you can see there is no del by del t term here there is no unsteady term time derivative is not there there is no delta t in this equation there is delta t is there but the delta t is analogous to minus of k. So this is analogous to if you convert into a differential form it is basically del square p is equals to some source term which is of Poisson equation and it is analogous to steady state heat conduction with volumetric heat generation now what is that volumetric heat generation analogy or source here is it is mass imbalance because this up star note that you have to multiply by rho delta y and what is up star into rho delta y this gives you I am showing you here all this star u star and v star it can be expressed into this now what is this term this is the discrete form of the continuity equation using the predicted velocity and in fact this is the negative of mass imbalance note that it is negative it is negative of mass imbalance and I would like to remind you that I told you that we should have expression pressure correction as negative of mass imbalance let me go back and show to you pressure correction is negative of mass imbalance if the mass imbalance is positive 4 minus 2 which is plus 2 pressure correction is negative when the mass imbalance is here 6 minus 2 sorry 2 minus 6 because we are talking of outflow minus inflow 2 minus 6 which is minus 4 then the pressure correction is plus 1.5. So the pressure correction mass imbalance should be related with the negative sign exactly similar relationship is shown here pressure correction as a function of negative of mass imbalance this is mass imbalance coming from the predicted normal velocities. So with this I had shown you all the equations I started with the philosophical discussions and showed you the equations showed you the equations for velocity correction mass flux correction as a function of pressure correction as a function of mass imbalance. So this is called a semi explicit method this is the formulation for the semi explicit method where I had given you shown you expressions. So this is what is called as formulation now let us go to the formulation of semi implicit before moving on to solution methodology that this equation because we do not have more we have more than one equation what should be the methodology what should be the procedure to solve this equation one by one. So before discussing solution methodology let us discuss what is semi implicit method what I am discussing here is basically simple method the pressure is taken fully implicitly as in case of explicit method and this is the original proposition to obtain u velocities and v velocities for the next time level. Now in this original proposition what is the difference between this original proposition in semi implicit method with that of semi difference is how we treat the advection term and diffusion term. Note that the diffusion term and advection term have a superscript n plus 1 which indicates that to calculate net y momentum or x momentum lost by the fluid in the control volume and net viscous forces acting in the x direction and y direction is calculated using the velocity of the next time level and you know that in implicit you get a system of algebraic equation where there are more than one unknown. So at each time step when you want to solve this equation you have to solve iteratively now here again there is a starting trouble you may say that I do not know the pressure of the next time level I do not know the velocities of the next time level. So let us use the pressures of the previous time level and let us predict the velocity and when we are predicting that velocity let us use the predicted velocity of the neighbors although they are not known because that is what happen in a implicit method because in this case note that you will get a linear algebraic equation for up star u w star u capital E star u capital N star u capital S star and in each equation there will be 5 unknowns if it is an interior point. So there are 25 equations there will be there are 25 interior points you will have 25 equations and in each equation there will be not more than 5 unknowns. So if you want to predict this velocities you have to solve this equation iteratively similarly for v velocities this is in semi implicit method note in that in semi explicit method this advection and diffusion were calculated using the velocities of previous time level. So in this in a semi explicit method you get there are 25 grid points you get 25 equation for u velocities 25 equations for v velocities and in each equation there is only one unknown. Now if you want to get an equation for mass flux correction as a function of pressure correction that is what I had shown earlier for from the momentum equation. So you just subtract the original proposition with this predictor step you know what is predictor corrector method prediction plus correction is equal to the final solution. So this is the final this is the original proposition this is equals to the prediction plus correction if you want to get an expression for correction you have to subtract this 2 equation. So when you subtract this 2 equation note that in pure explicit sorry semi explicit method when you subtract the original proposition let me show you note that in semi explicit method when you subtract this 2 to get a corrector step this an is same in both original proposition as well as prediction. So it cancels down completely diffusion also cancels down but in a semi implicit method when you subtract original proposition and the predicted equation for predicted step whether this 2 will cancel down this is a p at n plus 1 this is a p at star which uses predicted velocity in a coupled sense. So when you subtract this 2 this is non-zero and in fact this involves velocity corrections of the neighbor and you get an algebraic equation like this actually in CFD we want to get an algebraic equation and we want to get an algebraic equation which is sparse in nature. If there are million points you will get million equations but in an implicit method in each equation there will be not more than 5 unknowns. So million minus 5 in each row will be 0 so there are many zeros so that is what is called as sparse linear algebraic equation. If you want to solve this equation it will be like million equation having million unknowns for each point there will be almost no zeros. So your coefficient matrix will be completely filled with numbers till this step let me tell you all the steps which we are doing is fully implicit. So why it is called a semi implicit because at this step we neglect let us say it is like z is equals to x plus y and then we are saying z is approximately equal to y. So we are neglecting the this term let us say x so when you neglect the velocity correction of the neighbors and why we are doing there is a purpose because in CFD we want to have sparse matrix and we solve iteratively. Note that at this step it becomes semi implicit. Before this step it is fully implicit this we did for u velocity we have to do this for v velocity control volume also shown here. Here again we consider the velocity corrections of the neighbors as negligible note that there is an equal sign here and there is an approximate sign here and at this step it becomes semi implicit. So when you substitute the so this is the equation for normal velocity correction as a function of pressure correction when you substitute into the continuity equation the way I had shown earlier you get a linear algebraic equation this is basically the discretization of Poisson equation del square p prime is equals to minus of mass imbalance. So this results into pressure correction equation pressure correction equation is an equation to obtain pressure correction is a function of mass imbalance. Here I can note the negative sign negative of mass imbalance. So this completes the formulation before I go into the solution methodology let us go into the implementation details. Now actually we have three different control volumes. Now when we want to predict u velocity so what is the idea actually if you want to know how to proceed go back to that animation which I had shown for pressure what was the first thing which came in the animation we take a control volume we use the pressures of the previous time level. Then what is the first thing which came into that animation the predicted mass. What is the second thing mass imbalance what was the third thing pressure correction is a function of mass imbalance. What is the fourth thing mass flux correction is a function of pressure correction or a weighted value of mass flux to be proceed in a solution algorithm or our procedures that way. So let me start with the first step that is prediction of u velocity. Just to connect to you we are taking an example problem which was given to you in yesterday lab session lead driven cavity flow which is just a square cavity whose three walls are stationary left right and bottom top wall is like a conveyor belt which is moving from left to right with a constant velocity u naught. So what will happen as this is moving the fluid which is close to it will try to moves along with it at the corner it cannot go out. So it has to turn this will create a clockwise vortex the flow will generate a clockwise vortex near to that corner there will be an anticlockwise vortex small anticlockwise vortex. So this is your staggered grid which I had I had shown you animation for this in the beginning of this topic. So I will not show you an animation here what I will do is with this slide I want to show you the how to predict the u velocity for the u in the u control volume. Same thing I will show you to an animation it will be your feel better feel of what is happening. So let us come to the animation this is the square domain which we are using for lead driven cavity flow problem this is the lead this is the velocity with which lead is moving. Now you see that this vertical lines are not equi-space near to the boundary why because here I want to show you control volume for u velocity and you know that the control volume for u velocity are staggered in horizontal direction by delta x by 2 whereas the horizontal line will be equi-space because u velocity is not staggered in the y direction and where are your u control volumes this are the running indices note that we are not using different running indices for pressure velocities and we have used a convention where we are making sure that we use a single running indices. In fact in the formulations which I had given you in the slides although there are certain subscript in denoting the neighbors or the cell center but when you program you have to convert into appropriate ij. So this is the interior grid point for u control volume how many are they 4 in the x direction 5 in the y direction this are the boundary grid points left bottom right and top note that on the left and right boundary there is a half control volume but do you need to solve in that half control volume do you apply a conservation laws in that half control volume what is the grid point coming in that half control volume boundary grid point and I had emphasized in my earlier lecture that for boundary grid points we do not apply conservation laws what we apply is boundary condition we apply conservation laws in the full control volume and this is the typical computational stencils this is the interior control volume I am also showing you the running indices for the interior control volume. Now what I will show you is that as you know when you apply law of conservation of mass in this control volume you will have to calculate the fluxes mass flux x momentum flux viscous stresses at the faces. Here I am showing you the face centers at which you will calculate the fluxes in the x direction and I am showing you the running indices in the i direction which is which is different from the interior node but in the j direction it is same. So I am not repeating that note that this yellow squares are 5 in the x direction and 5 in the y direction. Now what are this squares this are the this pink square corresponds to the points at which we calculate the fluxes in the y direction mass flux in the y direction x momentum flux in the y direction viscous stresses here we will calculate always in the x direction whether it is a vertical surface or a horizontal surface for vertical surface it will be normal stress and for horizontal surface it will be a shear stress. These are the boundary points and these are the running indices. Note that for vertical and horizontal I am showing only one of for vertical I am showing you i running indices for horizontal I am showing j. The other running indices is same as what is given for interior control and this is the control volume where we are trying to conserve the momentum. Note that I will again repeat that when you want to calculate x momentum flux in this control volume x momentum flux also involved mass flux and note that although in continuity equation mass flux grid point was sitting at the faces of the pressure control volume. But in this control in case for to calculate mass flux you have to do linear interpolation. In my previous slides I had shown you expression m e u m w u I have used a subscript u. This is the convention which we are following. So, at this points we calculate the mass flux momentum flux and viscous stresses in the x direction and on this face again we calculate the mass flux in the y direction x momentum flux in the y direction and the viscous stress which is the shear stress in the x direction. Similar thing we do for y momentum now here you see that vertical lines are equi-spaced, but horizontal lines near to the boundary there is half control volume and these are the control volumes for v velocity interior control volume. In for you it was 4 in the x direction 5 in the y direction here it is 5 in the x direction 4 in the y direction. I would like to mention here that this pictures are very good although they seem to be that there are very few points, but if you develop a program to begin with it is just that you have to use some loop. So, the time x whether you give a 7 or whether you give a 7 million the code should be able to handle, but to decide about the running indices it is better to use the 7 points and develop the code and in my formulations I had given various specific neighboring information which you can convert into appropriate i's and j's. I will be showing you the pseudo codes also this is to help you as far as possible through equations, figures, animations, pseudo codes. So, that you do not struggle and start teaching the CFD course. This is the boundary grid point in the half control volume where we do not use the law of conservation of wire momentum. We have used boundary conditions this is at the top this is on the left boundary this is on the right boundary. Next the flux calculations. So, this are the interior control volumes this are the running indices or this yellow inverted triangle is the these are the points where we calculate the flux in the x direction. This is the boundary at the boundary on this green inverted triangle we calculate the mass flux in the y direction momentum flux in the y direction discuss stresses in the y direction and the running indices is this. This is a convention which we are following as far as running indices is concerned. So, x momentum y momentum now I am showing you mass conservation. Now, the control volume here is equi-spaced now here in this case there will be no half control volume. These are the pressure grid points these are the phases where you calculate mass flux in the x direction. These are the phases where you will calculate mass flux in the y direction. Note that this mass flux is only used in the mass conservation this is not used in the momentum conservation. There we have to do linear interpolations. If you look into the temperature the temperature control volume is same as that of a scalar control volume which is pressure. So, these are the grid point for temperature these are the running indices. You calculate heat flux and enthalpy flux in the x direction at this phase center. You calculate conduction heat flux and enthalpy flux in the y direction in this red circles. So, with this I had completed the showing you animation of the different control volumes and the different phase centers which are involved.