The following content is provided under a Creative
Commons license.
Your support will help MIT OpenCourseWare continue to
offer high quality educational resources for free.
To make a donation or view additional materials from
hundreds of MIT courses, visit MIT OpenCourseWare at
ocw.mit.edu.
PROFESSOR: Ladies and gentlemen, welcome to
Lecture Number 10.
In this lecture, I would like to discuss a solution of
dynamic equilibrium equations.
And in particular, I would like to talk about the direct
integration solution of the equilibrium equations in
dynamic analysis.
These, of course, are the equilibrium equations that we
derived already earlier.
M as a mass matrix, C as a damping matrix, K as a
stiffness matrix, U as a displacement vector, and here
we have the velocities and the
accelerations at the null points.
R, of course, being the load vector and this load vector is
now time dependent.
I also should draw your attention on the fact that I
would like to talk about direct
integration of these equations.
We will talk about the modes of opposition procedures in
the next lecture.
By direct integration, we mean that we integrate these
equations directly without--
that is, without--
a transformation of these equations into a different
form prior to the integration.
In modes of opposition analysis, now you will see in
the next lecture, we're actually transforming these
equations first, and then we integrate the response.
So the procedures that I'd like to talk about are
explicit, implicit integration procedures that directly
operate on these equations without a transformation and
integrate, solve these equations.
I also would like to talk about some computational
considerations with regard to these explicit and implicit
integration procedures, the number of operations involved,
the costs involved in performing these
integrations and so on.
And then there is an important consideration, namely the
selection of the solution time step, delta-t.
We will see that when we integrate the dynamic
equilibrium equations, we have to select the time step
delta-t with which we perform the integration.
And this time step has to be selected in order to preserve
stability of the integration and to preserve of course or
to obtain accuracy in the integration.
Finally, I'd like to talk about some modeling
considerations pertaining to the solution of the dynamic
equilibrium equations.
Well, here on this view graph, I'm showing once again the
dynamically equilibrium equations here, noticing once
more that R, the load vector, is now a function of time.
The same equations here have been rewritten once here,
where FI denotes inertia forces, now time
dependent, of course.
The mass matrix is constant, but the accelerations of
course depend on time.
Then this part here, which are the damping forces, has been
expressed here as FD a function of time again.
The C matrix is constant with time, but the velocities, of
course, vary with time.
So here we have damping forces.
It's a nodes of the finite element system.
And these are the elastic forces, KU, which are really
due to the internal element stresses.
And on the right hand side, we have the external forces.
So what we're saying is that we want to solve this force
equilibrium over all times t.
The inertia forces plus the damping forces plus the
elastic forces, it's the nodes of the finite element system,
of course, calculated in the virtual work sense, the way
we've been talking about in the earlier lectures.
The sum of these forces here must be equal to the
externally applied nodal point logs.
Well, the procedure that we will be
following is the following.
We will consider time steps delta-t a part at which we
want to calculate the response of the system.
In order to do so, of course, we have to
define our loads first.
And here, I have schematic shown once the description of
the load Ri, that is, the load at degree of freedom i as a
function of time t.
Notice that we are having here this distribution of that load
component as a function of time.
The distribution is quite simple.
In a computer program, of course, all we would have to
put in is this point, this intensity here, that load
intensity, at a specified time.
That load intensity, at a specified time, and that load
intensity at a specified time.
In between, these intensities, the computer program would
simply interpolate linearly.
This is the load at the component i, displacement
component i.
This is the load at displacement component j.
Notice here, of course, we would have in general a
different curve at degree of freedom j, then the curve is
used at degree of freedom i.
I show here time steps, discrete time steps, delta-t1,
delta-t2, delta-t3.
These of course apply for the complete system, for the i
component and the j component of displacements or loads, and
for all other displacements and loads also.
Notice that they can vary with time.
In many analysis, however, we simply pick one time step and
keep it constant all along.
The objective is the following.
In the direct integration of the equations of equilibrium,
the objective is knowing the initial conditions.
In other words, knowing the displacement at time 0, the
velocities at time 0.
We want to calculate the displacements, velocities, and
accelerations at these discrete time points.
We are matching here with the time interval delta-t1.
And here we are changing the time interval to delta-t2.
Knowing all the velocities, accelerations, and
displacement at this time, we want to calculate the same
quantities at the discrete time increments of
delta-t2, and so on.
So basically what we are saying is that knowing the
initial conditions, we want to calculate in a step-by-step
procedure, the displacements, velocities, and accelerations
at these discrete time points.
Of course, the real problem then reduces to the following.
Given the initial conditions at the particular time t, we
might call this time t, and by initial conditions, I mean we
know the displacements, velocities, accelerations, we
know all the quantities at time t.
Our objective is to calculate the unknown quantities at time
t plus delta-t1 in this particular case.
So this point is t plus delta-t1.
In general, or in many cases, we keep a constant time step,
and so we will simply talk about t plus
delta-t from now onwards.
We should remember, however, the delta-t might change doing
the response solution.
I repeat once more, this is an important point, given the
initial conditions or given the conditions at time t, our
objective is to calculate the conditions, displacements,
velocities, accelerations, at time t plus delta-t.
In the explicit integration, we proceed in
the solution as follows.
We use the equations of dynamic equilibrium at time t
to obtain the solution at time t plus delta-t.
This is important.
In explicit integration, we use the condition, the
equilibrium equations, the equilibrium equation, the
dynamic equilibrium equation, and by that I mean, just to
refresh your memory, this equation-- we use that
equation at time t in explicit integration to obtain the
solution at time t plus delta-t.
Once we have the solution at time t plus delta-t, of course
we can march ahead in the same way again.
In implicit integration, however, we are looking at the
equilibrium equations at time t plus delta-t t to obtain the
solution at time t plus delta-t.
In other words, if you look at once more at the equations
here, at the equilibrium equations, these equations
would be applied at time t plus delta-t to obtain the
solution at time t plus delta-t.
This is indeed the fundamental difference between explicit
and implicit integration.
In explicit integration, we are looking at equilibrium
equations at time t to obtain the solution
at time t plus delta-t.
And in implicit integration, we are looking at the
equilibrium equations at time t plus delta-t to attain the
solution at time t plus delta-t.
In the following now, I would like to show to you an
explicit and an implicit integration scheme.
The explicit integration scheme that is abundantly used
in practice is the central difference method.
I will later on then discuss how this method
is effectively used.
In the central difference method, we are applying, as I
said earlier, the equilibrium equations at time t, because
it's an explicit integration scheme.
The two additional questions that we're using to solve for
the unknown displacements at time t plus delta-t, is
firstly an equation for the accelerations at time t.
And notice that this equation involves as the unknown, as
the only unknown, the displacement
at time t plus delta-t.
In addition, of course, we are using an equation for the
velocities at time t and again this equation only involves as
the only unknown the displacement
at time t plus delta-t.
Substituting these two equations into the equilibrium
equation at time t gives us one equation in the unknown t,
displacement at time t plus delta-t.
We can simply solve for it, as shown, now on
the next view graph.
You can see here the equation that is obtained by
substituting the relations for the acceleration and for the
velocities in the equilibrium equation at time t.
The unknown appears here.
This is the displacement at time t plus
delta-t, which is unknown.
And there is a coefficient matrix here, that involves the
mass matrix and the damping matrix.
And on the right hand side, of course, we are calculating a
load vector.
And that load vector, notice once again, corresponds to the
conditions at time t.
That's why it's an explicit integration scheme.
The difference operator, the essential difference operator
that I've been used here, simple operators which you
have seen on the previous slide, that involve only the
displacement at time t minus delta-t and the
displacement at time t.
And notice that here are these two
displacement vectors appearing.
Notice that this equation is cheaply--
or can be used in a very cheap way to calculate the
displacement, the unknown displacement at time t plus
delta-t, if the mass matrix is a diagonal matrix.
And, say, the damping is neglected.
In other words, this term not being there, this being a
diagonal mass matrix, then the solution is very trivial.
In fact, all we need to do is calculate the right hand side
load vector and divide each component of that load vector
by this coefficient here which, if m once again is a
diagonal mass matrix, is simply Mii
divided by delta-t squared.
And in fact, this is really how was the procedure is
usually applied, namely, with a diagonal mass
matrix and no damping.
The damping usually is not taken into account, or at
most, we have diagonal dampers also, and the mass matrix
being a diagonal matrix, then the solution is
very cheaply performed.
There's one further important point.
This load vector of course is constructed in the usual way.
This multiplication here is quite cheap if the mass matrix
is a diagonal mass matrix.
K times tU.
Well, K here is of course a banded matrix.
It's the stiffness matrix of the system, constant.
Well, there are some operations involved here.
And if we look at this operation here in detail, we
see here how it can be effectively performed.
We have a stiffness matrix times the
displacement vector here.
And of course, the stiffness matrix remember is obtained by
summing the element contributions.
So you can write this K matrix in this way.
Since tU is independent of the elements, we can also take it
into the summation sign.
And if we look at this product closely, we find that is
nothing else than tFm.
It's a nodal point vector of loads or elements forces that
correspond--
element elastic forces-- that correspond to the current
displacements.
And this tF vector is much more cheaply calculated than
the Km times tU.
And in fact, this is the way how we actually perform the
calculation of this product here.
Now notice the following important further point.
If we were to leave this K times tU the way it stands
here, we would have to assemble a global stiffness
matrix of the system.
However, if we write that product this way, we never
need to assemble a global stiffness matrix
of the total system.
In fact, all we need to do is calculate element forces and
add these elements forces via the direct stiffness procedure
the way I've been discussing earlier.
So this is a very effective way of calculating K times tU.
We never need to assemble a global stiffness matrix.
In fact, this is one of the main reasons why the central
difference method is effective.
It's this reason that we can calculate the KtU part, this
times that, very effectively via this procedure.
That is the first reason.
And the second reason is that if M is diagonal, we really
never solve equations.
The equation solving here is so trivial, that we might just
say that we do not solve equations.
There is no factorization involved.
There's no factorization of a stiffness matrix or another
coefficient matrix involved, because if M is diagonal, and
as I pointed out, we simply divide the right hand side by
Mii divided by delta-t squared.
Notice that if we start the integration process, we know
of course the initial conditions.
We know tU.
And we do not know--
usually directly the t minus delta-tU.
That is the displacement prior to the solution.
And to obtain those, we can use this formula, which is
really the central difference method applied at time t,
where we however now involve the accelerations
corresponding to time 0.
We involve the accelerations corresponding to time zero
right here, and the velocities corresponding to times 0, and
the displacements corresponding to time zero.
And of course, these quantities must all be given
as the initial conditions to the solution.
I pointed out already that we're using this technique
mostly with the lump mass matrix.
In addition I should also point out that we're using it
primarily with lower order elements,
with lower order elements.
The reason for that is the following, namely, one of
stability and accuracy of the central difference method.
The method can only be applied if delta-t, the time step that
I'm talking about here, is smaller than a
critical time step.
And that critical time step is given by the smallest period
in the system divided by pi, 3.14 et cetera.
The smallest period in the system can be
of course very small.
And therefore the critical time step that we're talking
about here can be very small.
The smallest period in the system, just to give you an
example, would be extremely small if you have a bar
element, and you have some very short
truss element in there.
The shorter the truss element, the stiffer of course it is.
AU over L is the stiffness of a truss element, and
therefore, a very small period would be in that
assemblage of elements.
Therefore, what is required really is to have a mesh that
is more or less uniform, so that there is no artificially
small time step required in the integration of the
equilibrium equations.
Because there is this condition that
delta-t must be smaller--
we should really say must be smaller or equal--
to delta-t critical, but to be on the conservative side, we
are selecting delta-t smaller than delta-t critical in
practical analysis of complex meshes.
Because this condition has to be satisfied, we're talking
about a conditionally stable scheme.
It is only stable when delta-t is smaller than that value.
If we were to select delta-t larger than delta-t critical,
we would see that after just a few time steps, the solution--
it gets out of bounds.
And after 10 or 20 or 30 time step, in a computer program,
you would all of a sudden notice
overflow of the numbers.
So the solution wouldn't make any sense whatsoever anymore.
And this is an absolute condition that
does have to be satisfied.
There's no way you can use a time step larger than delta-t
critical, and still get an acceptable, reasonable
solution to your equations of motion.
In practice, of course, we have to estimate, in a
conservative way, a time step.
And that means we have to have an estimate
also of delta-t critical.
For continuum elements--
by continuum elements, I mean plain stress, plain strain,
axis a matrix, three-dimensional element, the
way we have been discussing it earlier.
A good way of obtaining time step delta-t is to use this
formula here, delta-t being delta-L over C. C is the wave
velocity through the system.
[UNINTELLIGIBLE]
is divided by rho, square root out of it gives us C. And
delta-L is an element length that has to be defined.
In fact, it's the smallest distance between nodes.
The smallest distance between nodes if we talk about low
order elements.
What I mean by that is the following.
If we have a mesh here that looks like that.
And I deliberately show elements that are not equal in
size, then delta-E in this particular case would be this
distance here.
That is our delta-L, sorry, our delta-L is
that distance there.
This, of course, as I have here in the heading, assumes
that we have low-order elements.
That means only corner nodes in the elements.
For higher-order elements--
let me put down here one higher-order element, where we
have these nodes now, say an 8 node element.
In this particular case, we take the smallest distance
between nodes, again, just as in this case.
And in this particular case, it would be, say, this
distance here.
That is our smallest distance between nodes.
However, we have to also divide this distance by a
relative stiffness factor.
And why does that relative stiffness factor
come into the picture?
Well, it does come into the picture, because remember that
the stiffness at an interior node is larger than the
stiffness at a corner for the element.
The amount it is larger depends on the particular
element used.
But since there's a lot of stiffness here, we have to put
here a relative stiffness factor into the calculation of
delta-L, the effective lengths that we're using in this
particular formula.
Now, for this element here, for example, the relative
stiffness factor would be conservatively
just equal to 4.
But for other elements, and particular also when there's
node shifting involved and so on, the relative stiffness
factor has to be calculated, has to be estimated.
And we should always remember that it should be estimated
here in a conservative way, because a smaller delta-L gave
us a smaller delta-t.
And then we are then conservatively applying the
fact that delta-t shall be smaller than the delta-t
critical, value that I have been referring hereto.
The method is used mainly for wave propagation analysis.
The reason being that we have to use a fairly large number
of time steps to predict the solution.
Large, because a time step delta-t is usually very small.
And in wave propagation analysis, many modes of the
system are being excited.
So we really want to pick up a large number of
modes in the mesh.
And the small time step integrating accurately many
modes in the mesh goes together with the fact that we
want to predict many modes in the mesh or a large number of
modes in the mesh, and therefore the method is
effectively used for wave propagation analysis.
The number of operations are proportional to the number of
elements that we are using and the number of time steps.
Of course, this is an approximate proportionality
here that I'm talking about.
This fact derives directly from the consideration that we
talked about here.
Notice that as we include more and more elements in the
solution process, of course, we have to sum over more
elements here.
And as we go more and more time steps, we have to apply
this equation as many times also more and more.
The important point is that we don't assemble a K matrix, so
the bandwidth off the system, the bandwidth of the system
does not appear in an operation count.
It's the number of elements that do appear in the
operation count.
The next procedure that I like to talk about
is the Newmark method.
And in the Newmark method, this being an implicit
integration scheme now, we proceed as follows.
We apply--
because it's an implicit integration scheme--
the equilibrium equation, the dynamic equilibrium equation
at time t plus delta-t.
Here it is given the additional formula that we're
using is this one here for the velocity at time t plus
delta-t, and this one here for the displacement of time t
plus delta-t.
Now notice that with this set of equations, with these three
equations, knowing the displacements and velocities
at time t and accelerations at time t, we have with these
three equations, three unknowns.
The velocity at time t plus delta-t, the acceleration at
time t plus delta-t, and the displacement
at time t plus delta-t.
These are the only three unknowns that appear.
There are three linearly independent equations.
And we can simply use them to
calculate these three unknowns.
Of course, we want to do so in a very effective way.
And I will show you how we proceed.
The important point is that with this scheme, we are using
an implicit approach, as I defined it earlier.
And when we do substitute from here into there, we will find
that on the left hand side, we have a coefficient matrix that
involves a bandwidth even if the mass matrix is diagonal.
And this coefficient matrix here is written
here down as K hat.
This coefficient matrix invoice this K matrix.
I should say maybe a few words to these formulae.
These formulae can be obtained by Taylor Series expansion,
around the conditions at times t plus delta-t.
They have been proposed originally by Newmark.
That's why it's called now the Newmark method.
Newmark in particular used the values of delta equals 1/2 and
alpha equal to 1/4.
And then it's a constant average acceleration scheme.
The constant average acceleration schemes use these
values as I pointed out already earlier.
The method is unconditionally stable.
This means that the time step delta-t does not to have to be
smaller than a certain value, as it was the case in the
central difference method.
The method is primarily was for analysis of structural
dynamic problems.
I will talk about why this is so just now.
And the number of operations that we're talking about here
are an initial factorization of the K hat matrix.
This is the K hat matrix, And that initial factorization
involves this K matrix, so the bandwidth of the K hat matrix
is the same as the bandwidth of the K matrix.
And therefore the operations here are one half NM squared.
And for each solution, in time, because this equation
here has to be solved for each time step, you're starting
with delta-t, and then we go to 2 delta-t, 3 delta-t, 4
delta-t, and so on.
And since this solution is to be formed for each time step,
however, remember the factorization is only done
once because this K hat matrix is independent of the time.
Therefore, we have a forward reduction and back
substitution in each time step.
And the number of operations in each time step are 2nm. t
is the number of time steps here.
t is the number of time steps here.
Well, let us look then at the accuracy
considerations that we have.
I mentioned already the method is unconditionally stable,
which means that the delta-t time step that we're using can
be any amount, and we will never have a blow up of the
solution, a blow up of the solution due to round-off
errors that I involved.
In the central difference method, of course, the
explicit time integration scheme, I mentioned that
delta-t has to be smaller than a critical time step.
So in the implicit Newmark integration scheme, what we
can do is we select our of time step based solely on
accuracy considerations.
Let us consider then these equations that we want to
solve, and I'm now leaving out the damping term in order to
obtain some understanding of how we would select the time
step, delta-t.
Of course, remember the larger the time step, the less normal
operations I involved in the solution.
So we want to select the large time step.
However, if we make the time step extremely large, then
surely we cannot expect any accuracy in the solution.
In fact, if delta-t is very large, we simply predict the
static solution.
This is one property of the Newmark method.
If you were to use a computer program with the Newmark
method, if you're interested in it, why don't you try, say,
for example, using ADINA, in which we have the Newmark
method, and just put delta-t very large, and all you would
be getting out is a static solution.
Well, since we're interested in the dynamic response, we
have to therefore select delta-t small enough to
integrate accurately the solution response.
And these are the equations that I want to now consider,
and I want to show you how can we select delta-t rationally
and to obtain the required accuracy.
The displacement vector here can be expressed in terms of
the mode shape vectors times xi values.
I will be talking about how these phi i vectors are
calculated in the next lectures.
This involves the solution of the eigenvalue problem
written down here.
The free vibration eigenvalue a problem corresponding to
this equation.
I like to defer the discussion at this point, because we will
be talking about it later on in detail.
All I like to say is that this transformation here, phi i
being a vector time independent and xi being a
scalar, dependent on time.
This transformation can be used in this solution here, or
can be used rather to transform this system of
equations into a different form.
And in that transformation, we use the fact that when we
calculate phi transpose K times 5--
notice these are capital phis.
I've crossed bars here, capital phis.
This is a lower phi, no crossed bars top and bottom.
This capital phi here stores all of the lowercase phi
vectors, and their n of such vectors.
If I use this capital phi and calculate this product here, I
get a diagonal matrix.
In fact, the diagonal matrix here stores on its diagonal
the frequency squared, the free
vibrational frequency squared.
Also phi transposed M phi gives us an identity matrix.
These, once again, the vectors in this matrix, capital phi,
are these phi i vectors which are the
eigenvectors of this problem.
And I will talk later on about how we calculate them.
Well, if we apply this transformation--
whoops, let me put it right once more-- if we apply this
transformation to this system of equations, and we use this
property right here, we obtain n decoupled
equations of this form.
X double ought i plus omega i squared, xi equals phi i
transposed times R. And they are in such equations, in
other words I run c from 1 to n.
In fact, if we use multiple position analysis, we actually
perform this transformation.
What I now want to do is simply refer to it
theoretically.
I do not want to transform it.
You're talking about direct integration of the solution of
equations, in which there is no transformation involved.
So I do not want to go through this transformation.
However, I want to refer to it theoretically.
And it gives me a tremendous insight into what happens.
So basically, what I can say theoretically is that the
direct step-by-step solution of this system of equations,
using the explicit integration scheme, the central difference
method, or the Newmark implicit integration scheme,
corresponds completely to the direct step-by-step solution
of n equations, of this form, where the displacements here,
U, are given in this way.
x being now a function of time, so if I differentiate
here, your double dot would be obtained by differentiating xi
also twice, with respect to time.
Once again, this is an important,
very important, fact.
I look at these equations, and I say that's a solution of
these set of coupled differential equations.
It's completely equivalent to the solution of n equations of
this form with this transformation used.
No dots here or dots there.
In other words, this equation holds for displacements
without the red dots and for the accelerations
with the red dots.
It is equivalent.
This solution is equivalent to that solution.
However, in direct integration, I never perform
actually the transformation.
All I know is that when I operate on this system of
equations, using the Newmark method or the central
difference method, I, in fact, operate on n
such decoupled equations.
That's important.
And this means that for an analysis on our direct
integration scheme, we can in fact analyze our direct
integration scheme by looking at these equations, instead of
these equations.
You see if I were to analyze my direct integration scheme
using these equations, I would have to deal with this huge
matrix K, the huge matrix M. All I now need to deal with,
in the analysis of the direct integration scheme, is with n
such simple one degree of freedom equations.
And in fact, I don't even need to deal with n such equations,
for the analysis purpose, I only need to deal with one
equation, where I make now omega a
variable and r a variable.
So my analysis on this equation, using this
particular time integration scheme that I want to analyze,
is completely giving me all the information that I need to
know when I use my direct integration scheme in the
solution of this equation.
For analysis purposes of my direct integration scheme, I
can simply analyze my direct integration scheme using this
equation and get all the relevant
information from that analysis.
This is an extremely important fact.
Once again, I do not actually do perform the transformation.
We don't do the transformation.
Only in theory for analysis purposes.
I'm looking at the solution of this equation using the
Newmark integration scheme.
And I get all the relevant information on the accuracy,
stability, et cetera of the integration scheme by looking
at the solution of this equation with that direct
integration scheme.
Now, even in the analysis of this
equation, of course, remember--
or in the use of this equation, remember that r is
available now.
And I have some choice what to put in there and so on.
But the important point is that what I want to identify
is the characteristics of the integration scheme.
And the characteristics of the integration scheme are really
quite well seen already if we look at a very simple case,
namely this case.
I set r deliberately equal to zero.
And I look at the case where I have initial condition x, the
displacements being one, the velocities being zero.
this means, of course, that the accelerations are minus
omega squared.
If I use the Newmark method, the one that I just talked to
you about, on that system of equations, then I can plot
errors that reside in that integration.
And what I've plotted here the percentage period elongation,
pe divided by t times 100, which is this amount here that
we see when we integrate the differential equation that I
showed you via as the Newmark method with different delta-t
over t values.
On the bottom here, we are plotting delta-t divided by t,
the period t, of course, being given for the particular
oscillator that we're looking at.
We are starting off with one.
The exact solution, I should take better now another color,
the exact solution would look something like that,
finishing off here.
That would be the exact solution.
And what we're seeing now is an amplitude decay.
It's already there.
It would come in right there.
That is the exact solution.
We'll stop here.
The period is of course this t here, as shown.
This would be the exact solution for one period.
What we are seeing as an amplitude decay, a drop in the
numerical solution, and period elongation in
the numerical solution.
And the percentage period elongation is plotted here.
We notice that the Newmark method corresponds to this
curve here.
So with delta-t over t being 0.10.
We have about--
let us just look it up here--
we have about 3%, 3% period elongation
for the Newmark method.
I'm also showing other techniques here, which are
also implicit techniques, the Wilson Theta method with Theta
1.4 has about 5% period elongation for delta-t over t
being 0.10.
Another point of interest is of course also is the
amplitude decay.
And that amplitude decay is plotted on this graph here,
percentage amplitude decay as a function of
delta-t over t again.
Delta-t over t.
And notice that the Newmark method does not appear.
In fact, the Newmark method has no amplitude decay.
The solution is right down here.
No amplitude decay for any delta-t over t, whereas the
Wilson methods and the Houbolt method is also shown here,
have the amplitude decay shown by these curves.
Notice that the Wilson method has less amplitude decay then
the Houbolt method.
The period elongation for the Houbolt method was also shown
on this view graph here.
So the Newmark method does have period elongation, but it
does not have any amplitude decay in this solution.
I should also mention this, of course, as one important
point, I forgot to mention that I'm looking here at the
Newmark method with delta 1/2 alpha equal to 1/4.
This is the constant average acceleration method.
Well, let us look then at how can we
actually perform the solution?
How can we select an appropriate time step delta-t?
Here we notice that the dynamic load factor that we
are knowing that we can derive for the solution of this
differential equation here--
here I included damping--
the dynamic load factor plotted along here as a
functional of p over omega.
p being the exciting free frequency, omega being the
frequency of the system.
The dynamic load factor looks as shown here.
If we have no damping, you're talking about this curve.
Now notice the important point.
The dynamic load factor, of course, is 1 for static
response, and it really gives us a magnification of the
static response.
In a static response, we would have this solution.
And if the p over omega value is say 2, we would have that
maximum solution in the dynamic response analysis.
If p over omega is equal to 0.5.
we would have, for zero damping, this maximum solution
in the dynamic response history.
Or rather this would be the maximum dynamic response that
you would see in the dynamic solution.
This is how the dynamic load factor is defined.
Of course, at p over omega equal to 1, we have resonance.
In practice, we have always some damping, and if you look
at this amount of damping, this would be the maximum with
p over omega equal to 1, we would have this as the maximum
response that is ever measured in the
dynamic response solution.
Well, if we look at this graph, we recognize of course
that as p over omega goes against 0, in other words as
the exciting frequency, as over the frequency of the
system goes against 0, we really talk about static
conditions.
As p over omega is very large, we are talking about no
response whatsoever.
And the important dynamic response really lies somewhat
in this region here with that being the tail end of it.
But the important dynamic response really lies in here,
that we have to include in the solution.
Out here we have basically a static response.
Let us look at two simple cases.
For p over omega 0.05, we can see that the static response
here shown in a dashed line, the dynamic response showing
as a solid line.
They're basically the same.
In other words, a dynamic response being equal to the
static response, practically equal to the static response.
We are down at this range here.
As another example, if we have p over omega equal to 3, we
are down over here.
Then the response, the static response is shown here in the
dashed line, and the dynamic response is shown
as the solid line.
Very little response here at all.
Because as p over omega goes to infinity, there is no
response for the system at all anymore.
Basically, what we're saying then is that the frequency of
the applied loading is so fast, that the
structure can;t just move.
It just cannot move at all.
And we don't have any response in the structure.
Well, the important point is that if we know the frequency
content of the loading, then what we have to do in the use
of the Newmark method is to calculate a delta-t small
enough to pick up the loading, of course, because that is
given in terms of this great steps or this great input
points to the computer program, and the time step has
to be small enough to integrate the dynamic response
of importance accurately.
The static response, the static response is contained
in the solution anyway, because the K matrix is on the
left hand side of the solution.
So the complete solution process therefore can be
summarized as follows.
Identify the frequencies contained in the loading, and
using a Fourier analysis, if necessary.
Choose a finite element mesh that accurately represents all
frequencies up to about four times the highest frequency
omega u contained in the loading.
The highest frequency contained in the loading is
really the cutoff frequency that will give us a time step.
And that then, the time step is in fact obtained via these
consideration.
Perform a direct integration analysis.
The time step delta-t for this solutions should equal about 1
over 20 Tu, where Tu is 2 pi over omega u, which is the
highest frequency contained in the loading.
What we are saying basically is that we want to select a
finite element mesh to represent the highest
frequency in the loading correctly, accurately, and
then we want to select a time step in the data integration
solution using the Newmark method that integrates the
highest frequency contained in the loading also accurately,
or the response in that highest frequency contained in
the loading also accurately.
And this then gives us basically a step-by-step
procedure for the selection of the finite element mesh and
the times step.
Notice that the mesh selection and the time step selection
both of course hinge on the structure and the particular
frequency content of the loading.
That's important.
This is the modeling procedure that we would use for a
structural vibration problem.
We would use this time step here in the implicit
integration, because we get all the accuracy we want.
In an explicit integration, we would have to use a
smaller time step.
The sentence down here "or be smaller for stability reasons"
really refers to the explicit time integration.
In the Newmark method, we would not use that condition.
We would simply use Tu, being 2 pi omega u and delta-t being
1 over 20 times Tu.
So in the Newmark method, I blank this out, and this is
here the delta-t that we would be using.
And just to summarize once the modeling of wave propagation
problem here, here what we would do is we assume that the
wavelength is Lw.
The total time for the wave to travel past a point is then tw
equal Lw over c, where c is the wave speed.
Assuming that n time steps are necessary to represent the
wave, we get our delta-t and the effective length, which I
talked about earlier, of the finite element would be then
given by this formula.
So the procedure here is somewhat different in that if
we have a wave that looks like this say, we now select a
certain number of times steps, delta-t, to pick up the wave
accurately.
And this being the length Lw, the time for the wave to
travel past the point is tw.
And delta-t then is the number of time points, or time steps,
that are used to represent the wave accurately.
And once we have the time step, we select the effective
length of an element.
And then, of course it depends on whether we have a low-order
element or a high-order element, that will give us
then the procedure of how to apply the
effective lengths criteria.
I like to finally now summarize the complete
solution process.
This is a table that I put together for you to show how
the complete solution process for explicit and implicit time
integration can be very effectively implemented in a
computer program.
We form a linear stiffness matrix K, the mass matrix M,
the damping matrix C, whichever applicable.
We calculate the following constants.
If we select the Newmark method, these are the
constants that you would use.
Let's not go into details.
These constants are given in the textbook.
They are given in the ADINA manual, if you are using the
ADINA computer program.
They are widely published.
We simply cast them in this form.
In the central difference method, we
would use these constants.
So there is some initialization involved.
Then in the central difference method, we also want to
calculate our first displacement, delta-t U, as
shown here.
And in the Newmark method, we don't need to do that, because
we simply march ahead from time t to time t plus delta-t.
You only need the condition at the conditions at time t to
obtain the solution at times t plus delta-t.
In the central difference method, remember we need two
previous values.
We then form an effective linear coefficient matrix.
In implicit time integration, this is
the coefficient matrix.
The stiffness matrix being there, as I
showed to you earlier.
The mass matrix and damping matrix being there.
If we look at this a0, we notice that a0 is 1 over
delta-t squared and there's a constant alpha here.
The Newmark constant.
But it's 1 over delta-t squared.
That is the important point.
So as delta-t gets larger and larger, larger and larger,
this effect goes against 0.
The same holds for the a1 constant.
As delta-t gets larger and larger, this effect
goes against 0.
And this is the reason why the Newmark method really goes
down or collapses to the static analysis, if we have
very large time steps.
In explicit integration, we use this coefficient matrix.
And now for each integration step, or for each time step,
we perform the following calculations.
We calculate the t plus delta-t R hat.
And implicit time integration, the Newmark method, in
explicit time integration, we calculate this vector.
Notice at tF here, no k times tU.
k times 2u is equal to tF.
In implicit time integration, having calculated for each
time step this load step, we now use the equations K hat
times t plus delta-t u, equals t plus delta-t R hat to
calculate t plus delta-t u.
Of course, we would triangulize the K hat matrix
prior to going into this loop.
We do it once and forward.
In explicit time integration, we use this matrix or this
vector here, I should say, this vector now, and we can
calculate directly for the t plus delta-t U also.
The two equations that I used are summarized here.
In implicit time integration, this is the
equation that we use.
Notice we only need to factorize this
K hat matrix once.
We do it outside the step-by-step integration, as I
pointed out earlier.
And then this gives us the increment in displacements.
In explicit time integration, this is the equation we use.
If M and C are diagonal, it's a trivial operation as I
pointed out earlier.
With this then, we have the new displacements.
And on the last view graph, I show you the formulae that we
are then using to calculate the accelerations in the
Newmark method from the incremental displacements and
variables that we knew already.
Having now calculated the accelerations, we can get the
velocities at time t plus delta-t.
And of course, we have already our displacement at time t
plus delta-t.
In the central difference method, we use our basic
equations that we already applied earlier in the
development of the method to calculate the velocities and
the accelerations at time t.
This completes the procedure that is effectively used in a
computer program for the integration of the dynamic
equilibrium equations.
I like to just now mention that I've placed heavy
emphasis on the Newmark method in an implicit integration.
Of course, there are other implicit integration formulae
available, the Wilson method, the Houbolt method, and so on.
The Newmark method has shown, has
proven to be quite effective.
And it's certainly a good example to demonstrate to you
how the actual procedure is implemented and how an
implicit integration scheme is being used.
Other implicit integration methods would be used in much
the same way.
This is all I wanted to mention to you, discuss with
you in this lecture.
Thank you very much for your attention.