 After discussing few issues about computational aspects like numerical integrations, model reduction and coupling techniques that is substructuring techniques towards the end of the last class we return to the theme of development of structural matrices for different types of elements and we started discussing issues related to modeling of plane elements, plane stress and plane strain elements, so we will continue with the discussion, so we can quickly recall we talked about plane stress models where we have a prismatic object carrying in plane loads there is no loading in the Z direction, there is no body force and no surface traction so that means these edges are free from surface tractions only these edges are loaded, so what we do is we write the boundary conditions on these edges which are traction free and these are the boundary conditions there is no normal stress and there are no shear stresses and on the edge which is carrying load the loading is independent of Z, we presume that the thickness of this object is small in relation to the lateral dimensions, so based on that we interpolate that these features that are strictly valid at the boundaries are true in the interior also, so that would mean we get sigma ZZ, ZX, ZY to be identically equal to 0 for all XYZ in the domain omega, and similarly this sigma XX, YY and XY which were independent of Z on the boundaries now become independent of Z throughout the interior, so this is the basic postulate of a plane stress model and from this we use the constitutive law and get the strains and then use strain displacement relations and arrive at model for displacement, so resulting from that analysis we saw that there will be 10 unknowns in this model, 3 stress components, 3 strain components and 3 displacement components, I think there will be one more component epsilon XY, now there will be 2 equilibrium equation, 4 stress strain relations and 4 strain displacement relations, so the number of unknowns and number of equations match and these are these equations and based on that we will be able to compute the strain energy in the body. Now another 2 dimensional model which again we briefly touched upon towards the end of the previous lecture is a plane strain model, so here again the geometry of the object under consideration is again prismatic, here the thickness of the object is much greater than the lateral dimension of the structure, so for a model to be treated as plane strain model in XY plane we need that the object should be prismatic along Z axis and it should be fixed at 2 ends Z equal to 0 and Z equal to L, against movement in Z direction and lateral dimension as already said is much smaller than the thickness, now the loading should be such that the surface traction and body forces are functions of X and Y alone and they are independent of Z and there are no body forces in Z direction and this is the assumption on geometry and loading. Now at Z equal to 0 and Z equal to L, so we have something like this, so this is the valley and we are assuming Z equal to 0 and Z equal to L, L is the thickness in this case, we have assumed that the displacement component W is 0, now body is symmetric since object is prismatic the body is symmetric about L by 2 and boundary condition is symmetric loading is symmetric about this Z equal to L by 2 plane, therefore we obtain the condition that W of X, Y L by 2 must also be equal to 0, now if you consider section between Z equal to 0 and L by 2 and we consider the plane Z equal to L by 4, this will be again a plane of symmetry for the section between Z equal to 0 and Z equal to L by 2, for this section again loading is symmetric about Z equal to L by 4, object is prismatic and displacements are 0 at the 2 ends and again by virtue of symmetry we get at L equal to 4 W must be 0, so this argument can be repeatedly, we can repeat this argument to reach the conclusion that W is 0 for all Z between 0 to L, since loading and body geometry do not change with respect to Z, we also postulate that U and V are functions of X and Y alone, so this is a basic plane strain model, postulates of plane strain model, as you see here the postulates are on displacements, whereas in plane stress model the postulates were on stresses. Now starting with these displacements, so the root will be now from displacement we will use strain-displacement relations and obtain strains, and from the known strains we will obtain stresses we will be completing the formulation. So now epsilon XX is dou U by dou X, and since U is function of only X and Y, this is epsilon XX is function of X and Y, similarly epsilon YY is function of X and Y, epsilon ZZ since W is 0, epsilon ZZ is 0, now the shear strains this will be nonzero, whereas the other two strains since A is independent of Z this is 0, and W is always 0 so this is 0, and similarly V is independent of Z therefore this is 0, and W is independent, W is 0 therefore this is 0, so we get 3 strains components to be 0, and now given the strains we can use the constitutive laws and obtain the stresses, so the corresponding to the shear strains epsilon XZ and epsilon YZ the corresponding stresses would be 0, but whereas all other stresses will be nonzero, and they are related through the constitutive law for isotropic hooky and material there will be 2 elastic constants, so we get the stresses in terms of the strains as shown here, so E is the Young's modulus, nu is the Poisson's ratio. Now based on this postulate we can show that this will be the equilibrium equation, the third equilibrium equation will be automatically satisfied, so in this model there will be 3 stress components which are not known, a 4 stress components which will not be known, and there will be 3 strain components and 2 displacement components, so there are 9 unknowns and there will be 9 equations, 2 equilibrium, 3 strain displacement relations and 4 constitutive laws, so the problem formulation thus gets completed. Now again we ask this question what would happen if we substitute this model into equations of 3 dimensional elasticity, will all the equations and boundary conditions be satisfied, it so happens that in plane strain model all the field equations will be satisfied, equilibrium equation will be satisfied, since we are starting from displacements and deriving the strains the compatibility equations will be satisfied, the constitutive laws will be satisfied, but the problem will arise in satisfying boundary conditions at Z equal to 0 and L, in that sense this model also is an approximation. Now we can show that plane stress and plane strain models are mathematically equivalent, so to convert a plane stress model into a plane strain model we need to adjust the value of mu in that model, so similarly we can convert plane strain model to plane stress model by making these substitutions, both plane stress and plane strain models are approximations, in plane stress model we have difficulty in satisfying a few compatibility equations and in plane strain models we have difficulty in satisfying boundary conditions on Z equal to 0 and L, so this is the elements of theory of plane stress and plane strain model based on a linear theory of linear elasticity. Now to be able to develop the finite element model we need to write the expressions for strain energy and kinetic energy, so to be able to do that we reorganize the notations somewhat and write now the stress as a 3 cross 1 vector, two normal stresses and one shear stress, and similarly strain is written as 3 cross 1 vector and the constitutive law is now expressed in terms of this 3 cross 1 stress vector and 3 cross 1 strain vector, D is the matrix of constitutive law which is 3 cross 1, the strain displacement relations are given by this form since W is 0 this gets simplified like this. Now the expression for strain energy we have seen already this is sigma transpose epsilon DV0, V0 is a volume and if 2H is a thickness this will be the expression for this integral over area and if I now use for sigma the constitutive law I get this as epsilon transpose D epsilon DA, similarly kinetic energy is half integral over the volume rho U dot square plus V dot square and over an area integral this becomes this, so we know now the relationship between all these quantities and now we should be able to develop finite element models for these continuum elements. So we will consider the finite element model for plane stress continuum, in the discussion of beams all the structural elements are modeled as line elements therefore there were no questions about geometry of the element, but when it comes to two dimensional elements one need to consider the geometry of the element also, so we can have thus a triangular element, a rectangular element, a quadrilateral element or a triangular element with corrugated boundaries or a quadrilateral element with corrugated boundaries. Now these problems represent additional difficulties where we need to deal with the geometry of the element, now if you carefully look at this the strain energy and kinetic energy will, strain energy will lead to the stiffness matrix, kinetic energy will lead to the mass matrix, so this integral over the domain you can anticipate that that would be simple for geometries like this, but when it comes to geometries like this the integration exercise becomes difficult and at that stage we will adopt numerical integration, so we will make a substitution for coordinates and what we will do is we will map regions like this into, regions like this into unit squares and integral will be on those squares, okay, so we will come to that but let us start with simpler you know geometries, so let us start our discussion with linear triangular plane stress element. Now this element is characterized by three nodes 1, 2, 3 and the edges are 1, 2, 2, 3 and 3, 1, the nodes have coordinate X1, Y1, X2, Y2 and X3, Y3, this is the coordinate system, U is the displacement along X direction, Y, V is the displacement along Y direction, the displacement field within the element are the two displacement components U and V, they vary within the element as functions of X, Y and T, there are three independent variables, two spatial variables and one time variable, now what we wish to do is we want to express the variation of these field variables in the interior of the element in terms of values of these field variables at these nodes, so to be able to do that we begin by postulating that the displacement field can be approximated as a linear function of X and Y, so this alpha 1, alpha 2, alpha 3 are generalized coordinates which need to be determined which are functions of time. Now I have been assuming that thickness is 2H, but from this slide onwards the thickness of the element will be H to follow standard notations, now so I start with U equal to alpha 1 plus alpha 2X plus alpha 3Y, now at X1, Y1 here U is given to be U1, so I impose that condition U1 is alpha 1 plus alpha 2X1 plus alpha 3Y1, similarly at X equal to X2 and Y equal to Y2 here, I know that U is U2, so that I condition I impose I will get alpha 1 plus alpha 2X2 plus alpha 3Y2, similarly I get the third equation. Now I will substitute that write this in the matrix form U1, U2, U3 this is known, these are unknowns, so this I will write it as matrix A into this unknown quantity, so this alphas can be expressed in terms of the nodal values of the field variables using this relation, where this A inverse is given you can compute that and we get this, where this quantity A naught turns out to be this. Now to formulate this what I will do is I will write U as before alpha 1 plus alpha 2X plus alpha 3Y, and I will rewrite it as a row matrix into a column vector and for alpha 1 alpha 2 alpha 3 based on this relation I will write it as A inverse U1, U2, U3, so this is what I am aiming to do, write the field variable U of XY, T in terms of the nodal values U1, U2, U3, and these functions serve as interpolation functions, so now I denote N1, N2, N3 as A inverse U1, U2, U3 so that this is 1X, Y, A, so that I can write U as N1, N2, N3, U1, U2, U3, thus I get U is this and similarly for V there are the other field variable, there are two field variables I get a similar expression in terms of V1, V2, V3, now I can combine the two states U and V into a single vector and write them in terms of the nodal values U1, U2, U3, V1, U2, U3 as shown here, this I will write it in a compact notation, I will write U and V the vector of UV as N into UE, that E means element, so UV is 2 cross 1 and UE is 6 cross 1 and N is 2 cross 6, so this element has 3 nodes and 6 degrees of freedom, okay, the degrees of freedom are U1, V1, U2, V2, U3, V3, now let's consider this equation and make some observations, now U can be written as N1X, Y, U1 of T plus N2X, Y, U2 of T, etc., now since alpha has been selected so as to satisfy the conditions at nodes 1, 2 and 3, we get that if I put now for N1X, Y, N1 at X1, Y1 must be 1, N2X1, Y1 must be 0, N3X1, Y1 must be 0, so you can write this for N1, N2, N3 and in general I can write Ni of XJ, YJ is delta IJ where delta IJ is a chronicle delta, now what is the implication of this, so to understand that we will consider an element A marked as 1 with edges 1, 2, 3 in some Cartesian and coordinate XYZ, I have shown here the shape function N1X, Y, now N1X, Y we have seen it varies linearly, now let us consider a point A lying on the edge 2, 3, now N1 of X, Y here is 0, so when I write the displacement field here U it is not affected by the nodal coordinate U1 because that will be multiplied by N1 therefore it is 0, so only U2, U3 contribute to the value of A, sorry that field variable at A, so for any point A on the edge 2, 3 displacement is not affected by U1 and V1, now let us add another element 2 to this, now by similar argument we can make the statement that now A belongs to both element 1 and 2 because it is on the common boundary, now if I look at U from the perspective of element 2 and look at U from the perspective of element 1 what happens, that is a question I am trying to answer. Now if you look at from the perspective of element 2 displacement field anywhere here is given in terms of nodal coordinates at 2, 3, 2, 3 and 4, now a point on this edge is not affected by value at node 4, similarly this point displacement here as a viewed from element 1 is not affected by U1, so when you are looking at this point the only coordinates that affect its value of the field variable here are U2 and U3, when V2 and V3, now by the way we assemble the nodal coordinates at 2 and 3 will be common for both elements 1 and 2, that would mean that the displacement field will be continuous across this edge, there won't be any gaps or discontinuities here, that is the consequence of the way we have constructed the shape function. So upon deformation the nodes will be displaced, so all the 4 nodes will get displaced but not withstanding that there would be no gaps along the line 2, 3 in the displacement function, now we have now therefore got the displacement field in terms of nodal degrees of freedom, now let us look at energies, kinetic energy is H by 2 evaluated over the area integral rho U dot square plus V dot square, now this I will write it in this form U dot V dot transpose U dot V dot and for U dot V dot I will write this, use this N into UE so that become UE transpose N transpose and for U dot V dot it is N into U dot E and integration is on DA naught, so this U dot of E is independent of area coordinate so that can be pulled out, so the integral will be over only this quantity N transpose and N, so I can write in this form and get the strain energy in the form of half U dot transpose ME U dot, and this ME which is a quantity inside this bracket is a element mass matrix, so to be able to derive the elements of ME I need to carry out this integration. Now N and we have already determined N, so I can substitute that they are polynomials in this case and we can carry out this integration and we get this mass matrix, in this case I am able to evaluate this integral in closed form, so I get this mass matrix this is the consistent element mass matrix, so it is 6 by 6 because there are 6 degrees of freedom this is a 6 by 6 symmetric mass matrix, the components of these are the thickness of the element area of cross section and the mass density. Now how about strain energy, the stress 3 cross 1 stress vector has elements sigma xx, yy, sigma xy and this is a strain also, now sigma is related to epsilon through this D matrix and the strain is related to displacement through this relation, so this I write it as now for U and V as we have just seen I write it as N into U and I write it as BUE, where B is this matrix, this product, this product. Now I have now sigma transpose epsilon DV0 as the, sigma transpose epsilon as the integrand, so I will write this, this is over volume, this over the area it becomes this, and sigma transpose I will write it as epsilon transpose into D, D is a symmetric matrix, so D transpose is D and I get epsilon here, and this epsilon now I will write in terms of the shape function and this B matrix this leading to U transpose, B transpose DB, UE, DA. Now again this U is a function of only time, it is independent of special variable so it can be pulled out and I get this matrix as H into area integral, integral over the area B transpose DB DA0, so this quantity is called the element stiffness matrix. Here again to evaluate elements of this we have to explore the elements of this matrix, and in fact we can show that in this particular example B into this N, this matrix this into N will be this, and this is independent of X and Y, so consequently what happens this integration can be done easily and we can pull out that and we get KE, HA0 B transpose DB, so B is independent of X and Y, and hence this element is called constant strain triangle CST, now epsilon is given by BUE, sigma is given by DB, so and this is constant over the element. Now so thus we have derived now the elements stiffness matrix and mass matrix, now the steps which follow development of structural matrices, that is beyond this step that we transformation to global coordinates, assembling of different matrices which constitute the built up structure and imposition of boundary conditions and evaluation of nodal forces leading to the final equilibrium equations would be similar to those described in the context of analysis of frames, so we need not have to revisit that problem, there will be some variations in implementing those steps but conceptually the details remains similar. Now we have dealt with triangular element and we can now look at rectangular elements, so I consider a rectangular element of dimensions 2A and 2B, I will consider a problem of an element having 4 nodes 1, 2, 3, 4, so these are 4-noded rectangular element with 8 degrees of freedom, because at every node we have 2 unknowns U and V, so it is 8 degrees of freedom. Now what I do is we introduce new coordinate system called natural coordinates, psi is X by A and eta is Y by B, the consequence of making this substitution is that this rectangular region now gets mapped to a square of dimensions sides 2, so the originate 0, 0 and the vertices are at plus minus 1, okay, X equal to plus minus 1 and Y equal to plus minus 1. Now the advantage of this is as we saw we are going to evaluate elements of mass and stiffness matrices by carrying out the quadrature over this area and by making this substitution the evaluation of those integrals become simple and systematic. So again we start with displacement field we assume that it is alpha 1 plus alpha 2X plus alpha 3Y plus alpha 4 XY, similarly V is a constant plus linear functions in X and Y and quadratic terms in XY, I will discuss why it should be chosen in this form subsequently, so the field variable is U of, suppose you consider field variable U of XY, T, it needs to be interpolated in terms of values at 1, 2, 3 and 4, so there are 4 constants which has to be used for interpolation and we select 1X, Y and XY as the basis functions to construct those in the interpolants. Now how do we evaluate this alpha I, we again follow the same procedure, we use the fact that at the nodal value the XI YI it is UI of T and we can obtain alpha in terms of UI for I equal to 1, 2, 3, 4, similarly for V we again follow the same logic and we get this, so conceptually this involves certain details but conceptually it is simple, as I already said it is advantageous to introduce these transformations, so upon making these transformations we write the field variables in terms of in XI and eta coordinate as shown here, so U is approximated in terms of the nodal values U1, U2, U3, U4 in terms of these interpolation functions, N1, N2, N3, N4, this is function of XI and eta, similarly V is interpolated using the same interpolation function, now the function, this interpolation function should be such that they have this Kronecker delta property, NJ of XI, eta must be delta IJ, so we can show that this N1, N2, N3, N4 can be obtained in this particular form, so you can quickly see that at node 1 the coordinates are minus 1, minus 1, so this will be 1 and at all other values N1 will be 0, that means at all other vertices it will be 0, similarly you can verify that this condition is satisfied for this 4 function and all these 4 functions can be compactly written in a single notation as shown here, NJ of XI, eta is written like this, for I equal to J, J equal to 1, 2, 3, 4 and XI and eta J takes values plus minus 1 at the vertices. Now again we can examine the nature of these interpolation functions, suppose if I consider an element 1, 2, 3, 4 with nodes 1, 2, 3, 4 and consider the function interpolation function N1, at node 1 its value is 1 and we can see on edge 1, 2 where eta is minus 1, N1 of psi, eta is half of 1 minus psi, XI, so similarly 3, 4 I get this, 2, 3 I get this, 4, 1 I get this, so that means displace again the argument is suppose if I add one more element here and consider the value of suppose I consider a point A and if I were to add one more element the question of what the point here belongs to this element as well as this element, so would the displacement be continuous across this edge? If you examine these relations you conclude that displacement should be continuous across the element boundary. Now what are the nodal degrees of freedom? U1, V1, U2, V2, U3, V3, U4, V4, so all that I assemble in a single vector 8 cross 1 vector and I define this matrix of interpolation functions as shown here, so that I will be able to write UV as N into UE as before, so the structure of this representation remains the same but the dimensions of these matrices vary, this is 8 cross 1 and this will be 2 cross 8. Now again for kinetic energy I will use this relation rho U dot transpose U dot DV dot, and by using this relation I get this expression and now to evaluate this integral first I will do over the thickness area it becomes area integral, this is rho A BH and integration is from minus 1 to plus 1 it is N transpose N D psi T eta. Now we have seen this N1, N2, N3, N4 are polynomials in XI and eta, therefore in the evaluation of this integration is not difficult, so we can carry out, so ME gets represented in this form and the typical element ME IJ is in the form NI XI eta into NJ XI eta DXI D eta, so now for NI and NJ I can write the expression that we have already derived and I can show that this integrand gets separated into one integral over XI another integral over eta and this are simple polynomials or simple limits, so these can be evaluated exactly in a straight forward manner, so this particular example doesn't pose any difficulty, so far we have not encountered any difficulty in evaluating these integrals. Consequently I will get now the 8 by 8 consistent mass matrix for the element, it is symmetric and the quantities that enter the formulations of mass matrix are the density and dimensions of the element and this is free from any other parameters. Similarly strain energy we consider V as this is the expression for strain energy in terms of strains and the constitutive matrix of constitutive coefficients, so again if I use this representation as before I get, now of course the dimensions of these quantities will vary but the mathematically the structure of this you know the form will remain the same, so the integrand now becomes you know UE transpose B transpose DB UE and the stiffness matrix is obtained as ABH integral minus 1 to plus 1 B transpose DB D psi D eta, now this is the element stiffness matrix, so here again you can show that in this particular case we will be still be able to evaluate this integral in closed form, okay, so the B matrix in this case for a simplification I have given it here which can be used to evaluate that integral. So D transpose DB you can see that this will be quadratic functions if you put that and you will be able to evaluate the elements of stiffness matrix, so in summary we have the stiffness matrix and the mass matrix and as I was mentioning in cases considered so far it has been possible to evaluate these integrals in closed form, if the element geometry is not simple like triangle or rectangle as we have been considering we need to resolve to numerical integrations to evaluate these integrals, so this leads to questions I mean what is on a isoparametric formulation and we will see that in due course. Now before we proceed further with element development we can consider few examples, the first example that I have taken is on a shear wall like this dimensions and elastic properties are shown here, thickness is about 0.2286 meters and you can see that the conditions for plane stress modeling prevails and the problem is to determine the first 5 natural frequencies and normal modes and we can also use Timoshenko beam theory for this and arrive at the natural frequencies from a beam theory and we can see how these results compare. What we are going to do now is we will develop different models for this using triangle elements and using rectangular elements with different mesh densities and see how the results behave. So the first discretization is as shown here we use triangular plane stress elements, in this model there are 32 elements and 48 degrees of freedom and you can formulate the K and M matrices and if you evaluate the natural frequencies, the first 5 frequencies are listed here and based on the examination of normal modes which I will be showing shortly the second and the fifth mode in this particular model turn out to be axial deformation model, that means oscillations are predominantly in this direction otherwise they are flexure. So how do the mode shapes look like? This is the first mode shape, flexure, the second one is axial deformation, this is again flexure, flexure, axial deformation, flexure. Now how do we know that this is acceptable? So to proceed further we have to refine the mesh, so the thing that I do is I will refine the mesh as shown here. In any mesh refinement problem the idea of refinement should be such that the nodes that we have produced in this model should again be remain as nodes in the refined model, so if you see carefully I have all the original nodes in the first model are retained here, but the number of elements have increased, so in this case there are 128 elements and there are 160 degrees of freedom and if I use again the free vibration, if I perform the free vibration analysis I get the natural frequencies as shown here, in this model the second mode become flexural model and third model becomes axial model. Here these frequencies have not yet converged with respect to the element size, so the errors in different natural frequencies need not be the same, so if the second mode and third mode have different amounts of errors and they are converging at different rates the order in which they appear will also switch, so that's what happens here, you can see here model the second natural frequency is axial model whereas here it becomes third mode. Now how do the mode shape look like? You can see here this is the first mode, this is second mode, this is axial mode, third mode, this is the fourth mode, these are the higher modes and I get two higher order bending modes and axial deformation modes, it is difficult to discern the patterns of axial deformation there will be nodes and things like that that is not easy to determine because the deformations are shown in the direct, the deformed geometry is also shown in the same direction so it will not be easy to delineate whether there is a node and things like that. Now if we summarize these three models also somewhere I indicated here if we do a Timoshenko beam theory and compute the natural frequencies exactly, there is no finite element model I get the natural frequencies to be this, the first two modes are flexural modes, third one is axial, fourth one is flexural and fifth one is again axial. So if we tabulate this for model 1, if I assume that since this is an analytical model and natural frequencies have been computed exactly and also the geometry of this structure is such that Timoshenko beam theory is likely to be valid, we can treat that as more exact than these two, if we agree with that postulate we can see that the first model has 6.39, the next one is 5.35 they seem to be moving towards this value of 4.97, here the second mode is axial mode and that should be compared with this mode, this is 31.99, 31.95 and this is actually 31.944, so what do we see? We will come to that point again, now the second mode this 26.391 it is 27.58 and that has to be compared with 32.19, so 32.19 and 27.58 are approximations to 26.39, similarly now the higher modes you can compare and the thing that we can observe is the axial modes are captured well, for example 31.99, 31.95 are good approximation to 31.944, whereas 32.19, 27.58 and 26.391 are relatively poorer approximations, similarly first mode 6.39 and 5.35 are approximations to 4.97, so the model of the story is in triangular elements the axial modes are more accurately captured than bending modes, is not quite a satisfactory situation because this model is already fairly complicated, in spite of it we are not getting good results. Now let's solve this problem using rectangular elements, so let's start with this discretization, again everything remains same, for this discretization it has 160 degrees of freedom and it is same as this, okay, so now the model that I get is 5.03 which is for approximation to 4.97, 26.22, 26.39, 31, 31.95, 61, 62, 96, 95, so here in addition to axial modes being captured well with rectangular elements where we are also able to capture bending modes with higher level of satisfaction. So the mode shapes are shown here, the ordering matches with the order that is observed in Timoshenko beam, so from this one could conclude that rectangular elements behave better than triangular elements at least for the examples that we have considered. Another example that we can consider is that of a triangular wedge which can serve as a model for an earth dam and one can analyze this triangular wedge as a shear beam and you can since the variation is linear we can that is area of cross section variation is linear along the height we will be able to get an exact solution in terms of Bessel's functions that is available in the literature. So this example now if we analyze using plane stress elements, plane strain element this will be now a plane strain problem because this is an embankment problem, so how does our models behave? So the first scheme of modeling is shown here and here we have 42 elements and 42 degrees of freedom and these nodes will be fixed at the bottom and the first few frequencies in Hertz we obtain as 1.25, 2.67, 3.70, 4.0 and 4.1 and from the shear beam model analytically computed natural frequencies are shown here. How do the mode shapes look like? These are the mode shapes. So now this should give us an idea within an element how well we are approximating the deformation fields, okay. Now the same problem can be tackled using of course refined triangular mesh and so on and so forth but we can also use combination of rectangular and triangular elements, so an illustration for that this is not a very good mesh, it is not something that I recommend for you to use but it serves to illustrate how triangular and rectangular elements can be used in the same model. So for this part we have used triangular elements, similarly triangular elements here, triangular elements here and for this part we have used rectangular elements. Now again we compute natural frequencies and these natural frequencies turn out to be this. So the mode shapes are obtained as shown here, you can see the shear deformation type of behavior here and this is mode 1, this is mode 2, this is mode 3 and this is mode 4. Now we can do this problem more systematically by refining the mesh in an orderly way, maybe at some point I will show how that can be done and we can understand how these natural frequencies behave with respect to the refinement in mesh, but the idea of showing these illustrations here is to simply illustrate that the two elements that we have developed can be used to solve this problem, and mind you for this type of geometry entire machine cannot be done with only rectangular element, we will end up using triangular elements to capture these details at the corners at least, the remaining part can be meshed with rectangular elements, but when it comes to these corners we will end up requiring triangular elements. Now this shows us that there will be situations where simple geometries like rectangle and triangle may not suffice, so we need to deal with more complicated geometries such as this. So how do we proceed? So what we do is the basic idea of considering different geometries is we have some Cartesian coordinates X and Y, so a point here has coordinates X and Y, now I will introduce a mapping so that this region gets mapped to a square with vertices at plus minus 1 as shown here, so this I call it as XI and eta, now a point here gets mapped to a point somewhere here, so we can assume that X is function of XI and eta, Y is function of XI and eta. Now the idea is that we will now represent this in terms of the nodal values, okay, and use the same interpolation function that we have used for deriving the shape function, for representing the displacement field, okay. So what did we do for displacement field? We wrote the nodal values and this interpolation functions, okay. Now I will use the same interpolation functions and represent these coordinates in XI eta plane, so that means the geometry of the element is also approximated using the same shape functions which were used in representing the field variables, so this formulation is known as isoparametric formulation. Now if we go back and look at the structure of the stiffness and mass matrices, we have got stiffness and mass matrices in this form, so these integrals will now be in the form of, will now be with respect to the natural coordinates XI and eta, so this function because of the nature of transformations involved will no longer be simple functions which admit close form solutions, so what we do is we use numerical integration schemes to derive the elements of KE and ME, that specifically what we use is the Gauss quadrature rules are used to evaluate these integral, so what we will do in the next classes, we will introduce the isoparametric formulation and derive the mass and stiffness matrices for quadrilateral element and revisit some of these problems and see how best we can, what would happen if I use quadrilateral elements and solve the problem, so with this we will conclude today's lecture.