 We have been discussing issues related to analysis of equations of motion, we are discussing alternative damping models, so we will continue with that discussion and then we will take up two special topics in this lecture, one is on dynamic stiffness matrices and other one is on transfer matrices. So we have been, we have classified damping models as viscous and structural and each one is either classical or non-classical and we have now developed the expressions for frequency response function for all these alternative damping models. There is one way of doing that is to simply invert the so called dynamic stiffness matrix that is if you have this equation of motion it can be expressed in frequency domain and then we get the matrix of receptance functions as inverse of this matrix, this can be evaluated for every omega and similarly the system is structurally damped, this is the expression for the frequency response function matrix, so if we are adopting this approach of directly inverting the dynamic stiffness matrix then there is no need to worry too much about handling damping matrix, so this procedure is universally valid no matter whether damping is viscous or structural or classical or non-classical etc. But we, this is a brute force method it doesn't provide too many insights it can be computationally very demanding, so we would like to develop a kind of an uncoupled approach using natural coordinates and if we have to do that then issues of modeling damping gets closely associated with how you compute frequency response function and we have gone through that exercise and we have derived the frequency response function for the different types of damping models, this is for classical viscous damping, this is for classical structural damping, and this is for non-classical viscous damping, and this is for non-classical structural damping, so in each of this case there is a different interpretation for the natural coordinates, but the common feature is that there is a orthogonality relation satisfied by the structural matrices, I mean the model matrices with respect to the structural matrices which enable uncoupling of equations of motion. Now we pose two questions in our discussion on non-proportional dam systems, the first one was what is the mathematical framework to uncouple the equation of motion, this we have now answered. Now there is another question that needs some discussion, this question is this, are there any simplifications possible so that the damping remains classical and yet at the same time we take into account the fact that the structure is made up of subsystems with different material, so how do we approach this? So we consider a system made up of NS number of subsystems, we will assume that each of this subsystem is made up of different materials and there are therefore damping matrix itself can be envisaged as being made up of NS number of contributions from each one of these subsystems. Now CI here is a contribution to the structural damping matrix from the ith subsystem, now what we will do is let us specify damping ratio as a function of the material, in fact this is what is done in some design practices, for example the ASME section 3 boiler and pressure vessel code specify the following damping ratios in the context of earthquake response analysis of industrial structures, if the structure is equipment it is 2 percent damping, piping system 5 percent damping, welded steel structure 2 percent damping and so on and so forth, so this way of describing damping that is modal ratios as a function of material is an alternative way of specifying damping, now. So let us do, try to see how we can handle these details, so again let C be made up of this CI's as shown here and let Phi be the undamped modal matrix so that Phi transpose M Phi is identity and Phi transpose K Phi is the matrix of diagonal matrix of eigenvalues. Now each of the subsystem we assume that damping can be represented as a proportional damping, Raleigh's proportional damping modal for each of the subsystem, so CI is alpha i MI plus beta i KI, so C therefore can be written in this form. Now if you now consider Phi transpose C Phi we can obtain we use this and we will be able to write this but it is important to note that this Phi transpose MI Phi will not be diagonal because MI is simply a contribution to, CI is a contribution to the total damping matrix from the i-th subsystem, similarly MI and KI are contributions to the respective mass and stiffness matrices from the i-th subsystem to the corresponding global mass and stiffness matrices respectively. So let's consider the following, let, we will assume that beta i's are all 0 that is the first case that means for each of the subsystem damping is mass proportional, now if you recall if C were to be alpha M in Raleigh's proportional damping modal the damping ratio for the j-th mode is given by alpha by 2 omega j, now alpha is the only free parameter, now let us consider Phi transpose C Phi this I will write it as 2 eta j omega j, I am assuming that this is diagonal so it is 2 eta j omega j and this itself is made up of alpha i Phi j transpose MI Phi j, now if I divide by 2 omega j I will get eta j as this ratio, now I can take the omega j to the numerator and rewrite this with understanding that Phi j transpose M Phi j is 1 because Phi transpose M Phi is identity matrix this is 1, but we can see that the numerator I am getting Phi j transpose MI Phi j, so this can be, this has an interpretation of being the contribution to the total kinetic energy in the j-th mode by the i-th subsystem and similarly this is denominator can be interpreted as total kinetic energy in the j-th mode, okay, and this eta ij is the ratio, so eta j now therefore I will be able to express in terms of kinetic energies stored in different subsystems weighted by this factors eta ij and divided by the total kinetic energy in the j-th mode. If we now assume that eta ij, actually eta ij is a function of mode as well as material, now if I assume that for a given material damping is constant for all modes then what happens eta ij becomes independent of j, if it is for example if I have a system made up of say RCC and steel and soil and things like that suppose RCC is 5% steel is 3%, soil is say 10% damping ratios, so for all modes that is true, so we will assume that eta ij can be approximated as eta i, so I get an equivalent damping for the j-th mode of the global system in terms of the material damping and kinetic energy in different parts of the structure. Now the similar argument can be extended if I now take alpha i as 0 and write the damping matrix to be in terms of only stiffness matrix, so again if you recall if damping matrix in Rale's model if damping matrix is proportional to stiffness the j-th damping ratio is beta omega j by 2, okay let's go through this calculation now phi j transpose C phi j is 2 eta j omega j and this should be equal to summation beta i phi j transpose K i phi j, so I will again write expression for eta j by dividing by 2 omega j I get this, now I will rewrite this in a slightly different form as shown here and for omega j square I will write phi j transpose K phi j, okay, so I am multiplying and dividing by certain constants so that the expression mathematically is not altered but I will be able to interpret terms in terms of energies. Now so eta j in this case will be again contribution to the total potential energy in the j-th mode by the i-th subsystem that is what I will interpret this phi j transpose K phi j as, and the denominator which is phi j transpose K phi j is a total potential energy in the j-th mode, there is a factor of half of course. Now if you assume that modal ratio for i-th subsystem is independent of frequency that is mode count J then eta ij can be written as eta i, so the equivalent damping for the j-th global normal mode in terms of damping ratios for different materials is thus given by this ratio, so what you need to do is for a given mode shape you need to compute the strain energy in the system, the total strain energy will appear here in the denominator and for each subsystem made up of a particular material you should compute what is the contribution to the total strain energy, multiply that by the damping ratio as relevant to that material, so sum it and divide by this total potential energy you will get the equivalent damping for the j-th subsystem, the j-th mode for the global system. A generalization of this of course would be phi transpose C phi is taken as you know alpha i phi transpose MI phi plus beta i phi transpose K i phi. Now again if you recall if damping matrix is as per the Rayleigh's proportional damping model alpha M plus beta K eta J which is the damping ratio in the j-th mode is given in terms of alpha and beta and omega J through this expression. Now we will construct now phi J transpose C phi J which is 2 eta J omega J, now on the right hand side I will have this expression. Now eta J will become now sum of these two terms, so we cannot now interpret these terms in terms of energies, so if you want to use this model the user has to specify the ratios alpha and beta i for a different subsystem not the damping ratios for different materials, but of course it can be related to the damping ratios of the material if one does free vibration analysis of individual subsystems uncoupled from the remaining subsystems, so in this model alpha and beta i need to be specified for all subsystems, this can be accomplished in terms of known damping ratios for each material by performing a separate free vibration analysis for each of the subsystems in the uncoupled states, so this requires additional effort, okay. Now before we leave this topic we can consider one more mathematical complexity, this is what would happen if structural matrices are asymmetric, we have been dealing with even when we dealt with non-proportional damping we configured our equation of motion in such a way that all the structural matrices were symmetric. Now it may so happen that either we choose to configure our equations in a way that we don't end up handling unsymmetric matrices or alternately there may be physical problems, for example as in case of gyroscopic systems where the structural matrices could be asymmetric. Now I will not get into the genesis of asymmetry in structural matrices but I will outline what one could do to uncoupled equations of motion if one encounters asymmetric matrices. So what I will do is to illustrate this I will start with traditional equation of motion MU double dot plus CU dot plus KU equal to F of T, and for purpose of illustration we will assume that all these structural matrices MCK are all symmetric, okay, they are not asymmetric, but now what I will do is I will introduce a new coordinate system wherein in the transform coordinate system the matrices become asymmetric. So suppose if I do introduce a vector X in terms of YI and Y double I as U and U dot, so if YI dot will be U dot U dot is Y double I, similarly Y dot double I is U double dot that I can derive from the equation of motion. Now if I rearrange these terms in a set of 2N first order equations I get this equation, so I can write this as X dot equal to AX plus F of T, now this A matrix here is a 2N cross 2N matrix you can see here this is a A matrix and this is clearly not symmetric, okay, so now the question we are asking is how can we uncouple the equation of motion if one encounter this type of equation, so now let us start with the usual approach we wish to introduce a transformation X is equal to TZ such that after transformation we get a set of 2N uncoupled first order equations, now the question is how to find capital T, as we have been doing till now suppose if you start with the free vibration problem X dot equal to AX and seek the solution in the form Phi E raise to ST, this is fine because this is a set of linear differential equations with constant coefficient so exponentials are always solutions, so for some value of S this could be a solution, so Phi and S are unknown, so now I will substitute into the governing equation and this I get this eigenvalue problem, that means this type of solution is admissible provided S and Phi satisfy this equation, where A is the structural matrix relevant to the problem, so this is an algebraic eigenvalue problem, A is 2N cross 2N, so let S1, S2, S2N be the eigenvalues and associated eigenvectors be Phi 1, Phi 2, Phi 2N, now let us blindly follow the procedure that we have used earlier to deduce the orthogonal relation, so we will just follow the procedure as we have been doing and see where we get stuck if at all, so let us consider Jth and Kth eigenpairs, so the governing equations are A Phi J is SJ Phi J and A Phi K is SK Phi K, so name these equations as 1 and 2, now pre-multiply 1 by Phi K transpose, so I get Phi K transpose A Phi J, SJ Phi K transpose Phi J, now similarly pre-multiply equation 2 by Phi J transpose, so I will get this equation, now let us transpose both sides of equation 4, so I get Phi K transpose A transpose Phi J is equal to SK scalar that remains as it is Phi K transpose Phi J, now at this stage we should notice that A transpose is not A, so we now have the fact that A transpose is not A, so we get stuck at this step if we blindly follow the procedure that we have used earlier, so what we are doing earlier we will subtract equations say 3 and 5, the left hand side, right hand sides are the same but the left hand sides are not the same, so we will not get the required orthogonality relation or in other words Phi K's are not orthogonal to A, okay, so we can't use eigenvectors of A to uncouple the equations of motion, so how do we proceed? So what we do is we define B as A transpose, and we note that A and B would have the same eigenvalues but different eigenvectors, okay, so if A is asymmetric B is A transpose the eigenvalues of A and B would be the same, now consider the eigenvalue problem now B psi is equal to SI, that is now since S is the eigenvalue of A it should be eigenvalue of B also, so we can return that here, psi is the eigenvector associated with matrix B which is A transpose, now again let S1, S2, S2N be the eigenvalues and these be the eigenvectors SI1, SI2, 2N, now let us again consider Jth and Kth eigenpairs I get these two equations, that is what I do is the Jth eigenpair is SJ Phi J that is the eigenpair for the first system of equations, that is A Phi is S Phi, SK psi K is the eigenpair corresponding to the B matrix, B psi is equal to SSI, so the eigenvalue problems are A Phi J is SJ Phi J, B psi K is SK psi K, so name this as 1 and 2, pre-multiply 1 by psi K transpose and pre-multiply 2 by phi J transpose, if I do that I will get these two equations, so psi K transpose A Phi J is SJ psi K transpose Phi J, similarly Phi J transpose B psi K is SK Phi J transpose psi K, name it as 3 and 4, now let us transpose both sides of equation 4, so what I will get psi K transpose B transpose Phi J transpose, transpose of Phi J transpose is Phi J itself, so on the right hand side it is SK which is scalar psi K transpose Phi J, now B transpose is A, so therefore I can write it as psi K transpose A Phi J, right, so because B is A transpose, call this as equation 5, now if I subtract 3 and 5 I get SJ minus SK psi K transpose Phi J, so that would mean psi K transpose Phi J is 0 for K naught equal to J, and psi K transpose A Phi J is 0 for K naught equal, so this is the orthogonality relations that the eigenvectors of A and A transpose satisfy, so we need to therefore solve two eigenvalue problems, first of all each eigenvalue problem is of size 2n by 2n, not only that we have to solve two of them, so we can select the normalization constant so that psi K transpose Phi J is a chronicle delta, so if I now denote capital psi as the matrix of eigenvectors of the second system and matrix phi as a matrix of eigenvectors of the first system I have psi transpose phi is I, and psi transpose A phi is diagonal of SI, now how does this help us to uncouple the equations, so let us start X dot equal to AX plus F of T, with X of 0 is X naught, now let us make the transformation X equal to phi Z, I get this equation, now instead of pre-multiplying it by phi transpose I will pre-multiply by psi transpose, so I will get psi transpose phi Z dot is psi transpose A phi Z plus psi transpose F, now we know psi transpose phi is identity matrix, psi transpose A phi is a diagonal matrix with eigenvalues appearing on the diagonals, so I get therefore a set of uncoupled 2n uncoupled equation, psi transpose F is a generalized force, so I will get this equation for K equal to 1 to 2n, so these are first order equations which can be solved, and initial conditions so to find initial conditions I use this relation I pre-multiply by psi transpose I get this, that is after making the transformation X equal to phi Z I pre-multiply by psi transpose I get, so equations are uncoupled, so the point is therefore we can uncouple equations of motion in a very broad class of linear time invariant system problems, and that always helps us to represent solutions in a somewhat revealing way as a summation, so we can understand the behavior of the resulting transfer functions by considering contribution from each of these terms in the summation. Now we now move on to a special topic, we know now that FEM is after all an approximate method we are using now models from FEM that is the models resulting from application of finite element models in calculation of frequency response functions. Now there are approximations in this, what are the approximation? There is special discretization and errors due to interpolating field variables within an element using polynomials, okay, that given domain is represented approximately through a mesh that we adopt for making the finite element model, and within an element we are using the nodal coordinates and interpolating the field variables using interpolation function, so both contribute to errors, and if you are using a model superposition method to derive the frequency response function there is questions about errors resulting from model truncation. So this would mean that a finite element model has limits on spatio-temporal scales within which the solutions remain acceptable, right, a given finite element model cannot be used to resolve you know the field variables in space, there is a limit on that and in frequency or in time. So now this leads us to ask the question can frequency response function of skeletal structures be evaluated exactly, that is within the framework of beam theory used, if so these solutions can serve as benchmarks to validate approximate solutions. Now indeed such approach is possible and that is the issue that we will be briefly visiting that takes us to discussion of dynamic stiffness and transfer matrix methods. So what we will do is we will start by considering an actually vibrating rod, in this discussion the forcing functions are all harmonic and we are focusing on steady state behavior, that means system is damped and their steady state exist and we are focusing on steady state behavior. So now let us consider an actually vibrating rod with field variable U of X, T, A, E, M and L are the properties axial rigidity mass per unit length and the span of the rod, and let's assume that there are harmonic actions P1E rise to omega T and P2E rise to omega T resulting in displacements delta 1E rise to omega T and delta 2E rise to omega T, so what is the condition under which all this can coexist, that's a question. So we can write now the equation for the axial vibration and see what we get, so we will assume that the rod is homogeneous, so the governing equation will be AE dou square U by dou X square plus a structural damping term plus a viscous damping term which are proportional to the strains that is equal to the inertial term MU double dot plus a viscous damping term and a structural damping term which are proportional to velocities. So these two damping terms are proportional to structural velocities and these two damping terms are proportional to structural strain rates, okay, strain rate in time. Now what are the boundary conditions? At X equal to 0 this is a displacement, delta 1E rise to omega T, at X equal to L this is a displacement, and we are assuming that these coexist therefore the other boundary conditions on forces are this, this is axial thrust which is P1E rise to omega T and at X equal to L it is P2E rise to omega, so under what conditions can the displacement delta 1 and delta 2E rise to omega T can coexist with forces P1E rise to omega T and P2E rise to omega, that's a question. So we will assume now the solution to be of the form phi of X E rise to omega T, why? Because the system is linear it is being driven harmonically, system is linear time invariant driven harmonically and we are interested in steady state behavior, so such systems always display a harmonic response behavior at a driving frequency, so the solution is of this form. Now substitute this into the governing equation, the derivative with respect to time will be taken care of and I am left with an ordinary differential equation in space which is phi double prime lambda square phi equal to 0 where this parameter lambda square is given by this ratio which is independent of spatial coordinate and time, but you must notice that there is frequency of driving which is present here. Now this is a second order linear ordinary differential equation, so these are the A cos lambda X plus B sine lambda X are the solutions, so now let us start with the boundary conditions on displacement phi of 0 is delta 1, phi of L is delta 2, so now if phi of 0 is delta 1 you look here phi of 0 is A that must be equal to A, phi of L is delta 2 therefore A cos lambda L plus B sine lambda L is this, so therefore delta 1 and delta 2 must be equal to 1 0 cos lambda L sine lambda L AB, so AB therefore can be written as inverse of this into delta 1 delta 2, okay. Now let's look at the questions on the forces, now AE phi prime of X is this from the assumed solution, therefore at X equal to 0 I get AE phi prime of 0 is minus P1 which is P1 AE B lambda, at X equal to L it is P2 that is AE phi prime L is P2, so if you substitute that I get this equation, so now P1 and P2 are related to A and B through this matrix equation which is actually assembled from this equation and this equation, so that is this, so consequently now I have A and B in terms of delta 1 and delta 2 that's what I have derived, and I am expecting delta 1 delta 2 P1 P2 to coexist therefore I should get same AB, so substitute for AB I get this, so now this matrix that is AE into this into this inverse of this, this matrix is called a dynamic stiffness matrix for the rod element, so if you carry out this simplification you can show that the dynamic matrix is given in terms of cotangent and cosecant functions which is given here, so therefore I get AE lambda into this matrix is delta into P, so it is something like KX equal to P, in equilibrium equations in static is of the form, static problems is of the form KX equal to P, so now this is something like D delta equal to P where D is the dynamic stiffness matrix, so it depends on not only elastic properties it also depends on inertial and damping properties and the driving frequency, so that is why the word dynamic stiffness matrix. Now this lambda is a complex valued number due to the presence of damping, if you look at expression for lambda this I is imaginary square root of minus 1, if damping terms are 0 that means H2, C2, H1, C1 are 0 then lambda square is M omega square by AE which is real valued, but due to presence of damping the lambda is complex value, consequently the dynamic stiffness matrix is also complex value, now you can see that D is symmetric but it is not Hermitian, now if you do now if the same problem were to be analyzed using a finite element method with one element I would get this as the dynamic stiffness matrix, where M is this, K is this, and C and D bar to be suitably chosen, now we will see that quickly as an example, so let's consider a rod which is fixed at the left end and driven harmonically on the right side by a force P raised to I omega T, so and let the system have, let the system damping be made up of both structural and viscous damping therefore lambda will be given by this, so this is the equilibrium equation, at X equal to 0 delta 1 is 0, and at X equal to L P2 is 1, so now the unknowns are the reactions at X equal to 0 which is P1, and displacement at the free end which is delta 2, so if we solve these equations I will get delta 2 as tan lambda L by AE lambda, and the reaction as minus secant lambda. Now the same problem if I were to solve using the finite element method that we have developed with linear interpolation functions I will get the similar, somewhat similar equation this we have done already, this is a mass matrix, this is a stiffness matrix, and we will introduce damping in the simplified equation, so we have got this expression if you recall we have tackled this problem, so as T tends to infinity the displacement at the free end this will be the amplitude, and at T tends to infinity as a function of space this will be the solution. Now instead of using one term approximation right I will now introduce an additional node, and do now this system will have two degrees of freedom that is U2 and U3, and this end is fixed, so a two degree approximation I can get a governing equation to be this, and a reaction is this, so I can solve that problem, this problem also can be solved. Now let us assume some numerical values and plot the responses, the blue curve that you see here is a result from dynamic stiffness matrix using technometric functions, so this is actually the solution I have got here tan lambda L by A lambda, and this minus secant lambda L, so this is displacement blue line is this, the red line is finite element solution with one mode, that is one degree of freedom, the black one is finite element model with two elements, so you can see here that with one element the response, acceptability of the response ends somewhere here, over this frequency range the one degree of freedom model is acceptable. For a two degree of freedom model actually it qualitatively captures the second peak but still the answers are not good, but as far as first mode is concerned it is lot better than results from a one mode approximation, that is a single degree freedom approximation. Now how about the reactions, so this blue line is again the solution from dynamic stiffness matrix method using technometric functions, and red one is one element model that is single degree freedom model, black one is two degree of freedom model, again the qualitative features observed are similar to what we observed in displacement, the two degree of freedom model qualitatively is acceptable at least in the sense that it has two peaks, but quantitatively it's not acceptable even for the behavior near the second mode, whereas first mode it seems to be better than the solution from the one mode approximation. How about in space, suppose I fix the frequency and plot this function, this function, displacement as a function of space for a fixed value of frequency, and that I can find A and B and plot from the dynamic stiffness matrix function approach also, so I have not named the frequencies here what has been done is these lines I have selected few frequencies you can see that as frequency is varied the wavelengths of oscillations in space also increases, higher in the frequency region there are more waves within the rod, so this is a depiction of how the system responds in space for different values of frequency. So if you are using one mode approximation you will be getting a straight line as a solution in a classical FEM, okay, so that is not acceptable, if you cross the frequency and go to higher frequencies where we can expect to get more oscillations in space. Now the earlier calculation was for one chosen damping model only the velocity dependent viscous damping model was included, now I have changed some of the damping models and qualitatively I am showing the spectrum of the displacement, amplitude spectrum of the displacement, this is for presence of both structural and viscous damping, this is the reaction, this is the force, so we could play with this and see how different parameters of the problem will influence the solution, so the good thing about these solutions is that they are exact, whereas you can easily see here the finite element approximations you know aspire to serve as approximations to these curve, so they have limited range of validity in frequency and also in space. Now I have discussed this for actually vibrating rod the same exercise can be done for an Euler-Bernoulli beam, so Euler-Bernoulli beam again I consider you know one translation and rotation at each node there are two nodes and these are the model parameters, so again I assume harmonic displacements and harmonic bending moment and shear force at the ends and ask the question when they can coexist, so this is a governing equation, so again this is strain rate dependent structural damping, this is strain rate dependent viscous damping inertial term, velocity dependent viscous damping, velocity dependent structural damping, now these are the boundary conditions on displacements and on applied forces bending moments and shear forces, so again the question we are asking is under what conditions can the displacements delta K E raise to omega T and the forces coexist, so how do we tackle? Again we consider the governing partial differential equation and seek the solution in this form, where psi of X is a unknown special function in time it is harmonic, because again we are considering steady state behavior system is linear driven harmonically response will be in steady state harmonic at the driving frequency, so now you substitute into this I get a fourth order ordinary linear differential equation and the parameter lambda now has a different meaning as displayed here is somewhat similar to what we got in the problem of analysis of axial vibration but a few details are different, so this is a fourth order equation therefore the solution will be in terms of sine and cosine functions and sine H and cosine H functions, so A, B, C, D are the constants of integration. Now I can now impose the boundary conditions on displacements, so this is the equation for displacement and this is for the gradient, so at X equal to 0 psi of 0 is delta 1 I get that must therefore be equal to B plus D, psi prime of 0 is delta 2 therefore that must be A lambda plus C lambda, psi of L is delta 3 and it will involve all the terms, psi prime of L is delta 4 it will again involve all this, so in the matrix form the deltas and these A, B, C, D are related to each other through this 4 by 4 matrix, so if I am interested in A, B, C, D in terms of delta 1, delta 2, delta 3, delta 4 I need to invert this matrix. Now how about the boundary conditions on forces? So I need the second and third derivatives so we can differentiate psi of X once twice and thrice I get all these terms, now at X equal to 0 the shear force is P1 by EI that must be equal to this, the bending moment EI psi double prime of 0 is given by, similarly at X equal to L term coming from shear force and term coming from bending moment, so I can write these equations relating A, B, C, D to P1, P2, P3, P4 through this 4 by 4 matrix, we already have A, B, C, D to be equal to the inverse of the earlier matrix multiplied by the displacements, so for A, B, C, D if I write this I will be able to relate deltas to the piece, so that is this equation, so this matrix into the inverse of this matrix into deltas, so this product of these two matrices multiplied by EI is the structural dynamic stiffness matrix, so that you can simplify and show that it is of this form, okay, I use the notation S lower case S and per case S lower case C and capital C and they have these meanings, again the dynamic stiffness elements of dynamic stiffness matrix are complex valued because of presence of damping and the damping, the dynamic stiffness matrix is symmetric, now we can again quickly consider a simple example, a property cantilever in which there is an applied bending moment, E raise to I omega T at the support, so this is the nomenclature for degrees of freedom and these are the boundary conditions. So now if I use the classical one element, finite element approximation using Hermitian polynomials I will get this equation and this we have done earlier a similar problem, so the steady state amplitude of rotation here that is U4 is given by this expression, okay, this is a one degree freedom approximation to the vibrating beam. Now on the other hand using the dynamic stiffness matrix with trigonometric functions this is the equilibrium equation, delta 4 is unknown the reactions R1, R2, R3 unknowns, so there are 4 equations in 4 unknowns, the equation for delta 4 will take me to this expression and delta 4 will be this, this is exact, there is no error here as far as spatial discretization is concerned and there is no model superposition either, so let's plot this, so the red one is a result from the dynamic stiffness matrix analysis, it continues to be valid over a wide frequency range whereas the blue one is a solution from one mode approximation from FEM, qualitatively it is fine, it is picking up this first mode but there is error, I think even this error we have computed in our earlier calculations it was, that number is there in our earlier lecture, so there is a, even the location of the peak is not correct but qualitatively it's okay, at least up to some frequency range. So this is the result for slope for the same exercise, so as far as comparison is concerned I have only these two comparisons. Now for purpose of illustration we can return to one of the examples that we have tackled, a beam with a hinge and a roller, inclined roller we have solved this problem, now what we will do is we will assume that this system is driven harmonically along one of the coordinates F1 of T is raised to I omega T, and how this system behaves, so if you recall we have made a 6 degree of freedom approximation for this and I am not going to give all the details I leave this as an exercise, so you can form the 6 by 6 dynamic stiffness matrix using the dynamic stiffness matrix that we have derived just now, and go back to the finite element mass and stiffness matrices, and this is the result that I have got from using the dynamic stiffness matrix approach, for displacement that is amplitude response of U1 of T and U6 of T, okay, so I leave it as an exercise for you to compare these results from 6 degree of freedom finite element model and this approximation. Now we can also notice that, I mean we can expect that there will be a relationship between dynamic stiffness matrix and the classical elastic stiffness and consistent mass matrix obtained using linear interpolation function for say actually vibrating rod, so the dynamic stiffness matrix the equilibrium in frequency domain was obtained earlier as this, now this is the dynamic stiffness matrix, now if we use Taylor's expansion around omega equal to 0 and expand these terms in this matrix you can show that this AE lambda cot lambda L which is this term if you expand the first term will be AE by lambda, second term will be M omega square L by 3 and some higher order terms, similarly the other term also can be expanded and we get this, now if I assemble them in matrix form this matrix can be expanded the first term will be the elastic stiffness matrix, the second one will be the consistent mass matrix and higher order terms are of order omega to the power 4, this is what we are ignoring and that's why we are getting all the approximation, of course it's not a serious dead end for application of finite element method because we can always increase number of elements and increase the size of the problem and get better approximations, similarly for Euler Bernoulli beam if one were to do the Taylor's expansion we will get the, if the dynamic stiffness matrix is expanded in frequency at omega equal to 0 the first term will be the elastic stiffness matrix, the second will be the consistent mass matrix and then we have higher order terms, okay. So this is an interesting property that explains the plots that we have obtained till now, so we can make some remarks now, so the dynamic stiffness matrix approach provides exact frequency response functions for skeletal structures that is within the framework of beam and rod theories that we adopt, no errors due to modal truncation or errors due to polynomial interpolation functions occur in these calculations, so therefore these solutions can form basis on which approximate solutions can be validated, that is the, that is its use that I perceive in this context of finite element course. Now it has great flexibility in treatment of damping that again adds to its value, then the problems with the dynamic stiffness matrix is extension to 2 or 3 dimensions, for example plate problems or shell problems etc., is not straightforward and often it's simply not possible for various boundary conditions and geometries and so on and so forth, also this method ceases to be applicable if system behaves nonlinearly, because the basic fact that we are assuming is that a harmonically driven system in steady state responds at a driving frequency harmonically and that property is valid only for linear time invariant systems, so moment the system nonlinearity is included in the modeling this approach becomes invalid. Now again if structural elements are inhomogeneous specially varying say elastic constants and so on and so forth, it is possible in principle to develop the elements of dynamic stiffness matrices by solving a set of initial value problems, but that would bring in some approximations. Alternatively after going into the frequency domain a finite element like solutions can be developed for such system, by such systems I mean inhomogeneous systems linear for linear time invariant problems, after removing time we get a special problem that can be tackled using finite element method, but the shape functions then would be functions of that lambda which is a frequency. So such approaches are on a spectral finite element methods where the trial functions are functions of the driving frequencies, the trial functions are functions of driving frequencies that helps the treatment of dynamic problems as if they are static problems. This could also be trigonometric functions derived from the homogeneous elements, but the computations would become TDS, but the advantage would be that the shape functions will adapt to the frequency of excitation, we need not refine the mesh if frequency increases. Now before we close this discussion there is one at another alternative approach which is in principle does whatever dynamic stiffness matrix does and that approach is known as transfer matrix method. So to describe that we again return to the problem of axially vibrating rod, now what we do is at every station X equal to 0 I call the displacement and the force as a state vector. So we state vector at the left end is delta P L and at the right end is delta P R and I want to relate the system states here to the system states here through a transfer matrix, so that is what I want to do, so how do I do that? So we can do that by various alternative roots but since we already discussed dynamic stiffness matrix approach we can derive it from dynamic stiffness matrix, so we consider the equilibrium equation for the system in the steady state in the frequency domain which is given by this. Now we rewrite this equation, the first equation will be AE lambda cot lambda L into delta 1 minus Corsac, this P1 I will take to the right hand side and the Corsac term I will take it to the right hand side, that means what I am doing is on the left hand side I am collecting state terms involving states at left end, similarly the second equation I will get in this form. So now delta 1 P1 which is a system state at the left station is related to delta 2 P2 which is a system state at the right station through this equation, so the transfer matrix obviously will be you have to get delta 2 P2 in terms of delta 1 P1 and this product is the transfer matrix, so that is T is this, so if you multiply this we get this as a transfer matrix. Similarly for Euler-Bernoulli beam we can do a similar exercise, the system states would be now consisting of translation, slope, bending moment and shear force, okay, so the state vector is 4 cross 1, so the transfer matrix will be 454, how do we get that? So again this is the equilibrium equation in frequency domain, we expand and rewrite all this and on one side we collect terms involving states at one station and on the other side terms involving states on the right side, if you do that I get this matrix equation, 4 by 4 matrix into delta 1, delta 2 P1, P2 which is translation, slope, shear force and bending moment at X equal to 0 is related to the similar state at X equal to L through this equation. Now my objective is to relate what happens on the right side, right hand side to the left hand side and the transfer matrix is given by inverse of this into this, okay, so now I carry out, this can be multiplied I am not going to give all the details, we can code it up and see how what we get. Now a simple example we can start with an study of an undamped single span beam, so this is the beam, so you look at now left side delta 1 which is translation is 0, delta 2 which is a slope which is nonzero, unknown, P1 which is a shear force which is unknown, P2 equal to 0 which is the bending moment which is 0, at X equal to L translation and slope are 0 but there is bending moment and shear force, okay, if a system is vibrating harmonically this will be the equation, it may do so under free vibration conditions, okay, under some initial conditions. Now we expand this and we get this equation, okay, now delta 2L and P1L which is the unknowns at the left end is given by this and for non-trivial solution the determinant of this equation must be 0, so this gives an equation for omega which are the natural frequencies of the system, okay, so using this approach even we can compute natural frequencies. Now the same calculation can be done if beam is inhomogeneous and made up of say N sections, so in each sections assume that flexural rigidity and mass per unit length are changing and lengths could be changing, okay, so I can use the transfer matrix now I start with X equal to 0 multiplied by T1 I will get transfer system state here and that system state multiplied by T2 will give me system state here and so on and so forth, so when I reach the right end this will be the equation, okay, so again if you write these equations the characteristic equation will be given by in terms of the system transfer function which is product of T1, T2, T3 up to TN, okay, so this will help us to find the natural frequency, this is again exact, although the beam is inhomogeneous this formulation is exact, so the beam is of course the beam is made up of piecewise uniform sections, so we can get the characteristic equation and tackle this problem and how do we get more ships, suppose I have finished solving this, so at the end I would have got this non-trivial solutions for this, so again the eigenvectors are nonunique to an extent of a constant multiple, so one of the degree of freedom I can assign an arbitrary value, right, so what I do I write that equation for the first segment and start with a known solution, so I will get solution here that I propagate, and when I reach here if all my calculations are 0 the delta 1 and delta 2 must be 0, that's a cross check, so you get a displacement pattern and that can be normalized later on to the desired normalization criteria, so these are the steps, you fix omega at one of the natural frequencies you start with delta 2 equal to 1 at station 0 which is left hand side and you get this therefore you get P1 0 is equal to this and you can also get the reaction at the right end, although it is not needed in this calculation at this stage you can get P3N and 4N, N is the last station, this will help us to validate our calculations after we reach the, by multiplying the various matrices when we reach the right end we should get the same forces as we got here, that shows that our calculations are right, now determine the state vectors for I equal to 1 to N using the transfer matrices, upon reaching the station N the states arrived at must match the ones obtained in step 1, okay, at this stage modal deflections are available only at N plus 1 junctions, how do you get more shapes within an element? You form now a transfer matrix which takes system state from this end to some point within wherever you want the more shape, so that you can repeat at whatever resolution you want and you will get you can trace the exact motion that allows for the actual spatial inhomogeneity, later on you can evaluate the normalization constant and scale the more shapes as is needed, now once the more shapes are obtained you can do a frequency response function calculation in terms of modal superposition, but it so happens that if you are not interested in calculating more shapes you can directly do the response analysis without computing the more shapes, so how do you do that? To explain that let's consider the same problem at one end I apply now a harmonic moment, T raise to I omega T, I want to analyze the response of the system, the system has N inhomogeneous, N piecewise uniform beam sections which are all you know not identical, so again I write the equilibrium equation in terms of transfer matrices now I get this as the equation, so delta 3 and delta 4 at X equal to L are 0, P3, P4 are unknowns, at X equal to 0 delta 1 is 0 but delta 2 and P1 are unknowns, but there is an applied bending moment whose amplitude is 1, so that will sit here, so this you have to use now to propagate the you know the states, so first what I will do is we need to solve for delta 2 0 and P1 0, for that I can use the equation this equation in the first, I mean this equation and I will get delta 2 0 and P1 0, then subsequently the reaction at the other end is also obtained using this, now you want at any junction K the response you use this, up to K you multiply the transfer matrices you will get the response at any Kth junction, okay, and if you are interested in this response within a subsection you have to replace the transfer function matrix now has to be written up to that value of X, not up to the end of the element but within an element whatever distance you want to move that X you have to put in your formulation and get this, both this dynamic stiffness matrix and transfer matrix methods are special techniques they have the potential to provide exact solutions, and if you since finite element method invariably leads to approximate solutions it is nice to have an alternate solution where the solutions are exact, so skeletal structures can be handled in an exact manner in frequency domain using dynamic stiffness matrix and transfer matrix also can do many things exactly, so this is a point that I was trying to make. Now in the next lecture we will now move on to certain modeling issues we have been talking about analysis of equation of motion, we will halt that discussion for a while now, we will go to the discussion of grids and three-dimensional skeletal structures in the next lecture, so this lecture concludes here.