 We have been discussing energy methods for structural stability analysis, so in the previous lecture we introduced the axioms, if we consider a system with n generalized coordinates and if we focus attention on statically loaded structures there are two axioms, first axioms is a stationary value of the total potential energy with respect to the generalized coordinates is necessary and sufficient condition for the equilibrium state of the system. This provides the condition for establishing the equilibrium, the axiom 2 helps us to examine the stability of that equilibrium, so a complete relative minimum of the total potential energy with respect to the generalized coordinates is necessary and sufficient for the stability of an equilibrium state of the system. So based on this we analyzed a few problems and basically we considered simple systems to bring out the central ideas of the analysis, so we considered the beam column problem and in this case when p equal to, for small values of p, y equal to 0 is the equilibrium position and as p is increased we find that a neighboring equilibrium position becomes possible and that we discussed as buckling due to bifurcation of equilibrium. So the load deflection path has this feature, there are infinite number of buckling loads and the lowest one obviously is of the fundamental interest, so as the load is increased at some point when p coincides with the first critical or the first buckling load the zero position becomes unstable and two new branches emerge. In another type of problem what is known as limit load buckling example for which is shown here there are two rigid lengths supported through springs as shown here and it is loaded by a force P and as the force P is increased the deflection of the system and the load increase together and at some point both of them reach a critical stage when a faraway equilibrium position becomes possible and this is known as limit load buckling. Here for we start with here at p equal to 0 this will be the deflection and as load is increased it will go the load deflection path rises and at this point it snaps and reaches this position and this is a load deflection path and at some value for certain values of P the load deflection path is multi-valued but although it is multi-valued there are certain portion of this load deflection path is not physically realizable because those positions are unstable. We also started discussing about approximate methods for stability analysis wherein we considered the expression for the V that is total strain energy as is shown here and we postulated the displacement field in terms of trial functions where A is a generalized coordinate and we evaluated V as a function of A and then applied the two axioms to get the critical load. Here we started with the total potential energy whereas in Galerkin's method we could as well start with the governing differential equation and assume the solution in this form where phi N of X are need to be admissible and we substitute this into the governing equation and we get a residue and we demand that this residue is orthogonal to the trial functions and we get a set of capital N number of equations and this leads to the equilibrium equation as shown here and for non-trivial solutions we need to do an eigenvalue analysis to determine the critical load. Now the question now that we should consider is how to generalize these two more a difficult class of problem, for example if beam is multi-span or we consider a plate with a hole with complex geometry and a hole and things like that or a building frame which carries lateral loads and as well maybe there could be vertical loads as well, how do we perform the stability analysis of this structure, for example if we consider this type of problems suppose this is the loading pattern that is acting on the structure we would like to know if the keeping the entire loading pattern the same we increase the magnitude of these loads, so at some point we ask the question can the structure lose stability, okay, so that is the type of questions we wish to answer. So we use to tackle this problem the finite element method, in the Rayleigh-Ritz and Galerkin's method we constructed trial functions which are globally valid and as we have seen in the discussion as in the discussion of vibration problems constructing trial functions which are globally valid for these type of built up structures is not possible and what we do is we divide the domain of the structure into subdomains and within a subdomain we interpolate the field variables in terms of the nodal values of the field variables and that is that finite element method and that provides a powerful alternative to tackle problems of stability analysis as well, so as we will see shortly this leads to equilibrium equation of this form where J is the stability matrix and P into J is the geometric stiffness matrix and K minus P J is the total stiffness matrix and we can either solve the response of the structure by allowing for effect of geometric stiffness or find out the critical conditions where the structure loses stability we can do a static analysis or we can do a dynamic analysis, so we will see some of this as we go along. So let's start with as we did when we discussed vibration problems Euler-Bernoulli beam problems and then we will generalize to more complicated elements, so the total potential energy for Euler-Bernoulli beam we have derived it is in this form, so we considered Euler-Bernoulli beam element carrying an axial load P it has two nodes and it has two degrees of freedom displacement and rotation at each degree of freedom and we assume the field variable displacement field in this form, where UIs are the nodal values and phi I of X are the trial functions, so as we have seen in the context of vibration problems we, you know the cubic polynomials are appropriate to handle this problem and these are the polynomials that we have used in the dealing with vibration problems, so see it's very clear there are two nodes and there are two degrees of freedom because the functional here the highest derivative is 2 therefore the nodal degrees of freedom are V and DV by DX there are two nodes so we need a polynomial with 4 coefficient that would mean we need to use a cubic polynomial. So we substitute now that into the expression for the total potential energy and we get now the strain energy or total potential energy in terms of the nodal coordinates, so these are once we substitute for this these being polynomials they can easily be evaluated and for equilibrium we have dou U by dou UJ equal to 0 and with that we get this equation and we define K IJ as EI phi double prime, phi J double prime X DX and J IJ as this, for equilibrium this must be equal to 0 for I equal to 1, 2, 3, 4 and this I can arrange in a matrix form and this K is the elastic stiffness matrix that we have seen earlier and this is a new thing J is the stability matrix and P is the axial load, J IJ if you notice is independent of elastic constants and loads on the system so that is why it called, that's why the name geometric stiffness. So we can evaluate those integrals and I get the elastic stiffness matrix this we have seen earlier and this is a new matrix that we have derived. Now to examine the critical condition we construct a Hessian and this is the Hessian is a matrix of these elements which is KJK-PJK and this is a Hessian matrix, for equilibrium I need this equation to be satisfied and for the critical state the determinant of the Hessian must be 0. So we can make few observations, K is the elastic stiffness matrix which has been derived earlier, KG is equal to PJ is the geometric stiffness matrix, it is called as a consistent geometric stiffness matrix since the shape functions which have been used in evaluating K have also been used in evaluating this matrix, the same shape functions have been used, J is called the stability matrix. Now KT which is K-PJ is called the total stiffness matrix, so it modifies the elastic stiffness matrix by contribution from the geometric stiffness term. Now KG which is P into J is independent of EI of the beam and depends upon only the length of the beam and hence the name geometric stiffness. Now we can also see the structure of U, so if we write this and expand we can see that U using KIJ notations KIJ and JIJ this can be written as half U transpose K-PJ into U, now based on this we can consider few examples so that we get an idea on how to use this, so we start with a proper cantilever, it has 1 degree of freedom namely U4 of T, U1, U2 and U3 are 0, so this is one of the simplest problem that we can think of, so for equilibrium I have 1, 2, 3 as 0 and U4 is to be determined and this is again 0, 0, 0, U4 and I get the equilibrium equation to be given by this, that would mean U4 equal to 0 is the equilibrium state, is this stable? Now you differentiate this with respect to U4 and demand that that derivative is greater than 0, so for that to happen P critical, P should be less than some number and the critical value is 30 EI by L square, exact Euler buckling load for this case is 20 point, around 20.1997 EI by L square and the buckling mode shape according to this theory is given by 5 4 of X which is shown here, so this is a very simple application of the idea of elastic stiffness and geometric stiffness matrix. We can now consider a single span simply supported beam discretized as a single element, so there are 4 nodal degrees of freedom, 2 of which will be 0, U1 and U3 will be 0, the equilibrium equation now in the, this is K into U and in U I have zeros here and this is a geometric stiffness matrix, now extracting the equations for U2 and U4 I get the equilibrium equation as this, now the critical condition is given by determinant of SCN matrix to be 0 and that leads to this condition, so this is this, finding this condition in fact is equivalent to performing an eigenvalue analysis of K and KG, that means I have, we need to perform this eigenvalue analysis, so we can do that and we get the equation to be this and from this I get lambda to be PI square by 30, I mean by substituting for lambda this number and carrying lambda will be the eigenvalue parameter, I get the modal matrix to be this and the two eigenvalues to be this, so the smallest of these eigenvalues is helpful in determining the lowest critical load and using that we arrive at the solution that critical load is 12EI by L square, whereas the exact one is 9.869EI by L square and this is the mode shape, okay, so we can see how the mode shape look like, exact one is a sine, half sine wave and the blue one is the approximate using these two cubic polynomials. So this is a model for, simple model for critical condition, simple model for establishing critical condition for a single span, how to improve upon the results, so in the first example suppose if we introduce an additional node, okay, so now what happens we have two elements and procedure for assembling matrices and transforming into global coordinates, etc. is quite similar to what we have done in dealing with vibration problem, so that story need not be repeated now. So we have two elements and the boundary conditions will be U1 and U2 are 0 here and U5 is 0, so degrees of freedom are 3, 4 and 6, so for the first element I formulate the elastic stiffness and geometric stiffness matrix, for second element I form this and now using the fact that 3, 4, 6 are the degrees of freedom we can assemble them, either you can form the A matrices for the two elements and assemble them through matrix operations, or if you are working with pen and paper you can assemble them by observations, so we get the 3 by 3 structure, reduced structure elastic stiffness matrix and reduced structure geometric stiffness matrix. For equilibrium it is KU is equal to PJU and a critical condition is this, so this is a equilibrium equation and again by introducing lambda as PL square by 30EI now I get a 3 by 3 eigenvalue problem associated with the 3 by 3 matrix and the minimum eigenvalue if you calculate will turn out to be 0.1726 from which I get P critical as this, exact we have seen to be this, and now we have improved on the earlier solution, note that in this case I have taken element length to be L, whereas in the first example the span of the beam was taken as L, therefore we need to do a slight adjustment in comparing the results. Similarly for a simply supported beam if we now improve upon accuracy we have now degrees of freedom will be U2, U3, U4 and U6, so 2, 3, 4, 6 are the degrees of freedom, so once we assemble the structure stiffness matrix, elastic stiffness matrix and geometric stiffness matrix we get these two matrices, and hence the eigenvalue problem to be solved is now 4 by 4 because there are 4 degrees of freedom and lambda which is now this parameter P by 3EI, if I solve this problem I get these are the 4 eigenvalues and similarly by using the smallest of these eigenvalues the P critical using this method turns out to be 9.438EI by L square, L is 2L here, and this compares much better with the exact solution because we have increased the degree of freedom. Now how do we tackle problems such as this, for example there is an axial load P here and an axial load 3P here, these are externally applied, so what we do is first we analyze the structure for finding the membrane forces which are the axial thrust here, so by simple logic we can see that the reaction at this end will be 4P and the axial thrust diagram will be this, so again I will divide this beam into two elements and while forming the geometric stiffness for element 1 I will use 4P as axial force, whereas when I form the geometric stiffness for the second element P will be the axial load, so that is the only difference, so for first element the geometric stiffness is 4P by 30L, for second element 1 it is P by 30L, all other details remain the same as far as elastic stiffness matrix is concerned, so we get the eigenvalue problem like this and lambda is again obtained as this and P critical in this case turns out to be this. Now what we are doing here is if we see when I increase P the two loads are increased in same proportion, for example if this load was Q and this was P then what I will do is I will increase P by alpha P and Q also by alpha Q, that means as I increase the load all the loads increase by the same factor, there could be transverse loads and other loads as well, so I will find out the value of alpha at which critical condition is reached and that value of alpha is known as load factor, so in problems of elastic stability analysis we aim to find the load factor, by that what we mean is for a given structural configuration carrying certain specially distributed load keeping the loading pattern the same I will increase the entire loading pattern by a constant amount, so that is alpha, so the question we ask is how far alpha should be increased before the structure reaches a critical state, so this example is illustrative of that procedure that we will follow. Now this is the axial thrust diagram of the structure when it is on the verge of buckling, so this is a critical loading condition as far as axial loads are concerned, so some more examples suppose we have encastered beams I take two elements and introduce a node at the mid-span and here again we can perform the analysis and I get 40EI by 2L square as against 39.4784, so these are some simple examples that you can run through and convince yourself that these steps are correct. Now one more example, suppose there is a cantilever beam which is mounted on a spring and it is load P is applied, so we discretize there are two elements now, the first element is a beam and the second element is a spring and the first element is discretized using single element with two nodes and U1 and U2 are 0 and I have basically degrees of freedom 3 and 4, 3 is also the degree of freedom shared by K, so I construct the elastic stiffness matrix and the geometric stiffness matrix and when we assemble I use a notation alpha as KL cube by EI and I get this equation, so clearly the critical condition, critical loads will be function of K, so that parametric study has been done when alpha equal to 0 it becomes a cantilever beam and when alpha tends to infinity it becomes a property cantilever, so for K equal to 0 and K equal to infinity we know the exact solutions, so here when alpha equal to 0 this is the exact solution and when alpha tends to infinity this is exact solution, so now as alpha is varied for alpha equal to 0 I get to this factor as 2.486 that has to be compared with 2.4674 and as this number is increased this model converges to 29.993 this factor whereas it should really converge to 20.1997, but this difference is basically due to the simplicity of this model, this can be improved by introducing more number of elements and discretizing the structure with finer detail, okay. Let's make some remarks now based on what we have seen, if you recall when we found natural frequencies we solved this eigenvalue problem, this was associated with the structure elastic stiffness matrix and the mass matrix, whereas to determine buckling loads we are again doing an eigenvalue analysis and one of the matrix here is the same that is elastic stiffness matrix, there is a new matrix KG which is a geometric stiffness matrix. The eigenvalue here where the natural frequencies whereas here they correspond to the load factors or the critical loading conditions there will be N such situation in N degree of freedom model and lambda N eigenvalues provide that. Now one important difference that we have to notice between these two problems is P critical is a function of the applied load, the same structure can have different P critical for different applied loads, see for example we considered this problem because there is a load 3P the same property cantilever will have a different natural frequency than when suppose this 3P was not there only P was there you will get another load factor. In this aspect the buckling load is different from natural frequency of a structure whereas natural frequency of the structure unless you include geometric stiffness effects also into that is not dependent on applied loads. So load factor is a constant factor by which a given loading pattern needs to be increased so that the structure would be on the verge of buckling. If axial stresses are tensile the stiffnesses of the structure increases this effect is called stress stiffening, so in a more general problem this will automatically get accounted for based on the sign of the axial thrust, so that automatically will be included in the analysis. So a small exercise would be to improve upon the result that we obtained for this cantilever beam supported on the spring by introducing one more element and you can redo this exercise for different values of K and see whether we get improved solutions. Now after this point we can take a traditional route and for example compute geometric stiffness matrix for a 3D beam, so 3D beam has bending along two planes twisting as we have seen earlier that we can derive geometric stiffness for that, then we can move on to a plate and other elements, but instead of taking each of these problems one by one we can since given that we already have experience in formulating these elements we can take a look at the general formulation as applied for a 3 dimensional solid and then we can specialize it to specific structures of interest, so with that in mind we will first start with the preliminaries as you have seen it is important to include non-linear strain displacement relations in our formulations for stability analysis, so to be able to do that in a systematic manner we will quickly review the notion of green lagrange strain tensor, so we will consider a Cartesian coordinate system X1, X2, X3 in which there is an object which is supported and loaded by surface structures and body forces in certain manner and our interest is this is the structural configuration before deformation and this is the structural configuration after deformation. We are focusing our attention on a line element PQ which is in the undeformed configuration and it will occupy a position P star Q star in the deformed configuration, so this length is DS and this length is DS star, so we want to examine properties of these two line elements which will help us to define the strain, so we will consider P coordinates is X1, X2, X3, Q which is X1 plus DX1 plus X2 plus DX2 and X3 plus DX3, similarly P star and Q star have these coordinates, PQ length is DS, P star Q star is DS star, now the displacement field is clearly given by X1 star is X1 plus U1, X2 star is X2 plus U2 and X3 star is X3 plus U3, so these relations define the displacements U1, U2 and U3, so this is a standard definition of displacement field. Now components of DS are DX1, DX2, DX3, components of DS star, DX1 star, DX2 star, DX3 star, similarly components of PP star, this is a displacement field, P moves to P star and that is U1, U2, U3, so if you now find length, square of the length DS square this is given by DX1 square plus DX2 square plus DX3 square, similarly length DS star square is given by this, now you find the difference I get this equation, now each point in the deformed configuration is a function of the initial coordinates, so I have XI star as XI star of X1, X2, X3, so we are selecting X1, X2, X3 as a coordinate system, so the coordinate system is with respect to the undeformed geometry and that's what we are doing, so we can differentiate this with respect to, we can find the DX1 star, DX2 star, DX3 star using this relation and I get DX star as J into DX, and if I now write a equation for X1 star which is X1 plus U1, X2 star is X2 plus U2 so on and so forth, using that if I evaluate elements of J I get them to be in terms of displacement, gradients of displacement field to be this. So now we will look at the change in this quantity DS star square minus DS square, so I can write in this form DS star square minus DS square by 2, I can put it in this form, so what I have to do is I have for DX star in this expression should be written in terms of U1, U2, U3 using this relation, so it involves quite a bit of an algebra at the end of which I will be able to write this quantity DS star square minus DS square by 2 in terms of quantities known as epsilon IJ and DXI DXJ, and if I divide by DS I get epsilon IJ in the summand we have DXI by DS, DXJ by DS and this is nothing but the direction cosines of the line segment so I call them as NINJ, so this is epsilon IJ NINJ, what are epsilon IJ that has to be synthesized by actually performing these calculations and it turns out that these epsilon's are given by these relations. So you immediately recognize that the quantities shown in the black are the equations represent the linear strain displacement relations, whereas the quantities that are shown in the red are the contributions due to nonlinear effects and this is a complete set of nonlinear strain displacement relations, so we call this quantity epsilon IJ as strain components and the quantity DS star square minus DS square by 2DS square lambda PQ is called the strain of the line segment PQ, and that is given by epsilon IJ NINJ, so this quantity epsilon IJ that is the matrix epsilon IJ is called the Green Lagrange strain tensor at P, it is a tensor therefore the rules of coordinate transformation is epsilon prime in the transform coordinate system is C epsilon C transpose. Now if DS is equal to DS star the magnification factor lambda PQ is 0 and all epsilon IJs will be 0, so in tensorial notation we can write this in this form where repeated index implies summations and the indices run from 1, 2, 3. Now if you now take, let us focus on this expression lambda PQ NJ, now if I take N1 to be 0, 1 and N2 to be 0 and N3 to be 0, on the right hand side I will be left with only one term which is epsilon 11, so what it means is if you are taking a line segment which is parallel to X axis here, and epsilon 11 represents the magnification of that line element due to deformation, that is the definition of strain that we are interested in, similarly 0, 1, 0 if you take it is epsilon 2, 2 that means a line segment along Y axis will get amplified by epsilon 2, similarly the line segment along Z axis gets amplified by this. Now we would also like to evaluate the final direction of the line segment, so direction cosine of PQ is N1, N2, N3 which is given by this, and direction cosine of P star, Q star is N1 star, N2 star, N3 star, so now we wish to evaluate N1 star, N2 star, N3 star in terms of strains and displacements, so again we use this relation XI star is this and if I differentiate with respect to SS star I get this and this is nothing but Ni star. Now lambda PQ is given by this and from this I can manipulate and obtain DS by DS star to be given by this, so DXI star by DS therefore is obtained in this form, by writing for displacements we can write the three equation DX1 star by DS, DX2 star by DS etc. in this form, and consequently by rearranging these terms I get the direction cosines of the line segment in the deformed configuration to be given by the direction cosines in the undeformed configuration through this matrix. This is important in the definition of shearing strain, so to define shearing strain what we do? We consider a point P and erect two line segments at an angle theta and after deformation P will move to P star, A will move to A star and B will move to B star and the subtended angle will be theta star, so some notations PA is DS1, PB is DS2, P star A star is DS1 star, P star B star is DS2 star, so direction cosines of PA is N1, N2, N3, PB is M1, M2, M3 and the star quantities are for P star, B star and P star A star, theta is an angle BAP and theta star is an angle B star, A star, P star. Now cos theta, using this I can write it as N1, M1 plus N2, M2 plus N3, M3 and cos theta star is given by this, so now using the relation between the direction cosines of the line segment in the deformed configuration with those in the undeformed configuration which we derived a while before for the two line segments I write this and if we now evaluate cos theta star in terms of U's and N1, N2, N3 it can be shown that this factor, there will be these factors so that if we absorb and take it to the cos theta star on the left hand side I get this quantity to be cos theta plus 2 epsilon i j N i M j. Now I define this quantity gamma AB as a shearing strain, so why that is shearing strain, for example if theta equal to 90 degrees that means the two line segments are at 90 degrees, one of the line segment is aligned with X axis, other one is with Y axis, if you compute gamma AB it turns out to be 2 epsilon 1 2, which is in the engineering this is gamma 1 2 or gamma XY if you use XYZ notation, so similarly if you take a line segment along X axis and a line segment along Z axis, gamma AB turns out to be 2 epsilon 1 3. Now if theta is small we can approximate theta star as pi by 2 and cos theta star as pi by 2 minus theta star, and in which case gamma AB turns out to be pi by 2 minus theta star which is the definition of the shearing strain for a small deformation. So the traditional definition of strain is recovered, so based on this we accept that gamma AB given by this quantity is a measure of shear strain, so this is some basics of Green Lagrange strain tensor. Now a few observations, all displacements are measured with respect to the stationary reference coordinate system, which is the original coordinate system for the structure in its undeformed configuration, this is true no matter how large the rotation of the structure, rotations of the structure are, now the Green Lagrange strain components are equal to 0 for all rigid body motions, the rigid body motions can be smaller large. The stationary coordinates are also called material coordinates, this approach to solving non-linear problems is called total Lagrangian approach, so there are alternative notations in tensorial notation with repeated indices all these relations can be compactly written in this form where repeated indices imply summation and comma implies differentiation, and in terms of XYZ coordinates we can write this expression instead of X1, X2, X3 I will write XYZ and I get this. Now with this brief introduction to non-linear strain displacement relations we will now return to the problem of a general formulation of stability in non-linear solids, so let us consider 3D solid, let an initial stress sigma naught exist in the body, okay, now assume that the body as the body strains the stress sigma naught remains constant, the work done as strain epsilon occurs is given by epsilon transpose into sigma naught, the sigma naught is a multi-axial state of stress, in the beam problem if you recall what this assumption means, so when we analyzed beam we found out the state of axial stress under this load P and P by A was the sigma naught, sigma XX was P by A and all other stress components were 0, so as the body strain, this strain this load did work on the deformation, that is what we are using in computing total potential, now in a 3 dimensional solid the generalization of that is there exist a multi-axial state of stress in the body and this remains constant as the body strain, thus here the axial load in the beam problem was not changing as a structure deformed, that was held constant, so the same assumption if you make it means sigma naught is constant as the body strains, so the work done is given by this, now stress components there are 6 stress components initial stress components, 3 normal stresses and 3 shear stresses, and the strain components are 6, so now we compute the work done this epsilon transpose into sigma naught by including the non-linear strain displacement relations, okay, suppose if you do that I will have to put all this now into this equation, sigma XX naught will get multiplied by epsilon XX which is entirely all these terms, not just the first term, all these terms, so it is the interaction between sigma XX naught and these terms that create the issues that are related to buckling or loss of stability. Similarly I get other terms, now I will split this work done into a component U1 and U sigma, the U1 component is on the linear part of the strain displacement relations, and the other part we call it as U sigma is the contribution that is crucial for the analysis of stability of the structure. Now this U1 if you consider what it does, there is an initial state of stress in a finite element analysis so how will you take that into account, it will become simply set of equivalent nodal forces, and this is what leads to that, so to see that let us consider U1 to be this, here I am retaining only the linear part of the strain displacement relations, next I will approximate UVW using finite element approximation, so there are say S degrees of freedom, N1 and these N's are the interpolation functions I represent like this, and by assembling UVW in a matrix N and in terms of nodal degrees of freedom I write in this form and the epsilon linear is given by this acting on NE this is BUE, so if you now consider U1 this will be UV transpose B transpose sigma naught DV, so this is nothing but UV transpose B transpose sigma naught DV, this is nothing but, now for sigma naught I will write D into epsilon naught where it is a matrix of elastic constants, this is essentially the equivalent nodal forces due to presence of an initial stress, so this becomes a load vector that has to be handled. Now how about the remaining terms, that is U sigma terms, so all the terms that are present here are gradients, dou U by dou X, dou U by dou W by dou X, etcetera, etcetera, so I declare this gradients of UVW with respect to XYZ, the nine coordinates at delta, so this delta is related to UVW through this matrix, okay, so this I call it as capital XI UVW, where capital XI is this a matrix, this is similar to the matrix that we formulate when we write the elastic you know stiffness or in this case the effect of initial stress on linear part of strain displacement relation, so now I have U sigma to be this, and I am now making the substitution delta is capital XI into U, and again UVW I will interpolate in terms of that of interpolation functions and nodal coordinates, and delta consequently becomes a capital XI into N into UV, this matrix capital XI into N I call it as capital G, so this is G into UV, this is exactly similar to B into UV, whereas B was this matrix of these operators delta I have here capital XI to be this, so consequently this complicated expression can be put neatly in this form, U sigma is half delta transpose into this matrix of initial stresses delta DV, S is a stress matrix due to initial stresses, so U sigma is this, now for delta I will now make this substitution G into UV and I will get U sigma as U E transpose K sigma UV, and K sigma is this matrix, this G transpose into this stress matrix into G DV, this is the geometric stiffness matrix in this general form, so when you have to apply this theory to specific structural elements you have to simply identify this G matrix, so for example if you are dealing with Euler-Bernoulli beam, Mortimer-Schenko beam, etc. you will know the assumed displacement form from which you will be able to derive the strain in displacement relations, and you will have to evaluate K sigma with that, so let us consider a simple example of a axially deforming bar, suppose I have axially deforming bar in XZ plane, so U and W are the displacement fields, so this is the undefined configuration and these are deformed configuration, so this point now will have U2W2 as the degrees of freedom, it is not necessary that this should deform in its own axis, along its own axis, initial stresses we will assume that in this case all components will be 0 except sigma XX naught, now U and W I approximate using, there are two nodes and linear interpolation function we will use, N1 is L-X by L, N2 is X by L, so this is the representation for the field variables, there are two field variables U and W, now G matrix in this case turns out to be this, and you can do this simple calculation and find out this integral of G transpose S into G ADX, it turns out that A is, the GK sigma matrix is given by this, so this is the geometric stiffness associated with axial deformation, where P is A into sigma XX naught, now if A is taken to be much greater than P we can show that this matrix reduces to this, so if you see text books, this K sigma for axial deforming bar sometimes will be given in this form, sometimes in this form, both are valid if you accept this assumption, otherwise this is a more general form which we can use. Now we will consider a general 2D beam element and this is a beam element carrying an axial load, and this is a 2-noded beam element with 3 degrees of freedom at each node, U1 and U4 are the axial degrees of freedom, and U2, U3, U5, U6 are the flexural degrees of freedom. Following the procedure that we have developed it can be shown that the geometric stiffness matrix for this type of element can be shown to be given by this matrix, I leave this as an exercise, the kind of context in which this element becomes useful is depicted here, for example if you consider a member like this it is subjected to combined action of axial deformation as well as bending, and when we formulate the geometric stiffness matrix for this type of situation we use this matrix, the other steps involving coordinate transformation assembly etc remains similar to what we have already done in previous occasions. The rules for coordinate transformation and assembly all that is similar, so this we have seen already I will flash these slides, you can go back to the earlier lectures for details, so this is how we transform the nodal degrees of freedom in global coordinates relate the global degrees of freedom with local degrees of freedom through this relation, and this is the transformation matrix and I write U bar as T naught transpose into U where T naught is this etc, so all this we have seen, now I can substitute that in the expression for total potential and the elastic stiffness matrix in the global coordinate system we already derived, similarly the stability matrix can also be derived using this formulation, so K bar and J bar are the element stiffness and stability matrices in the transform coordinate system, so this we need not again get into all the details, now the next question is how do we solve this type of problems, so now we are equipped ourselves with analysis of 2D beam element, so we can discretize this frame into different elements as shown here, and for each element I can use the 2D beam element that we have 2-noded beam element with 3 degrees of freedom per node for which I have the elastic stiffness matrix and the geometric stiffness matrix, so first thing what we have to do is we have to perform an elastic stress analysis and find out the axial thrust in the structure, so the first step would be don't include geometric stiffness and analyze this problem, you will find the initial stress that is what we have to first determine, once the initial stress is determined for every element I know the state of initial stress, then the question that we are asking is this loading pattern that exists will be now multiplied by a factor alpha, and as this alpha is increased the initial state of stress is not I mean assumed to be constant, so the question that we are asking is as the structured strains the initial stress won't change that is the assumption, we want to find out the value of alpha at which the structure as a whole loses its stability, so this alpha is known as load factor, so the procedure for this would be to first of course we have to assemble the K matrix and perform the elastic stress analysis and find out the initial state of stress, then you formulate KG matrix, then we have the eigenvalue problem to be solved and we evaluate this, those values actually this is alpha, so I will write this in this form, number alpha that needs to be increased so that this quantity becomes 0 is a critical value of alpha at which structure loses its stability, so what we will do is we will close the lecture, this lecture at this stage, in the next class we will consider some details about this. Now following this we need to consider how this formulation can be extended to three-dimensional situations, now one way to assume that is when we developed a general 3D beam element assuming that cross section is symmetric, what we did was we analyzed each of the mode of deformation separately, axial deformation separately, twisting along the longitudinal axis differently, bending in one plane and bending in other plane, these four entities were analyzed separately and they were all subsequently assembled, so there was no qualitatively new feature associated with the built up stiffness matrix, it was only that the amount of computational effort increased, there was no new concept that entered our formulation, but whereas when we come to doing the same problem for geometric stiffness analysis that is not so, it is not just that the axial loads interact with bending, the other interaction can also take place, the initial state of stress is multi-axial, so in the initial state of stress there are torques, there are axial stresses, there are two bending moments, so there are four forms of initial stress resultants in the system, they interact, so we need to include that, I mean that is fairly obvious that that is not going to be simple by looking at these expressions, there is complex interaction between sigma XX naught and these slopes, it is not just that the axial load interacts with bending, so depending on to what extent we are able to capture these interactions, there are various forms of non-linear analysis possible, and to be able to perform a reasonable stability analysis we should capture at least a few important interactions, and we will see some of these details in the next lecture, and that discussion on beams will be followed up by discussion on plates, and let's see how it goes, we will conclude this lecture at this stage.