welcome to the first part of the lesson on
numerical solutions to thermal field and fluid
flow in welding this lesson is part of the
nptel mooc on analysis and modelling of welding
my name is gandham phanikumar i am from the
department of metallurgy and materials engineering
iit madras so the scope of this lesson will
be a limited to a one of the two major classes
of methods for numerical solution of the nevestok
and generalized fourier heat conduction in
the equation that we have discussed and derived
in earlier parts of this course there are
two major classes of ah numerical solutions
the first class is based on the variation
formulation and these are the methods that
are used for f e m type of calculations and
many commercial softwares which use thermal
and distortion calculations usually adopt
this method and then there is a second class
of of ah numerical tools which are based on
the control volume method which can also be
related to the finite difference methods and
these are very popular where the fluid flow
and heat transfer are going to be considered
so we are going to look up the solutions using
the second method namely control volume method
mainly because that's a method i have used
to derive the equations and its also a method
where the flux balance is going to come out
very naturally and we are going to go through
the details in a depth that will give you
an appreciation of how much attention to details
is necessary before you can simulate these
processes using a computer
and the scope does not involve you to write
a program out of this course mainly because
that would take a lot of time and effort and
i would normally would like you to consider
that option once you have finished this course
and if you are interested to develop your
own program and the objective of numerical
solutions is basically to obtain thermal and
velocity profiles in the welding scenario
and we would like to have ah the results analyzed
by knowing what are all the details that going
to the simulation as you go ahead so we have
looked at analytical solutions earlier in
this course and we have seen that there are
a number of methods that are available including
for example the popular rosenthal solutions
the analytical solutions are very valuable
because they can give you answers to the thermal
field at every location in the domain if you
wish and which means that the solutions can
be obtained as smooth curves and you can see
those plots for example the violet line that
is shown in this plot could have come from
one such analytical solution however the limitations
are already discussed earlier analytical solutions
are subject to a number of limits for example
they do not take into account the variety
of heat sources that are possible in welding
and the variety of heat removal processes
that take place and also fluid flow completely
avoided in analytical solutions because it
is practically impossible to derive analytical
solutions in welding taking the fluid flow
into account
so therefore if want a more realistic solution
for welding then you must know you must go
to the ah numerical solution procedure and
usually when we pick the numerical solution
we would like to take them at discrete locations
within the domain the locations can be chosen
a priori or adapted to the solution that is
emerging in the domain but we must know that
we do not have the solution available at every
single location in the domain but at discrete
locations and this is illustrated in the plot
temperature verses distance across weldment
you can see a curve over which there are some
points that are given you could thing that
the curve may represent either an analytical
solution or the exact solution that is actually
prevent and the circles in blue are represent
in the numerical solution that we would like
to generate as part of our solution and the
difference between numerical analytical is
that the analytical solutions can give you
smooth curves while as numerical solutions
give so solutions for where the features more
comprehensive
the ah on going discussion in this lesson
and the following lesson to cover the numerical
solutions will be referred by ah this book
suhas patankar ah book on numerical heat transfer
and fluid flow this book is very important
in the area of the control volume method to
solve the fluid flow equations and i would
strongly recommend you to have a copy of this
if your planning to go further into this particular
subject and this book is also not very thick
so you should be able to go through that quite
soon and its also going to give you hands
on experience on how to write the program
because the expressions that are used are
readily programmable in a language such as
protran and there is indian addition also
available which is not very expensive from
on a books the outline of my lesson today
is going to cover the following aspects governing
equations we will just ah refer to the equations
that we have derived earlier and then convert
them to a form that is generic for numerical
solutions and then we will see how to discretize
the governing equation each of the terms will
discretized and we will see how we can use
it to write the differential terms as ah properties
of different locations and then we would see
how to interpolate various parameters at intermediate
locations and that's where a lot discussion
will come when we come to the addiction term
in the next lesson and then once we have interpolated
and written the differential terms as discretized
then we will able to obtain what are called
a set of linear equations which could then
be used to solve to obtain the solutions
so the solution method also will be discussed
and finally we will summarize our ah lesson
in the end
the governing equations we have seen earlier
the one i am writing here is for the thermal
field you can see that in this governing equation
we have got varies terms we have got already
gone through this so we have already seen
the governing equation that we have derived
in the earlier lessons and i am showing it
to here for the thermal field the first term
is the transient term the next three terms
are the addictive terms and then on the right
hand side you have the first three terms referring
to the diffusive term and the last term is
the source term and this source term essentially
will take into account the latent heat if
the phase change is being considered as part
of weld modelling and the expression is for
example delta h f where ah delta h f is the
latent heat of fusion and dou f l by dou t
where f l is a liquid fraction so this equation
ah is coming directly from the generalized
fourier heat conduction equation which we
are applying for welding and we also have
derived ah by analogy the fluid flow equation
the namely the new stoke equation applicable
for fluid flow in the fusion zone of the weldment
and here also the first term is the transient
term the next three terms on the lift hand
side are the addictive terms on the right
hand side we have the first three terms referring
to the diffusive process of a momentum diffusivity
given by mu by rho and then we have the last
term which is basically the source term and
this equation is written for the u component
of the velocity which means that we will have
two more such equations for the w component
and the v component and a set of these three
equations will then completely described the
fluid flow in the fusion zone
and the source term is here going to have
various terms depending upon the phenomena
we would normally have the pressure gradient
term where the ah gradient will be in the
direction of the component of the velocity
for which we are writing the equation and
then we will also have for example a body
force terms referring to for example thermal
buoyancy or solutal buoyancy so one example
is given here for thermal buoyancy rho bar
the average density g beta t which is the
expansion coefficient and then t minus t ref
where the t ref is the reference temperature
normally it is chosen as the melting point
so we have the body force terms coming in
as part of the source term s u and we also
have for example source terms coming to handle
the change of the phase from liquid to solid
so we have already discussed earlier that
we wanted write any equation that is valid
from the entire domain so that we can have
single domain equations and to handle that
to ensure that the velocity would go to zero
in the solid we wanted to use the porosity
approach because it is very elegant and here
is the a term that refer to ah from that approach
and is would be also coming as part of the
source term and here ah the epsilon i have
written as it is from that formulation it
is a same as the liquid fraction because that's
the fraction through which the liquid can
flow so if you notice these two equations
they are very similar and therefore we could
write them as very generic form here the generic
form is written with phi as the variable so
that this can take different ah forms for
different component to the velocity so you
could see that the phi component term referring
to the transient term the second three are
for the addictive and then on the right hand
side we have diffusive and the source terms
and for different ah fields ah thermal and
the fluid flow what are the forms that the
gamma and s will be taken is given in this
particular ah table
and we can write the same equation in one
dimension as follows so you could ah dropout
the terms for the other three other two directions
and you can see that this is a form that ah
can be att attempted to be solved numerically
and we will see it by term by term so that
we could get take first initially the diffusive
term and then later on the addictive term
and what we mean by discretization is a following
process we want to basically divide the domain
into several control volumes if you want to
take a one d domain essentially we want to
divide a line in two segments and if we to
take the two d domain we want to divide an
area into squares or rectangulars and then
if you are taking the three dimensional domain
we want to divide them into cubes or cuboids
and each of these elements ah will then be
analyzed for what is a flux that is going
through their phases and we want to also identify
at what locations with in each of this control
volumes is any variable phi known so is it
that the center of the control volume or is
it at the control of the phases of the control
volume is something that we need to decide
and we will see that both will be used for
different parameters as we go along and then
whenever we want to express the value of the
parameter phi at any location other than the
location where it is specified then we need
to interpolate and the way we interpolate
will be borrowed from the taylor series expansion
so essentially we will be knowing the var
variable phi at a given location and then
we will expand the phi as a smoothly varying
function around that location and then use
only the first order terms to see what would
be the value at a neighbouring location where
we need interpolate
and then we will then ah use this expansion
to derive a various slopes of this parameter
phi a first order slope and the second order
slope namely dou phi dou x and dou square
phi dou x square so this is how we will be
able to write the different terms of the one
d convection diffusion equation which we have
seen earlier so this is where we are showing
you how we do the discretization so let us
say that the variable phi is specified at
the location two and we want to find out what
would the value of the phi at location one
or at location three so what we can do is
that if it is known at location two then in
the vicinity location two we can expand the
variable phi as if its a smoothly varying
function using the taylors expansion the first
expression shows you how to do that so essentially
we can see that the value of the variable
at the location where it is specified that
is ah given as phi two minus delta x dou phi
by dou x evaluated at the location two plus
half into delta x square dou square by dou
x square evaluated the location two so everything
is evaluated at the location where the variable
is specified and you could do that to expand
in the left hand side and the ah right hand
side directions to show the distances as minus
delta x and plus delta x now once you have
these two equations then you can add them
and subtract them to get different terms and
that's how we have written the slope and the
ah slope of this slope namely dou phi by dou
x and dou square phi by dou x square evaluted
at the location two as a function of the value
of the parameter phi at the neighboring locations
and this the same as central difference method
and you can see that the slope is given by
three minus phi one by two times delta x if
you were to write it only with the value at
phi and phi one then you would call it as
forward difference and if you were to write
it in terms of phi two and phi three you would
cal it as backward difference so what we have
written here is the central difference method
which is known to have better accuracy then
the forward or backward difference methods
and the second differential is then written
in terms of this slopes at either ends of
the control volume phase so you can write
it as phi one plus phi three minus two phi
two by delta x square
so now that we have the differentials written
in terms of the value of the parameter phi
at discrete locations one to three etcetera
which are then given in the domain then we
can go further to see how the differential
equation is going to look like so the generic
equation is written here and then we want
then apply it for a specific case so to illustrate
we will choose the case of only heat conduction
so that we can just take terms on the right
hand side first so we are taking the diffusive
term and the source term and then we are basically
trying to solve this in a numerical manner
without considering the addictive and the
transient term so you can say that what we
attempting now is steady state one d heat
conduction ah with source term using a numerical
so let us see how the method would evolve
so we are basically looking at this equation
and we will see that p is used always to indicate
the location at which we want evaluate and
the control volume is having the faces on
the east side and the west side detonated
as e and w the vertical lines are showing
you the control volume faces and this is a
control volume over which we are going to
essentially integrate this particular differential
equation and when you integrate you already
have a differential so you could see that
integration of the first term will give you
k dou t by dou x at east face minus west face
that would be the integral of the first term
the second term will be integrated in this
manner and the source term could be constant
within the control volume or it may vary in
either case we would be taking an average
value with in the domain so that the second
term here for the source term integration
can be taken as an average value multiplied
by the delta x that is the width of the control
volume so we would do that now and while we
do that we also need to implicitly assume
how would the a parameter ah phi or in this
case temperature vary across the control volume
ah and beyond so what we are doing basically
is a piece wise linear variation which is
reasonable and this is different from another
assumption that is possible which is basically
to choose that the temperature is constant
in the entire control volume we are not choosing
that assumption because that would lead to
a sudden jump in the temperature across the
base of the control volume we are not doing
that what we are doing is that we are saying
that the temperature is assume to be specified
at the center of the control volume and across
to control volumes it is varying in a piece
wise linear form so that then we want a temperature
at intermediate locations we can do interpolation
linear interpolation and we can get the temperatures
at those locations
and if we do that then what we can then assume
is that dou t dou x can then b given by the
talyors series expansion that we had discussed
earlier so that a linear approximation is
possible and you could then look at it like
this the slope at e is given by ah capital
e minus capital p that is temperature at this
location minus temperature at this location
divided by the distance here which is delta
x e so you can see that the first term is
ah discretized to be just difference of the
temperatures divided by the distance similarly
the second term also has been discretized
and then the integral of the source term is
then taken as an average value of the source
term across entire control volume multiplied
by the width of the control volume now once
you have this then you can then proceed to
gather the terms and one small clarification
here in case this source term is a function
of temperature then necessarily you would
like to linearise it in this form s c plus
s p times t p ah the s p doesnt indicate evaluation
of the source term at the location p it only
indicates the coefficient of the source term
with respect to the temperature so this is
basically what is called as linearised source
term which means that if your source term
is going to be ah non linear function of temperature
we need to make it linear we have discussed
this briefly in one of the the lessons earlier
and here we are seeing the reason why we need
to get and usually you have the s p having
a negative value so that the temperature evaluation
will be ah stable and we would see that it
is helping in one four rules that we will
be talking about for the coefficients of the
discretized equation that we are writing and
when we gather the terms this is how it would
appear you would write it as on the left hand
side the temperature at location p t p and
it is depended upon a weighted average of
then temperature at the neighboring locations
and then there will be constant term that
will be added which takes into account the
source term and you see that it is a essentially
coming as a summation of neighboring location
temperatures weighted by a coefficient and
what are those coefficient it has been seen
from the previous slide that it comes as basically
a ratio of the thermal conductivity and the
distances and this way we can now see that
at any any (( )) location we can evaluate
what would be the temperature as a function
of the neighboring locations the coefficients
are expanded here you can that the coefficient
of t e and t w are very simply given by the
thermal conductivities and the distances on
either ends of the ah locations where we are
evaluating where as the coefficient a p is
given as summation of the all the coefficients
minus the source term which means that in
the case where you don not any source term
the coefficient a p will just be the summation
of the coefficients a e and a w which is again
a very important role it will help us validate
the discretization process for more complex
equations that will come across later on
and there are four rules that we would be
ah at ah we would be addressing ah when we
write this equations we must understand that
ah control volume method is based upon the
flux balance which means at the face of the
control volume the flux that is arriving should
be the same as what is leaving which means
that we must definitely use the same mathematical
expression to write the flux that is coming
in which is balanced by the flux that is going
out and ah the piece wise linear relationship
we have seen is going help in that in case
we chose for example piece wise parabolic
variation of temperature then it will not
work out to be the same and you may actually
have ah some numerical errors that will be
coming up so it is very important to ensure
that mathematically when we are write flux
at the control volume faces the expression
is looking identical for both faces both sides
of the face for two adjoining control volumes
and all the coefficients must be positive
you can just verify that from the expression
for the coefficient you can see that it is
given by thermal conductivity which will be
positive divide by distance which is also
going to be positive so which means that all
coefficients must be positive this is also
going to help in evaluating what would be
the ah problem in case the discretization
by using various expressions is going wrong
for a more complicated situation and we must
ensure that all the coefficient should be
positive and then the negative slope of the
linearization is also to be ensured and if
that is ensured then you can ensure also that
the sum of the neighboring coefficients will
be equal to the coefficient of that t p in
the absence of the source term and once we
follow this rule then numerical scheme is
guarantee to be stable and this is also a
check to see if there is a problem with this
scheme or with the solution process whenever
we get absorbed results out of the numerical
competition so we must first ensure that these
rules are followed and later then see where
could the problem b in case the results are
not coming out to be nice
the interpolation scheme here is expanded
here in a bit more and its not as if we can
take every parameter in the domain and then
make a linear interpolation for example take
the thermal conductivity let us say the thermal
conductivity we have the choice of using interpolations
using either linear or harmonic and here the
argument is as follows let us say that the
control volume to the left of p is having
a thermal conductivity which is very poor
compare to the control volume in the right
hand side then what would happen is that we
should not have much of heat flux coming from
the left to the location p however if you
were to do a an averaging of thermal conductivity
then the thermal conductivity at the location
w here ah would be average between the two
which means that if the thermal conductivity
on the left hand side is zero then you would
still have a finite thermal conductivity at
face w which means that would be some heat
flux that is going from the left to the right
which should not be possible in case the thermal
conductivity at w is zero and this kind of
a situation can be avoided if you use actually
harmonic mean and that where we are actually
showing you here the thermal conductivity
can be taken harmonic mean which means that
the way you ah interpolate the thermal conductivity
is to take this kind of an expression ok so
one must pay attention to evaluate in the
variables at different locations in the domain
before we directly use the interpolation which
is linear which is default for most of the
parameters and this can be done for any location
however when we come to the boundary we see
that we need to treat the boundary little
bit specially then reason being as follows
for any location we saw that the linear equation
we wrote out of the discretized form is going
to contain the neighboring values for example
for p we saw that the expression will have
the values involving w and e however when
you come to the boundary then we see that
for the note that is exactly on the boundary
you have an neighbor on the right hand side
but there is no neighbor on the left hand
side which means that we will not be able
to use the same expression which means that
we need to treat it separately and the way
we treat it is like this there are three ways
of doing that essentially we must have the
value of the boundary known ok the numerical
solutions requires that the value of any parameter
phi in this case temperature is known at the
boundary and very often the boundary temperature
may be specified for example if you have a
large plate that is being welded and that
domain is also a small portion of plate then
the boundary temperature may be given as the
ambient temperature itself so in such situations
there is no problem the boundary condition
is well define as a value boundary condition
however if that is not the situation and if
you have for example heat loss that is specified
at the boundary and if it is a constant heat
loss which is the second condition or for
example if you have ah ah heat loss that is
given as a function of the boundary temperature
then we need to do ah some more treatment
and we will see that in a movement and essentially
ah we have to now agree that we will not solve
the linear size linearized set of equations
on the boundary we will use only the interior
points for the solution and then we are going
to use the boundary values separately so in
other wards when you write a program you will
have two routines one routine that will be
solving for the ah values of the variable
t at all the interior points which will be
basically solution of the linear set of equations
and then you would have another routine where
you will be evaluating what would be the boundary
temperatures and that would be discussed now
this is now we can do there are three ways
of doing the boundary conditions one way is
the boundary temperature itself is specified
and there is very simple and it is also applicable
for welding where the domain is much smaller
then the actual welding plate
the flux at boundary is known and this is
possibly in a situation where you have a constant
heat plus that being either given or removed
from the boundary but most popular condition
in welding would be the third one where the
flux at the boundary is given in terms of
the boundary temperature and this is valid
for for example convicted heat loss on the
side wall and on the top it is also valid
for radiative heat loss on the surface and
in all the situations the boundary temperature
is playing a role and at the end of the analysis
we need to have t b known so that for ah solution
of the linearized set of equations for the
location two on wards you have the neighbors
known and the way we handle the second and
third condition is as follows when the constant
heat flux is given what we do is we balance
the heat flux at the boundary and obtain what
would be the value of the temperature at boundary
namely t b so the flux balance is given here
as follows the flux that is coming is q b
and the what the flux that is leaving is given
by the ah fourier heat condition first law
and that is given here and the source term
is evaluated at the boundary as follows then
once you gather these terms so from the discretized
from of the a flux balance then you would
see that again it is given as a b t b is equal
to a i t i plus b i is the interior point
at which the temperature will be known so
which means that you can use this equation
to solve for t b and all the coefficient are
also available so which means that we have
a way by which we can find out the value at
t b and once t b is known then our solution
process can start
in situations where the flux at the boundary
is given as a variable with respect to the
boundary temperature for example convict to
heat loss or heat gain then you would have
the expression little bit more involved and
you could write the flux balance as follows
you could write it in this manner and the
moment you then substitute what is the expression
for q b into this ah flux balance you would
see that the coefficient have changed their
form so you could compare with the previous
slide look at the value for b for example
you would see that the v b form has changed
and your seeing that the for field temperature
t infinity is also coming into the ah solution
of t b and you can use this expression to
evaluate what the t b and once we use these
methods essentially whether it is constant
heat flux or variable heat flux or constant
temperature we are finally achieving what
would be the boundary temperature calculated
separately once that is known then we can
go and try to solve the equations for the
interior points
and that then written as a set of linear equations
you see that the linear equations are written
always as ah weighted averages of the neighboring
temperature in this case for example a i b
i and c i are the weights for the averaging
and then a source term that is given in the
end as d i so all the ah equations are going
appear in the same form and these are going
to solve from location one to n where n is
number of control volumes that we have divided
the domain into but we are going from two
to n minus one because i is equal to one and
i is equal to n refer to the boundary and
these are then obtained from the boundary
condition separately and how are we going
to solve this equations if you were refer
write the these equations as matrix then you
would see that you have a matrix in which
only three diagonal rhos are filled and therefore
you could use ah for example a tri diagonal
matrix algorithm this algorithm essentially
is in to stages you have what is called the
forward substitution and then backward substitution
or backward evaluation so ah what we are going
to do is as follows essentially we are going
to start from one end of the domain and because
on left hand side we have the known value
then you can use the equation to find out
what is on the right hand side and then you
can then proceed to go further by one step
and again you have got the left hand side
value know we can find out what is a right
hand side so that is the way we are going
to do we are going to write this equation
where the p i q i are going to be evaluated
and then once we go through the entire set
of equations then you can start evaluating
from the boundary to get what be the temperature
and therefore by looking at ah the sequence
here you can see that from this linear set
of equations by substituting you can get the
values from n minus one on wards back wards
and you can then find out the values of t
n so directly we obtain the values and you
do not need to for example ah invert the matrix
to get the solution inversion of the matrix
is possible in case you have one d problem
you can definitely do that but then you have
a partially filled matrix then it is not efficient
to invert the matrix because it would be computationally
intensive and tri diagonal matrix algorithm
memory efficient and it can give you the answers
directly without having to have matrix inversion
so this is one algorithm that is freely available
also as a source code if you want to barrow
from the internet and you could also use these
equations to directly program it if you are
interested
and ah once we have this then we can see what
changes will be required in the coefficients
in case we are interested in unsteady conduction
which means that we are now adding one more
term and instead of the source term we are
now trying to add the unsteady conduction
and the ah the difference between unsteady
condition verses conduction at steady state
it is as follows we have basically a temporal
variation which means that you now have that
temperatures stored at different time steps
and at time t is equal to zero we essentially
have the initial condition that are being
stored and the previous time setp values are
all given with zero as the superscript so
that we can identify them and essentially
you have the same kind of a discretization
that is going to be applicable and when you
do that and gather terms this is how the equation
would look like and you have the variable
f that is being used here with the particular
intention this f essentially will tell you
whether your going to use for any temperature
at the neighboring location whether we are
going to use it from the previous time step
or from the current time step so f is going
to determine that and we are going to then
talk about that in a movement the coefficients
are then evaluated in the same manner as described
earlier this f for example if were to use
f is equal to zero and you can look at the
condition here if f where to be zero then
the coefficient of a e is not the current
time step temperature it is a previous tempted
temperature that which means that we are going
to use known values of neighbor temperatures
from the previous times to evaluate the current
location temperature at the current time step
and this would be called as an explicit scheme
the reason being that in the scheme you are
basically going to ah evaluate temperature
where everything on the right hand side is
known and if f is equal to one you would see
that the previous time step temperatures of
the neighboring locations are not used you
have only the current time step ah temperature
that are being used for the neighbors which
means that it is called as an implicit scheme
which means that all the variable on the right
hand side of the equation are also at the
same time step and are unknown it requires
that then we have to iterate this equation
to find out the values
so f is equal to zero would give you the explicits
scheme f is equal to one will give you the
implicit scheme and in this expression if
you were to use a f is equal to zero point
five then you go into a special scheme called
as crank nicolson scheme where we have a fifty
fifty as usage of the previous and the current
time step temperatures and this can lead to
some problems in some schemes but generally
it is a good mixer between the explicit and
implicit schemes and this equation can then
be a simplified for explicit and implicit
as follows and you see that in the explicit
scheme every thing on the right hand side
is available from the previous time step temperatures
and therefore you can directly obtain what
will the temperature at p in the current time
step in the implicit scheme everything on
the right hand side is actually also at the
current time step and therefore you have to
iterate and how do these three scheme compare
with each other so here we are plotting what
would be temperature at location p and you
can see from tem time t to t plus delta t
how would the variation be considered ah explicit
scheme would mean that until we reach the
current time step we are assuming that the
variable is changing taking the same value
and crank nicolsone scheme is actually something
in between and you would see that explicit
scheme is when f is equal to zero implicit
scheme when f is equal to one and crank nicolson
scheme is basically showing you the variation
across the two discrete time steps very smoothly
going and these three schemes have to be ah
discussed before we adopt one of them for
our solution and generally one would use implicit
scheme in this problems because you do not
have the time step problem ah whenever you
want to do a long time simulation
and the reason why that is happening is as
follows explicit time scheme is used whenever
you want the ah program to be very simple
that is the source code writing is a very
simple exercise and explicit scheme can give
you simple expressions and you do not have
any iteration of the loops therefore the programming
part will be quite easy however it is subject
to a numerical instability whenever you have
the time step that is too large what is the
upper limit of the time step that you can
choose is given by neumann stability criterion
that is given here the smallest time step
that you can choose below that you can choose
but to the maximum times that you can choose
is given as follows square divided by k by
rho c p which is basically the thermal diffusivity
which means basically we are looking at how
much time is available for the diffusion of
heat across a control volume of a width delta
x and that is a time that you can choose as
a maximum time step you can choose something
less then that and as you can see the delta
x it would be very fine at the center of the
weld pool where you are going to the heat
source which means that delta x that you are
going to use for calculations must be the
smallest grid spacing that will be used in
the simulation and the thermal diffusivity
is quite large and is going to denominator
which means that the delta t is going to be
quite small and this must be also kept in
mind whenever you are using ah further changes
in the heat source and we will come to that
in a movement implicits scheme would not have
any such limit you could use any time step
that you wish to use however ah if you use
a large time steps you may have to do more
iterations and ah the previous the time set
values can always be taken as initial guess
in the case of implicit schemes where as explicit
scheme you always have the initial condition
to start your calculations and so a time marching
scheme is available ah so some ah care has
to be taken in the choice of the time step
and this is where we are talking about for
example let us look at the grid spacing itself
you know that the heat source in the case
of welding is going to be uniform it is going
to be focused at a particular location so
wherever it is being applied they you must
have find grid spacing so that the variation
of the the heat source is captured properly
which means that the grid spacing should not
be uniform across the domain it must be fine
at the locations where the heat source is
applied and it can be course away from the
heat source
and when you then use the finest grid spacing
to find out according to the neumann stability
criterion what would be the largest time step
that you could choose then you must compare
that time step with the pulsing time step
which is used in the heat source so very often
the welding is used in a pulsed current mode
and then the pulsing is done at a particular
frequency and that would also set a particular
time step for the change of the heat source
so the time step you choose in the simulation
must be smaller than this pulse so you can
say that it could be smaller then either the
neumann stability criterion or the pulse time
step whichever smaller than ah there so it
means that the choice of time step is very
very important otherwise you may actually
loose out information with respect to the
change of heat source during the pulsing itself
and this pulsing can be as fast as fifty hertz
and fifty hertz would mean twenty milliseconds
is the time step that you would choose and
twenty milliseconds is a very small time if
you where to do a welding simulations for
several seconds or several ah minutes ah and
which would mean that the total number of
times sets for which you have to calculate
would run into several thousand or even tense
of thousands
and how do we go about programming programming
exercise is going to involve you to have arrays
to store the temperature at various locations
so you would have one array for the temperature
one array for the coefficient at the p location
that a p one array for the left side coefficient
and one array for the right side coefficient
in case you use a constant grid spacing then
you would not need an array to store the coefficients
of the linear equation however as i mentioned
to you in welding we normally have a non uniform
grid spacing so you definitely will have to
store the coefficients of the neighboring
temperature in arrays so would have as many
arrays as the number of neighbors and then
you would also have an array to store the
source term which is definitely going to be
location dependent so essentially you must
have a program in which you have multiple
arrays to store this variables and in case
you are going to use unsteady state solution
then you need also arrays to store the previous
time step values ah and you normally use the
subscript zero to indicate them and which
means that you would have more arrays that
would be required if you are going to go for
a unsteady state conduction problem
and normally you would have one sub routine
or a function where you would be setting the
initial values of temperature or the guess
values of the temperature to start with and
also to set the property values a t various
locations if you have location dependent the
properties then you may have to have arrays
for the properties also and if do not have
that problem then you can actually use only
just the ah constants that to use at in a
program in case you are introducing temperature
variation a properties then because temperature
is also varying at each location then you
may have to again use arrays to store the
property values and you have to update the
properties at every time step as and when
the temperatures are change and then you will
have to have separate routines or functions
to write the boundary conditions and the tri
diagonal matrix algorithm solution and you
can repeat these solutions at each time step
as many time steps as you would need and after
solution you must copy it to previous time
step and then update the time and then you
can keep repeating the cycle as many times
as a total amount of time for which the solution
is to be done and once the time is over then
you can write the output data and then use
the output data to visualize the temperature
field that is calculated and then you can
stop the execution once all the time steps
are completed
so you can see that basically there is a multiplication
of the number of operations so how big is
a problem of competition you can say that
it is as big as the number of time steps as
big as the number of grid points and as big
as for example a number of iterations you
need to do and what would one normally do
to do these programming very often that the
program are written in fortran language because
it is very intuitive to write a mathematical
expression directly you could also use c or
c plus plus languages ah (( )) python and
matlab are also being used because they can
be used for visualization within the same
environment and for one d problems surely
these are all alternative can be used and
initially values are generally read from a
file which means that a program that does
a numerical solution should also have the
capability to read the input files from a
file and the ah grid spaces are also being
specified as a file ok you could have initial
values and grid spaces coming as an input
into your program and outputs are then written
as either asci files or binary files asci
files are written when you want to have a
look at the output temperature numbers directly
yourself and if you are not interested in
that and directly use in them to plot then
you can have a faster i o when its done with
a binary file the input output operations
are always fast when you use a binary file
which means that when you go to three d simulation
this necessarily you will be using binary
input and output because the amount of data
will be writing in or reading in would be
quite large
and what kind of tools are used for plotting
you could use tecplot software or matlab software
or a free software such as g n u plot to use
the data and then make the plots out of them
so this is the overall scheme in which we
are doing it so we have ah done discretization
and we are in the linearized set of equations
and we are using an algorithm to solve them
and then that algorithm is then implemented
as a program and that program is then written
in some language such as fortran or c once
that program is written then you would compile
that program using popular compilers you have
for example the g n u compilers c r fortran
g c c or g fortran for example or you can
also use intel compilers if you are using
intel platform so on so once you compile then
an executable is created which can then be
run so that executable in the case of a unix
linux environmental will be called as a dot
out and they can if you if you execute that
then it would read the input parameters from
a file as you have specified and those parameters
can be changed for each execution and once
the program is executed it would be write
in the output data and that output data can
be a binary format or an asci format and that
output data format can be there visualized
using a software such as g n u plot or matlap
and now you have got the temperature variation
available so you can see now how the entire
layout is done so you start from a differential
equation you discretize that and then you
write a denest of equation then you find the
algorithm to solve them then you write a program
to implement that algorithm compile program
you execute the program and then you write
the output data out into a file and then you
use that file data file to visualize and then
you can start seeing the temperature variation
on the domain as intuitively plotted for example
using a color map for example ok so with that
we close the first part this lesson in the
second part we would take more detailed aspects
of the numerical implementation and we will
continue ah shortly
thank you
