So, today we will study stability of bioreactors which stressed up on it a little bit in the
last class, but we will start of fresh today.
So, if you go to this screen so, just hope such process is called stability analysis
of bioreactors by which we test the stability of the bioreactors. Now, the method that we
used is called the Liapunov method of stability analysis of chemical reactors. It was started
in 1955 by Bilous and Amundson, and the way to do it as I said, is to start to be able
to write the unsteady state equation given over here del C del t equals f unsteady state
equation of the system. Now, then first thing we do is actually find
the steady state of the that, we will see here we find the steady state f C s s , P
equals 0. So, that is the steady state so, why do we find the steady state because what
we want to find is that, we want to give some small perturbations into the system and small
perturbations are always available in the system because in a real system, you all know
that we discussed and we want to see the how is the system behaves with respect to this
small perturbation. It has the do the perturbation increase with
time do the decrease with time or die wave if they decrease with time and die die wave
then this system is stable, if they stay with time and increase over time then the system
is unstable. Now, what we have to do is we have to write everything in a vectorial form
over here as the we have written. So, that to retain the generality of the process or
the analysis and we write it as d d t of C vector equals f vector which is a function
of C vector mouse. You can see C vector and P vector. So, C vector
for the chemostat as we have discussed in the previous classes consist of two variables.
One is the substrate S and the other one is the cell X and the P vector is a vector of
parameter which consist of the following elements D, the delusion rate Y the yield ratio of
X over S, k s the Michaelis kinetic constant mu max is a maximum growth rate and S 0, is
a initial concentration of the substrate.So, these are the variable that are concerned
that P consist of now, x is define as a deviation from this steady state. Why?
Do we define deviation from the steady state is because we want to study how this deviation
would change? As you keep you know as you change the steady state, as you go away from
the steady state. If this deviation is going to increase then it then with time then, if
the system is unstable if this deviation x of t is going to decrease with time then the
system is stable. So, that is the basic definition that is what we are aiming towards that all
points of time. So, x of t why is that x of t because x is the deviation of the function
of time and you want to study the dynamics of it as a function of time. So, that is written
as C vector as a function of time minus this steady state.
So, what we did we realize in the last class is that we have to do a Taylor series expansion.
So, C t is written as C s s plus x. So, C s s we have to do a Taylor series expansion
of C t about the steady state C s s and so, what we do is we take the function f which
is the function of the steady state of C t. Now, C t is written as the steady state plus
the deviation and a function of the parameters. Now, if we are to do a scalar, Taylor series
expansion for a certain scalar g this is the way we write it g C s about the steady state
let us say, g C s C s s plus x is g C s that is zeroth order term, then you have the first
order term g del C evaluated at C s s. I reminded you the other day that this is not g del g
del C s s, but this is del g del C evaluated at C s s and then del two g del C square evaluated
as C s s time the x square over to vectorial plus the higher order terms, which be del
three g del C cubed x cubed over three vectorial and so on.
Now, for a linear stability analysis first question is that why do we do a linear stability
analysis? We do the linear stability analysis because we at this point are not interested
to find how the system behaves, far away from the equilibrium linear stability analysis
is valid. How valid for systems close to the equilibrium and we are trying to find out,
how the system is going to behave close to the equilibrium on the application of small
perturbation or noise. The reason we are not interested in how the system, behaves far
away from the equilibrium using a Taylor series expansion is because far away from the equilibrium
this Taylor series expansion is not valid. Why is that?
Hm What ?
Is that is a question is Taylor series expansion valid for all values effects no.
What you think Liza is it valid for all values of x, what kind of values of x is it valid
for what do you think is it valid for all values of x ?
Yes, it is Taylor series expansion is valid for all values of x then there is no constraint,
on x for Taylor series expansion to evaluate, but if you want to truncate it. Then and then
only it valid for small values of x, do you see why because if x is large for example,
then as you go to x square, x cubed, x to the power 4 these terms are dominate. So,
leaving that those terms out would be very stupid thing to do.
So, x you know if you see look at that series over here, if x is much smaller than one then,
what happens this term for example, x square over vectorial two is very small the next
term x cubed over vectorial three is even small, x to the x to the power four over vectorial
four is even smaller. So, you can leave all those terms out.
So, that is a point so, if you want Taylor series if you want to truncate, the Taylor
series at the first order and the second order then x has to be small, but otherwise this
is like a mantra you need to remember because most of the students do not get it, they all
say that Taylor series is valid for small values of x, which is absolutely wrong Taylor
series is valid for all values of x provided, you take into account up to infinite term.
A very large number of terms say for example, if you have large value of x you know reasonable
large value of say almost close to one and 0.9999 and you took into account 200 terms.
Then Taylor series is fine it is a pretty good approximation fine. So, here we are concerned
so. What is our concerned? Our concerned is small values of x, why is that because it
is a linear stability analysis, we are doing a, b. We are trying to figure out that how
the system behaves in the vicinity of the steady state, we are trying to figure out
how the system behaves in the vicinity of the steady state.
If you are trying to figure out how the system behaves in the vicinity of the steady state.
What is the deviation from the steady state x therefore, x has to be much smaller than,
what x has to be much smaller than what, just by saying that its small does not make any
sense whenever, you say something is small you have to compared it with something else.
So, x has to be much smaller than. (C s s)
C s s so, that is my criteria if I want to truncate it over the first after the first,
order term that is right here. So, g g C s plus del g del C c a valid at s s times x
then, x has to be much smaller than one much smaller than C s s. So, which is criteria
that will take and we want to do a linear stability analysis and I am trying to impress
upon you the fact why we are doing a linear stability analysis is because we want to qualitatively
understand which direction, the system is going to go see a linear stability analysis
will only give you a good qualitative understanding near the steady state.
So, we want to say we want to figure out that this is a steady state, I have a steady state
over here I perturb the system slightly. So, how is the system going to behave when I perturb
it, just after I give the perturbation not too far in time. I mean in the vicinity of
the steady state means, not too far away in time is small time from the time you gave
the perturbation. How the system is going to behave because morning shows a day right,
you know if you give for a certain perturbation, you see the system behaves in a certain way
you kind of sort of figure out from that, how the system is going to behave on the long
term and we will look at that you know in great detail.
So, this is the scalar now what would be the vector form of it. So, we want to for the
system for the reactor are the chemostat, you saw that it is a vectorial system right
need to consist of two many parameters, parameters and two variable the substrate and the and
the excuse me and the substrate and the cells. Now, what is the problem in this system ? Our
major problem is that you cannot, it is a coupled system over you cannot decouple, the
system remember if you could decouple then you can treat them as two scalar equation,
which you cannot right because look at the substrate equation go back and look at the
substrate equation I will show you, if you need to it is a here.
This equation for example, it the both the substrate equation, and the equation for the
cells; so look at the equation for the cell it has X over here and substrate concern through
mu right similarly, the substrate equation S over here and the cell comes in through
X. So, both these equations are coupled and there is no way to decouple then, if you decouple
these equations then you can solve individually for these equations using scalar expansion,
but there is no way to decouple them. So, what we need to do is kind of divide this
strategy of how to do the Taylor series expansion for a vectorial system. I am not sure if you
have done that, but even if you haven’t we will go through it and if there are questions
you can stop me, if you do not understand anything.
So, this is a Taylor series expansion what you see on the screen is a Taylor series expansion
for a vectorial system; so f here being vector f C s s plus x and A being matrix and x again
being a vector; so plus h o t higher order term which we are ignoring fine. Now, let
us first look at the look at the orders of these. So, f is a, n cross one vector both
f a is a n cross n and x is a n cross one vector of variables.
So, as a result this is of course, n cross one so, everything is matched out both in
the left and right hand side. The n cross one thing is kind of matched out right dimensionally
consistent fine. Now, what does A what is A composed out for to intuitive I will show
you in a minute, but I want you to into it. (Jacobian)
Jacobian another function why is that because if you go back here, look at the stuff. So,
this is for a scalar its del g del C.
So, you can you tell me what the see if f consist of let us, let me write it here if
f consist of f equals f 1, f 2 do you do you know how to write it this way. So, if you
are writing a vector and you do not want to waste your space in the vertical direction,
this is how you write f 1, f 2 transpose so and x is x 1, x 2 transpose and my A is going
to be what this is very simple case. Del f 1 del X 1
Del f 1 del X 2 Del f 2 del X 1
Del f 2 del X 2 So, this my A matrix and similarly, for n
cross one n cross n matrix you can write it. So, A is n cross n or in other words f is
say f 1 to f n transpose and x is x 1 to x n transpose. So, this would be del f 1 del
x 1 up to del f n del f 1 del f x del x n and this direction, it will be del f n del
x one and this direction it will be del f n del x n right, fine with everybody. So,
this is my Jacobian this is how we define the Jacobian.
Now, I can approximate as you can see in the screen, I can approximate this as this. Now,
you have to tell me what are the, what are the assumptions involved in that approximation.
Let us put lets go back to the screen here, what is the what is the approximation involved
in that, from here to here there is straight forward I mean, it is not a rocket science
A the higher order terms are ignored. And then f C s s P is a steady states so that
0. So, we just put that over here so, f C s s P is 0 and this is ignored. So, this is
approximately equal to approximation sign is not for f C s s P, f C s s P is equals
0, but for ignoring the because we ignored the higher order term that is why we get.
So, f C s s plus x , P equals A x fine and what was if you, if you let me remind you
what was f C s s plus x , P, that was d X d t. so, essentially what we have is d X d
t equals x is it right. Now, this is what I just said, so the magnitude
of the higher order term. If this is x is very small then the magnitude of this is going
to be small, you know magnitude of the higher order terms is going to be smaller than the
first, first order right this is what I explained in the scalar case. So, because higher order
X square X cube is obviously going to be larger. So, n over here the order of the matrices
is going to be same as the number of species that is there in the system fine, if you have
two species the just that you have in the chemostat.
So, most of chemostat equation will be a two gross two matrix and which also equals the
number of equation, in the system it has to equal the number of equation because for a
system which is less for a system, which is has to be complete the numbers equation has
to equal for well defined system, the number of equations has to equal the number of species
right is that clear up to this point.
So, next we go to the next step which is that if x is small then x, but this is also, something
I discussed x bar its much, much less than the steady state the deviation is much, much
less than the steady state. Now, somebody has to tell me that I discussed this just
now, with respect to the scalar concentration. Now, how does it fit into the vector form?
So, let us look at the screen so, if x bar just now I discussed that x for deviation
for this whole thing to be valid the deviation x has to be much smaller than, the steady
state concentration, but that was for a scalar for a single species. Now, what if I have
a vector what, how do I define this it could be that, you know concentration of one is
smaller the other one is not smaller, but how do I define this you know how do I satisfy
this condition.
Each of these are smaller no that is not the way we do it in mathematically, that is not
the way we do it what you are saying is X 1 is smaller much much smaller than C s s
1, X 2 and no is its very hard to you know it will be a lot of comparisons, involve then
how. So, what is the way we do it mathematically you know you should some of you might have
taken that advanced math course. So, you should know if you had when you have a vector and
you want to compare the magnitude the vectors. So, how do you do?
What determinant vector does not have determinant. What is determinant of vector, how do you
compare magnitude of vectors you want to go back to class twelve or what just before the.
How do you differ, how do you compare magnitudes of vectors 2, vectors I give you three you
know each of them say consist of 3 variables. How do you compare the magnitudes the 2 variables,
3 variables, 5 variables whatever anything larger than one variable, how do you compare
the magnitude?
Norms, yes by comparing the norms of the vectors so that, that is the way to compare and so
there are different kinds of norms, there are there is a two norms there is infinity
norms and so on.
So, we have the two norms for example, here it will be it will be x 1 square plus x 2
square plus you know x 3 square. Let us say if you have 3 square root of that then you
can the infinity norm you can the absolute norm. What is the absolute norm? May be absolute
is be the absolute values of this x 1 plus. So any of these norms the most popular is
this one norm and this is, say I believe you know how to calculate the length of a vector.
So, this is the this is defined as a length of the vector that norm and that is the most
usually, norm that you use to compare to vectors you can use other forms of norms also like
the absolute norm for example, the infinite norm and so on, but typically what you use
is this is this norm, where you know just compare the lengths of those lengths of those
vectors take the square of them and sum, sum them up and take the square root. So, if I
am to going back to the screen.
So, if I am to decide whether as my x is small, I just compared the length of x as compared
to the length of the other vector. So, when I have done that and when that turns out to
be small then, I can make this approximation that d x d t equals A x clear. So, where a
i j that is the elements of the matrix A are given by just as, I showed you know f del
f. So, this is a Jacobian and given by this. So, del a i j could be written in the compact
form as del f i del C j is that clear to everybody any question, you ask me please. So, this
is this is the compact form of writing a i j. So, a i j is the elements of the matrix
A. Now, what I want you to sort of you know tell me if not completely, but at least to
some extent what would be the solution of this if I am allowed to go back to this.
No, no. What is the solution of this equation? Let us not go that part d x d t equals A x
what would be the solutions?
Yeah,
How do you know that. So, this is my equation del x A what is the solution ?
Eigenvector into.
Some constant some constant alpha I times.
Exponentially mean.
Hm Minus.
Yeah
Into into
Values into t and what sigma.
Sigma of that so how do you, how did you get this?
How do you know this?
See there is a easy way of doing it, which is very straightforward to do it. Now, this
is the equation and the solution of this is if X 0, X equals X 0 at t equals 0. What is
the easy way? The easy way is X equals X 0 exponential A t. Now, you have been taught
how to break up A t, if A is know A is given then any f of A it could be written as a function
of the Eigenvectors, Eigenvectors and the Eigen values. You aware of that I think its
you know f of lambda f of lambda times the Eigenvalue, Eigenvectors also coming there.
And if you use that you can go back to that other option is let us look at the screen
other option is this is the solution. Other option is to actually breakdown A x, there
when you have the matrix in terms of it Eigen values and Eigenvectors and then you can decouple
the system and solve it, solve it in this way that is like, but just coupling the system
that is one way other way is to straight away go and solve this in the matricial form and
then convert the matricial form into this form, but whatever it is this is the final
solution that you have to remember this solution and now, this is a well known solution.
So, d a d d t of x equals excuse me a matrix of times x the solution of that would be x
is given as sigma times some constant alpha i times beta i, beta i being the Eigenfunction
times e to the lambda i t you said minus lambda i t, but that is fine, but that depends only
how you define you Eigenvalue. If you define your Eigenvalue as you know lambda t and its
one and if you define it as minus lambda t that is another, but that is all right. So,
if you define A x equals lambda x then it comes out to be lambda i t fine, that is that
is how it is and if you define A x equals minus lambda even is minus lambda i t, but
typically you define A x lambda, lambda x. So, let us not confuse it so given the definition
that A x equals lambda x then, this is your solution there are two ways of reaching the
solution. One is as I said to use the Eigen values and the Eigenvectors right, away and
decuple the system and solve it there the other one is to not decouple, it use a matricial
solution form first get the matricial solution and then use the Cayley-Hamilton and other
theorem to break up that matrix say. So, exponential this is something might want
to do in Matlab or mathematic has a good thing. So, if you give your exponential if you give
a matrix in three by three matrix and you give exponential of that matrix then mathematic
will immediately spit out what. So, you put exponential so you put a three by three matrix
A input it then you say exponential of a and then mathematic will immediately, spit out
what that exponential of A is in terms of just another matrix you see what I am saying.
So, basically it is like this.
So, A let us say is two by two; so a 1 1, a 1 2, a 2 1, a 2 2. Now, exponential of a
would be some other matrix b 1 1, b 1 2, b 2 1, b 2 2, but these are related these, these
are related through the Eigen values and the Eigenfunction, Eigen functions and Eigen values
of the previous matrix. Now, give you quick problem say A equals simple
A equals a 1, 0, 0, a 2 what would be my exponential of a how do you make exponential of a does
not say or whatever three by three whatever any idea no it would be e to by a 1, 0, 0,
e to by a 2 and this is valid for any value. So, if you have diagonal matrix then any exponential
or anything that you do on that is same as doing on its individual elements is it clear
Liza. If you have a diagonal if you have diagonal matrix anything any operation that you do
on the matrix is same as doing on the individual elements of the matrix fine. So, I hope you
agree with me on this and I hope you remember this also those of you do not that, the vector
something that we will need all the time that if d x d t equals x then x is given by this
form alpha i beta i into the lambda i t where lambda i is the Eigen values of the system,
given as you know x equals lambda x and beta i are the Eigen functions. Now, how what alpha
i is a certain constant how am I going to define?
How I going to obtain those constants is through defining the initial values. So, if my initial
value is given at x t equals 0.
Let us say so my solution over here is x t equals sigma alpha i beta i e to the power
lambda i t, i going from say 1 to n then my initial value is given x 0 equals x 0 fine.
So, I put my x 0, t equals sigma i equals i going from 1 to n alpha i beta i right.
So, this is how I define my now my Eigen functions are known. So, there has to be a unique there
has to be there have to be unique values of alpha i with such that this will work out
right. So, if your beta i has you know three elements
and with each of these three deductions. So, let us say if you have two by two matrix and
beta i has three, three elements each of these. So, let us say alpha 1 times beta 1 1 plus
alpha 2 times beta 2 1 should equal x 0 1 first element similarly, you will have alpha
1 times beta 1 2 plus alpha 2 times beta 2 2 should equal the second element of this
and you have. Now, two unknowns two equations alpha one and alpha 2 2 unknowns equations
you can solve for them and get the values fine. So, you just need to invert it.
using the Eigenvector. (That is transpose of the matrix
You can solve what alpha i. (Alpha a
Eigenvectors of the transpose let us call it.
Y i. It’s called Y i.
beta i. That is correct. So, you can do it that way
also, but any how this is a system of equations with alpha i is as unknown and these two as
known. So, you can solve them any which way.
So, the Eigenvectors are solved and you know. So, the Eigen vectors and Eigen values you
know how to solve them and I do not need to repeat them, but still I will just quickly
go through. So, this is the Eigen values are solved using determinant of a minus lambda
i equals 0 and then Eigenvectors that obtained by solving a minus lambda i times beta equals
0 fine. So, these are some basic things you learned
in matrix algebra, but these are things you should be able to do because when a problem
is given to you. You have to be able to actually find the lambda i then the beta i fine. So,
and just as I said the initial condition just as I showed is shown by that. So, alpha i
depends on the initial conditions so, alpha i times beta i equals X 0.
So, there are n for n by n nonsingular matrix, there is n Eigen values and n Eigen vector
n corresponding Eigen vectors. So, which means that the system is uniquely defined if there
are n Eigen values and n Eigenvectors in the system is uniquely defined if there is a,
if there is a less than n Eigen value n minus 1 Eigen values then what is what does it imply,
if there is. Unique solution is not possible solution is
possible, but unique solution is not possible and in other words the system is ill-defined
is it is an ill-defined problem. So, essentially what is happening is that you have n variables,
but n minus one equations like we said if the system has to be completely defined then
you have n invariables and n equation. Now, if you have if the matrix is singular
then you know it what it means is that, if the for example, is less than the order then
what it means is that the system is ill-defined or in other words the number of equations
that, you have from the system is less than the number of variables. When that happens
what are you suppose, to do you are actually not suppose to settle around with the algebra,
you are suppose to go back to the model and check if there is something wrong with the
model. You are missing out an equation with the model it is there is no point fiddling
around with the algebra or the mathematic when there is a when the system is nonsingular.
Now, these Eigen values now we come towards the stability to you. Now, we are kind of
walking towards the stability equation these Eigen values that you get are complex numbers
in general, they are complex numbers then they will have a imaginary part and a real
part for example, here any lambda i am sorry this should be lambda i excuse me is written
as a i plus minus b i, b i times i. So, i b being the imaginary part and a i being the
real part and now, we will look at what the effects of these imaginary and real parts
are on the system.
So, if my Eigen values let us see so basically what is happening is that. I need to have
the solution out here. So, the what is happening is this is my solution. So, let us write lambda
i as a i plus b i times t. Now, what happens is let me write this over here and then let
us see.
So, X equals sigma alpha i beta i e to the power just write it like this a i plus b i
times i times t i from 1 to infinity. So, if this is our basic thing and now all stability
results will depend on this. Now, you tell me that if my a i is what, what are the effects
that you think of a i and b i. So, look at this what we are interested in is how it varies
with time right how it varies with time. So, what are the two factors that are going to
decide how it varies with time obviously a i and b i fine.
Now, let us try and you know try to discuss what would be the effect of different values
of a i and b i. So, a i is the real part b i is the imaginary part. So, a i is a real
and b i is the imaginary part. So, you have to tell me that what could you think would
be the effects of these. (Real part should be negative)
What? Real part should be negative.
Real part should be negative for all kinds of stability. Now, the imaginary part may
or may not be negative. So, what do you think if it is negative fine if it is not negative,
if the imaginary part is not negative then what will happen?
System will be stable, but something else will happen.
(Oscillations) Oscillations see what does this imaginary
part, you know e to the power i theta i to the power i what does this imaginary part
is composed of when cos and sin. When you write e to the power i x, what is it composed
of cos and sin. So, if the real part say a is the real part and b is a b is a imaginary
part. So, here over here if a is the real part is negative then definitely means, that
the system is stable. If the real part is positive then the system is always unstable.
If the real part is 0, what is it called you know what is it called, if the real part is
0.
No, the system is called neutrally stable that is just stable slight perturbation from
this side, slight increase in a will lead to growth or growth of noise slight decrease
will lead to d k. So, it is called neutrally stable just stable so that, that are those
are my explanations on the real part. Now, let us look at the imaginary part b i. So,
what do you think are the effects of imaginary part?
If the imaginary part is 0 it means, that there are no oscillations in the system, if
the imaginary part is nonzero that is positive or negative then there could be oscillations
in the system. So, what you have to do when you looking at the total system you have to
couple the effects of the real part and the imaginary part.
So, what we did over here is in the let us go to the screen we divided into 5 different
parts. So, first case is where the real part is negative and the imaginary part is 0, this
system is completely stable fine. Next, is the real part is positive and the imaginary
part is 0, the system is completely unstable third is real part is 0, but the imaginary
part is nonzero positive or negative, what it means that there will be oscillation, there
will be cos and sin part as a result of which there will be oscillations and sustained oscillations.
sustained oscillations means that, these oscillations are not growing with time or decaying with
time they remain what they are, but there are sustained oscillations I will show you
pictures in a minute and you will appreciate this better next, is that with a negative,
negative real part you have the imaginary part nonzero which means that, that would
going to be oscillation, but these damped oscillation.
So, oscillations are large at the beginning and then the slowly decay and the last one
is the real part is positive and the imaginary part is nonzero which mean that the oscillations
start and then they grow with time. So ,out of this how many stable a is stable first
case is stable a number 1, number 2 unstable number three neutrally stable number 3 is
neutrally stable just stable number 4, stable number 5 is unstable. So, two of them stable
two of them unstable one of them neutrally stable.
Let’s look at the pictures. Now, so this is the case first case where it the real part
is negative and the imaginary part is 0, second case real part is positive imaginary part
is 0. So, first case as you see x decreases x is a deviation right. So, you started you
gave a little perturbation started with a little deviation at t equals 0 and the perturbation
or the deviation decreases with time and goes back to the steady state.
So, this is this is your steady state value at 0 is at the bottom of the x your steady
state value, what do you see and we increased it to this extend and it goes back to the
steady state value fine. So, that is why this is stable node and then the next one is a
is larger than 0 and b is a 0. So, this you start with a small deviation let us say whatever,
the deviation value is and it grows over time. So, the system even if it is a very small
you give a large so, that is a difference for a stable system you gave even whenever,
I mean sort of significant perturbation and still it goes back to 0, for unstable system
you gave a small perturbation and that blows up and these are the most dangerous ones as
far as biochemical reactors or any kind of reactors are concerned and I think I discussed
this the concept of hotspots here. In this class at one point of time and this
is what it is related to systems, where you have hotspots you give them a small perturbation
and they growth and let us look at other things then I will come back to the hotspot question.
This is the case where the vortex point it is called neutrally stable which, where the
real part is 0 and the unreal part that the imaginary part is nonzero and you have sustained
perturbations. So, whatever perturbation you gave remains,
it does not go away and does not increase either then, this is a stable one it is called
the stable focus. Stable node with a imaginary part or with a perturbation or with a oscillation.
So, you start with an osculation and you give a perturbation it starts oscillating and then
it decrease all of this in the control theory, I believe this starts decreasing and, but
it keep keeps oscillating went through the oscillation it goes to 0 it decays, but decays
in an oscillating way.
And the last one is an unstable focus which is that you start with a little give a little
perturbation it starts to oscillate and the oscillation start to increase. Now, this is
one I am coming back to the hotspot question, this is the one that actually is most seen
in hotpots. So, there are two possibilities of instability one is the this one this is
the example the second one here, where it grows uniformly and this one where it grows
the oscillation. Now, what is seeing in most reactor than what
we are studying biochemical reactors that a chemical and biochemical reactors. What
happens when the system goes unstable, it goes unstable through this through this unstable
focus the reason being that it typically that is an oscillatory part and if it is a especially
this happens, when there is a what kind of coordinates do you think this will happen,
what kind of systems will happen in circular or spherical systems in cylindrical or spherical
systems, why because you have the radial direction over here.
Let us see in this in this cylindrical system and these oscillations happens in theta. So,
the hotpots starts to rotate so, these oscillations will actually lead to rotation of the hotspots.
So, you have a hotspot over here and these are in practical system this is what happen.
When you study the systems theoretically, we kind of assume there is a hotspot over
here and there is another hotspot emerging and another hotspot emerging and so on.
But in a real system that does not happen, what happens is that these hotspots start
to travel you know because there is a imaginary wave, there is a cos and sin part of it and
these are oscillation, these hotspots start to travel and as they travel they gain in
size and then and then it will lead to a huge hotspots and then, whole explosion the reactor
blows up and that is how reactor you know in industry for example, all these accidence
in the plants that you hear of actually happen this way.
So, hotspot means a temperature a small region in the in the reactor, where the temperature
grows suddenly and then it starts to gather you know gather in space like a for example,
if you start started over here and then it will goes over here it becomes like this it
goes over here it becomes like this. So, this is our hotspot is generated small increase
in temperature is a major problem lesser concentration temperature is a major problem.
You can have concentrations spots also you know, but concentration hotspots, but temperature
hotspots is a major problem because it lead to the reactor blow up. So, you have you start
with a small hotspot over here and because there is a if there is no see, if there is
no unstable focus if it is just simply unstable or in other words, the imaginary part is 0
then what happens is that hotspot grows in one space. So, it just stays there and keeps
growing, keeps growing which is bad also, but this is lot worse. Why it lot worse? Because
it is it propagates you know it helps in propagation. So, it goes to another regional space and
spreads it over there with that is why it is if you, if it just increases in one space
you have an idea and you can cool and through cooling and stuffs like that, but then this
propagates and as it propagates at it continues to infect other regions to hotspot, where
regions wherever there are no hotspot it continues to infect therefore, the unstable focus so
unstable focus means, why if there is a focus of hotspot. So, this is my focus of hotspot
and that itself is unstable. So, there is a focus and the focus itself
is moving whereas, this is a stable focus you know the system is unstable overall, but
the focus that is there is not unstable that is a its unstable in one place is that clear
this mathematics and how it translates to the physics of the problem and these are very
real problem, which is not it is just not something that is of mathematical interest.
You know this is these are things that happen in the reactor and you will read in a papers
here and there every now and then, that reactors blow up and this is the reason why reactors
blow up because you have certain noises or perturbations in the system that lead to because
why is this, what is the root cause of all of this is because the nonlinearity of the
reaction and as we do a problem, you know just in the next few minutes you will see
that how the nonlinearity try start to infect the system.
So, because there are nonlinearity in the system these hotspots are generated to start
with it its a completely linear system nothing like that would ever happen, you have hotspot
it will go away you know you increase the introduce the perturbation, you can try this
with a linear system introduce a small perturbation to the system, it would be decayed. You know
it will damp out and go away, but with a non-linear system you introduce the perturbation. The
system it starts to grow there are possibility that it grows and what we here or today trying
to understand is under, what conditions that will happen? Why are you trying to do that?
Because as a engineers we want to predict a priory whether there is a possibility of
hotspot or not, you do not want to have the hotspot blow up into an accident and then
you know. So, you will get fire if that happens so, you want to take care of your job and
you want to figure out how the, whether the hotspot is going to happen or not and that
is why? We are doing all these analysis I am think some of you might have done similar
things.
So, as we said that a i plus b i t is the whole thing, you know b sorry a i times plus
b i times, I is a whole thing and what is a whole thing the Eigenvalue of the system
right. So, what we need to find out is what is the, what is the nature of that Eigenvalue
is it positive is the real part positive or negative and is a imaginary part positive
or negative, but what is a prime importance for us more than imaginary part, the imaginary
part only leads to oscillations. So, what is a prime importance to us is a
real part whether, the real part is positive negative or 0. If it’s positive its unstable,
if its negative its stable and if it is 0 its neutrally stable. So, we want to do that
now, but we want to little clever in what we want to do is we do not want to compute
all the Eigen values. Why do we not want to compute the Eigen values is because if you
have a system like a two by two system for example, or three by three system is still
possible to go and compute the Eigen values, but in industries you will network of reaction.
You do not have reactions where you just with three species you know and I showed you the
examples of this multiple reactions. What happens is when you start a reactor with two
or three species a whole network of reactions emerge because you know, you cannot control
once you put the thing in you do not have control over the entire system. So, you want
a going to b going to c or a going a reacting with b producing c, but then d e and e and
f and all will come up and it will lead to a network of reactions.
Now, when is there is a network of reaction, this matrix a that we are talking of will
be a large matrix it will not be a 2 by 2 matrix n is going to be large number maybe
10 or 20 or 50 there are you know most of the industries many of the industries in petroleum
industry for example, the number of n value of n is very, very large it could be 100,
200, 500. So, how I going to 5, 500 Eigen values each
time it is very hard. So, you do not want to find all the Eigen values and that is a
lot of work. So, what you want is a shortcut to figuring out so your aim is not to find
the Eigen values you know it is we are not doing a course applied mathematics over here,
that we our job is to find the Eigen values. Your aim is to find out whether the reactor
is going to remain stable or blow up and for that you do not need the Eigen values, but
what you need is the science of the real part of the Eigen values and that is what you concentrated
on. So, what all we need to do is we need to determine
if the real part of the Eigenvalue is positive or negative. Now, so this is your equation
for obtaining the Eigenvalue a minus lambda equals 0. Now, this equation you for finding
the Eigen values could be written as this right. So, you can just open this up and write,
write this in this form fine is that clear if it say n by n matrix, then your Eigenvalue
in this characteristic equation is called the characteristic equation for Eigenvalue.
So, characteristic equation for the Eigenvalue will have n cross. You know it is it could
e to be off of order n right is that clear why because is that clear Liza why or do you
want me to go through that.
So, if you have a i 1 1 minus lambda a 2 2 minus, minus lambda and so on. Say a n n minus
lambda then, what is the order of this. So, lambda, lambda, lambda, lambda to the n and
then the other terms this is the highest order and then the all the lower order terms will
be there is that clear completely. So, when you arrange these terms you collect terms
of lambda to the n and so on. Then divide the numerator and everything and
get it in the cleanest form possible, this is the characteristic equation you get where
lambda n is given by lambda n plus b 1 lambda n minus 1 plus b 2 to the lambda n minus 2
and so on, b n minus 1 to the lambda to times lambda plus b n this is how it is given. Now,
all you need to do is find and you know what is a criteria called do you remember?
Routh root criteria and this is the Hurwitz criteria, if you remember so of finding the
two different ways of finding what the signs of this a would be.
So, the Hurwitz states that all what this is what we do and you can do the routh criteria
also it is either it is mine they states that all roots of an equation will have negative
real roots. If one is here B 1 the first coefficient here is greater than 0 and what is the other
one and if you remember the determinant of this B 1, B 3, 1, B 2 that is you go this
way 1, B 1, B 2, B 3 like this is greater than 0. Now, if this is I n by n matrix if
it is not a 2 by 2 matrix.
Then, this is how it is 1, B 1, B 2, B 3, B 4, B 5 this to see how it is going, its
going zigzag like this; so, one this reduction, then this reduction and so on. Fine it is
going like this so, if I look here let us go back to the screen; so 1, B 1, then B 2
then B 3 then B 4 then and so on. And then you start at B 3 again here and then go on
in the same direction, you just move and then go on in the same direction and then up to
B n. So, luckily for our case you know for the chemostat case, we have how many components
what will be the order of this matrix. Two
Two yes typically two by two, but it could be three by three if you know multiple substrate
or if it is a inhibitor that is involved then it could be three by three, but typically
it is two by two. Why? Because even if there is multiple substrate its different even,
if there is a inhibitor. What we do? We get rid of the inhibitor out of equate problem.
So, we come up with two by two system. If there are multiple substrate then yes you
know if s 1 and s 2 are two substrate then, how many what would be the order of the equation
three. So, it would be three by three. So, either it would be two by two or it would
be a three by three matrix for us. So, this matrix is this either if it is a 2 by 2 then
this is going to be the form, if it is a three by three then this is going to be the form
of it. So, let us remember this and you can also
prepare for the three by three and I do not think you need the generalize form, but I
still gave you the generalize form. So, these two criteria have to be satisfied and they
are they are joined by the logical which means, that both of these criteria have to be satisfied.
So, let us now go and look at the stability of this chemostat fine.
So, this is the stability analysis that we did, and this is generalized this is for valid
for everything. Now, what we are going to do is we are going to try and apply it to
the chemostat. So, the first equation for the chemostat is what we have done remember
dx dt equals dilution rate time x 0 minus x plus mu times x this is mu the whole thing,
the initial calculation that we will do we will assume it to be mu. Why would we assume
it to be mu? We own put the whole thing mu max S over K s plus x we just assume, it to
be mu, why would you do that would assume it to be mu as a function of as, but why would
you do that, what is the growth rate of this system?
What is it, what is the ideal typical growth rate of the system? It is mu times x, what
is the growth of the system its mu times. So, I think there is a basic problem here
so, let me go to it. So what is a difference between mu times x, and just you know the
concept of mu; what is the fundamental difference, growth is mu times x so, is the growth rate
is what is seeing seen on the screen over here mu, mu max S x over K s plus K plus S
is that the growth rate at all time. Yes or no, I do not understand you not if.
So, I do not know why you could not you know we are not so, what is the growth rate? The
growth rate is mu times x at all times, but the specific growth rate is not the Monod
growth kinetics at all time see, there is fundamental difference growth rate is mu times
x is all time there is no confusion about that then. Next, question comes what model
do I use for the specific growth rate and I can use different kinds of model for specific
growth rate. So, as a result when I am going doing this
analysis I want to do it for mu. So, that I can plug in any kind of growth rates at
a later time and the model analysis would be still valid is that clear because I think
that is the major thing that, you are not differentiating between the growth rate as
a whole and the fact that you are using just a Monod growth kinetics. You can use modified
Monod growth kinetics, you can use you know the Malthusian kinetics or different types
of kinetics. So anyhow, so these are the coupled equations
that you have; and so my f 1; so I need to define my two functions in order to be able
to do all these calculations. So, what I need to do first I need to define my two functions
and then I need to obtain my Jacobian obtain my A using these two functions and then I
need to look at the Eigen values of a not specifically, but we need to look at it through
the characteristic equation fine. So, first step in this is therefore, looking finding
the functions so, this is my first function and this is my second function. Next thing,
I need to do is find the Jacobian. So, why would you work on your copies to find the
Jacobian so, these are my f 1 and f 2. So, I want you to find the A matrix.
This is minus mu over Y, you got plus minus or plus.
Minus and this is minus D can we see other fine. So, these are my calculations we will
stop here and continue with the rest in the next class.
So, I just quickly capture back say f 1 and f 2 are my functions, these are my Jacobian
and I take the derivates in each step and get this.