Now in this segment, we'll move on to looking at solving the transient heat conduction problem in the hallmark five template. If we look over here, now that we're solving the transient problem, we have to plan initial conditions. And I've actually set up most of this for you. What we need to do is we loop through all of our nodes, and we check the position of each node. And, depending on the position of the node we'll define the initial temperature at that node. Pretty straightforward. You'll just need to insert those values, in some cases it actually depends on the x value of the note itself, as a function okay? All that information is given in the homework, you just need to calculate the value. Now that takes care of d naught. Storing the values of the initial values of the temperatures in D_trans. But now we will also need the initial values for v. We need v naught. And so to do that, we will solve a matrix vector system, okay. Now remember, this problem is based on this equation. We have M times V plus K times D is equal to F. F in this problem is of course zero, but in general, it wouldn't and it could be changing with respect to time. Actually, all of these objects could be changing with respect to time, potentially. So, we can also put the initial values in this equation, and if we need to solve for V naught, well first we move KD over to the right hand side This becomes our right hand side vector. And then simply do a matrix inversion to solve for v not. Now in this problem the mass matrix and the stiffness matrix don't change with respect to time, which is why I can get away with calling the symbol system only once. Okay, so let's look in the code at how I'm doing this matrix inversion, and how I'm doing these matrix vector operations. Okay, I've done all of these for you, because you'll be doing several of these on your own within the solve trans function. Okay, so first, I set RHS equal to K times D trans. And to do that I used this dot vmult function. Now with sparse matrices it will do a matrix vector multiplication. I pass in RHS and D trans. And it's acting on my sparse matrix K. But what it does is it stores within RHS the values of K times D_trans, okay? In the next step, I just take RHS times equals minus 1. Essentially I'm just changing the sign on all the components in my RHS vector. Okay, so this point RHS is equal to minus K times D trans or minus K times D naught. And then finally, even though F is equal to 0 and in this case I've included it for generality, I do RHS add 1 times F When is the scalar coefficient to F. And so that takes RHS and adds on to it one times F, which is, now, gives us a right-hand side vector of F minus K times D trans, or F minus K times D naught. Okay, so now my right hand side vector is set up. In this case, I'm only doing a matrix inversion on M, so I don't have to do any matrix operations there. But I did, up here, copy M into the system matrix, okay? So that's the same as saying system matrix equals M. Now that all these data objects are set up, I can apply my boundary values. Now note that, since I'm solving for V trans in this case, I will apply boundary values of V. The matrix that I'm inverting is system matrix, the solution is v trans through an outside vector is RHS. I create this sparse direct umf pack A, which copies over from system matrix. And then I do my matrix inversion. Now it may be a little confusing. I do A dot vmult here, but you notice it's doing a matrix inversion and then multiplying it by the right hand side vector. It's just the way the functions were named and so you just need to remember that with the sparse direct GMF pack, sparse matrix dot vmult does the matrix inversion. However, when you're dealing with a sparse matrix and a vector dot v mult, just as a simple matrix vector multiplication without any matrix inversion. Okay? All right, so with that background, let's scroll down to the solve trans function. Okay so the first thing we do in solve trans, is apply the initial conditions. So that's the function that we are just looking at. And the next is to define delta_t. I have delta_t set to equal one and that will give you a stable algorithm for all alphas, alphas = 0,1/2, and 1. You're free to change that on your own however, when you turn in your assignment, leave it at delta t = 1, okay. I declared this vector D_tilde, which you remember we use as we're updating. From DN to DN+1. And then we go into our time loop. Now I'm going to point out here that we are using the V method here. Okay. And so we'll look over those equations as we go through this code. All right, so the first step is to find detailed. Let me write that out for us here. So detailed, n plus 1, again this is a global vector, is equal to d at time n plus delta t times 1 minus alpha v n, right. These are equations taken from the lectures so you can look it up there. So this is step one. Step two is to use d n plus one tilde to find v n plus one. And we have this system maybe, equations here. We have M plus alpha delta tK, and I'm not going to re-derive these during the lectures, but here's the result. Sorry, we haven't quite inverted that yet. going to delete that. Okay, so M + alpha delta t K times vn + 1 is equal to Fn + 1, which in our case is just 0. Minus k times d tilde n plus 1. Which of course you would then to solve for vn plus 1 you would do a matrix inversion. So we take m plus alpha delta tk, inverse times fn plus 1 minus k. Detailed m + 1. This will be referring to as our system matrix in the code. And this we'll set up as the right hand side, the RHS factor. Okay, once we have the n plus one, step three is to find d n plus one, which is d tilde n plus one, plus alpha delta t, v n plus one. Okay? Note in the code dn at this point is d trans, vn is also v trans. Here we update v trans to be Vn+1, and at this step we're updating d trans Sorry that's a capital V in our code. Updating D-trans to be dn plus 1. Okay. So we'll be following these three steps in the code itself. Okay, so let's go back to the code. You'll notice I've given a list here. Of how to perform these matrix vector operations. A lot of these we've already seen above in the apply initial conditions function. You can use this same .add function with the system matrix, like we did with the vector. Again dot v mult will do a matrix factor multiplication. If you want to set m equal to or system matrix equal to a vector or to a matrix excuse me, you would use copy from. However with a vector you can just say RHS equals f for example. We already saw this. I'm operator times equals with a vector to change the sign. Okay, so that should be enough for you to be able to do all of these matrix vector operations. Okay, so the first step again was to find D tilde. And again, at this point D trans it does hold the values of D so then nv trans holds the value of V sub N. Okay so once you found D trans you will then create system matrix, as we wrote on the board and RHS. Again as we wrote on the board. Once you've done that, we apply boundary values of V, again because V trans or VM + 1 is the solution vector we're solving for. So we apply boundary values of V using the system matrix that you've defined and the RHS factor. Okay, and that matrix inversion happens here with these three lines. Again using the sparse strict plume of pact matrix, A. Okay, once you've done that matrix inversion and V trans is now equal to V n+1, you can use D tilde and V m+1, V trans, to solve for D m+1. Okay and there you have the solution for this time step. The last thing to do is to output results. Now I have it set to output results every 1,000 time steps. You'll notice I'm using this modular operator. So every time the current time step index, divided by a 1,000 gives a remainder of 0. Then we'll output the transient results as a VTK file. And also store at the L2 norm of the results in this L2 norm results. Or, sorry, I'll calculate the L2 norm of the difference between the transient solution at this current time step and the steady state solution. I'll store that within this vector here, okay? Alright. So the last step, after outputting results, or the last function I should say, after the output results functions, is to calculate the L2 norm. This is going to be very similar to the L2 norm that you calculated in Homework Two, except of course you can now use dl2 basis functions. You will still need to interpolate to find U steady and U transient in exactly the same way you interpolated U H, you will find that element solution in the torched point. And you will construct the intergrand of you know, the square of u-trans minus U study. Okay, and then we we return the square route of that value to actually calculate the l2 norm itself. Alright, so that should be pretty straight forward with what you've seen before, and that concludes our discussion of the template for homework five.
