Advanced Fortran Programming : 019 : Solving LAE using Runge-Kutta 4 Method for LAE + GNUplot Video Lecture | Introduction to Fortran Programming (AdvancedLevel) - Database Management

hello guys my name is Aaron and welcome
to my channel this series is a series of
tutorials on advanced Fortran
programming now in this tutorial we'll
be using the same example that we used
last time to solving a linear advection
equation and this time we'll be used
we'll be solving this in a very
important method that will be useful for
by many people for time marching Carla's
rk4 method or the runge-kutta method no
people's parents at different man
different plays differently in different
places so I pronounced it as Runge Kutta
method okay but I mean first simplicity
I will call it RK 4 okay now I'll walk
you guys through this code over here it
may look a little big but it's simple
I'll just walk you guys through over
here now but those of you guys who don't
know what runge-kutta method is it's a
finite difference it's a final different
scheme based design it's a mathematical
technique and a finite difference scheme
designed for you know time matching time
matching solution try matching solutions
and this is known for its highly
pressure with known for high accuracy
and when compared to the other methods
available okay but RK for us a little a
little complicated to work with
especially the mathematical derivation
of it but computationally it's not that
difficult especially with photon you can
do this much very very easily okay I'll
walk you guys through what's going on
here okay
I have this program called a solver and
will dust it will solve the linear
advection equation just like all we did
last time okay then what we do orders
has that I just time it has the 4001
time steps and the thousand and one
space points okay so this entire domain
is actually divided house inand one
space points or else has thousand the
thousand gaps in between or thousand
intervals in between okay and there are
like thousand and one four thousand time
steps okay
the game ideas in spgs for space points
and using that I'm defining a different
let's see yeah hey I'm you and be
amusing that I'm defining a defining how
do I put and defining a Matic so I'm
defining a matrix over here which is
like four one which is like oneness
which ranges from oneness to four
thousand and one for the four thousand
and two over here I think it's just four
thousand four thousand one that's it I
guess I kind of that's a typo anyway and
then and then I have am at four thousand
one cross two thousand one matrix called
as you which will store all the values
of U and then I have real value which
whose diamond which is one dimensional
matrix a thousand and one entry to the
one dimensional matrix is mattresses
namely UD k1 k2 k3 and k4 this k1 k2 k3
k4 they are actually used for calculate
used for calculating the intermediate
some intermediate values intermediate
values in the raga could have method
okay for the also because you do not you
do not know I strongly recommend you
have a look at look at on he could have
method we could have method and in
Wikipedia some other places okay some
other places and then I set the
parameter parameters here c2 b2 which is
the real it's just the velocity of the
sine wave and if you guys went through
the last tutorial of mine you'll find
the same the domain I'm keeping it as
thousand meters just like the previous
tutorial setting PI to be the PI to big
value setting the pi value here I'm
setting the time domain setting that
setting that there are like thousands
thousand seconds over here now DX and
DTR the space entire space intervals at
time intervals accordingly now here what
I'm doing is that unlike last program
this time using the number of points I'm
deciding how much per Harmon how much
should be the space in between them and
this works out this will work out to
this will work out to be nice so DX will
be around see the domain is thousand and
the space points are thousand and one
so there will be like so DX is work one
meter and DT is thousand seconds later
by four thousand so it will be 0.25
seconds so this time will be ad will be
moving the way slowly at every part
every 0.25 seconds point of a seconds
and the distance will be 1 meter so this
will be this plot will be much longer
but anyway this is the same cfl criteria
we looked at last time okay and I'm not
signing the last values of SP to be the
car to be the values of L now this by
the X to this x value will range from 0
to thousand
due to this setting due to this kind of
an iteration with I okay this will vary
from 0 to 2000 okay now here I am
setting the initial in the first time
step value i won instead of starting
from the carnival in sub starting from 0
which I used last time and I don't want
to use that I thought of you know using
with the normal notation so I started
from 1 and it will be going from what
I'm setting me from the first point to
the last point to be the sine value of
the X sine value of sine value okay sign
well twice PI sine of L so that so
thereby it starts from 0 and ends with 0
so we'll have a big sine wave occupying
the entire thousand meter domain okay
and I'm setting these mattresses k1 k2
k3 k4 to be zeros okay and I have a
format trying a format with format
number formatting
ID number 10 okay which will print to
will print two floating numbers
separated by five spaces separated by
five spaces and then here I'm calling
the I'm starting the time nursing scheme
Runge kutte a time matching scheme what
I'm doing is that I'm calling this
subroutine sub underscore lae underscore
RK 4 okay
with such that I'm passing ck1 UD dx/dt
and SP now if you guys know how wrong he
could have met runge-kutta method works
what what do you need is that you need
the value K 1 and then the initial
velocity
initial velocity that DT D DX DT and C
are kind of needed that we got that will
be part that will be part of the
function okay and then the second step
what I would have to do is that I'm
passing I have to pass you do to pass
the velocity plus half of k1 to get the
value K to okay here I pass you D DX DT
and C to get the value k1 here I pass
you D plus half of K one-half of K 1 DX
and D T and C to get the value k 2
similarly here I pass you D plus half of
k2 and DT and DX and then C to get k3
and then here I pass he here a pass do
UD + half of k3 and then followed by DX
DT and C to get the value K for okay SP
stands for the number of spaced points
for further a scheme to work you don't
you don't need all this but the reason
why I'm passing this is because this
will allocate the matrix will it it will
determine the size of the matrix
honestly now when these two these four
subroutines get called I will get will
get for these matrices matrix values k1
k2 k3 k4 will be will be set to some
nonzero values and those values
okay will be used here so the the new
value with the current velocity at the
current time step current time step for
all the points from the fire from the
beginning to the end are from all the
points will be equal to the previous
velocity value plus the weighted average
of k1 k2 k3 and k4 the weighting is such
that the k2 resist twice to 80 meters of
2 and K 365 days of - ok K 1 and K 4 K 1
and K 4 receives one unit of weight is
so when you divide this it's 1 K 1 plus
2 K 2 plus 2 K 3 plus K 4 divided by 6
so this weighted average is actually
added to Kate duty and this you D was
actually set as this one and this is the
same thing we did last time there this
is this is the same thing with the
dashed in the last point of the funk's
last point of the function in the
previous time step is made as the first
point
point of the first point of the wave
with the current time step simple as
that now now what I'm doing to do is
that I'm going to and I start just after
opening a file called as data dot that
using this open for open and sub
function and using this with a unit
number one okay and I'm writing all the
values of J from 1 to SP you want ASP
okay with into the Phi into the spell
data that that is in pharma turn this
power button is given on the top okay
and then I'm writing the values of X X
and the current time current values of
UI okay obviously I'm ignoring the first
time the first time step but it doesn't
matter but anyway it doesn't matter much
I'm just saving all the values in
parallel in along the columns and then I
close the file over here and I close
close that loop to loop over here and I
close the unit file over here and then I
call the sister I call the system I call
the system function to invoke the gnew
plot function you know give me plot over
here
alright now let me tell you guys how
this subroutine works it's very it's not
that difficult
ok it's not it's not that very difficult
what I do is that what I do is that I
get the values over here and I set the
value as such SH is actually the shape
nothing much I set the value here and
using that shape SH or which is like SPR
which is the previous thing
okay I'm setting the value of U and K
over here okay you I'm getting you as
intern in I'm getting KS intent in out
because K value will change will change
okay after modification so our K value
will change and I'm getting CDX and DTS
just see DT has just you know values
that I will be using over here and then
what I'm doing is that this is this is
the important line
what I'm doing is you're doing here is
that I'm just following the Runge Kutta
scheme and right when you task
and I'm writing the right-hand side
function of the that will defining the
right-hand side function of the Runge
Kutta scheme in this line and this is
what you have to do and once this is
done and once this is done and I mean
I'm setting all the values of all the
values of K from 2 to 2 to the last
point the first point will not be will
not be modified ok because we are using
that there isn't it linking linking
thing where and we lose the last point
of the last point of the previous time
step as the first point of the next time
so so I'm doing that so I'm ignoring one
other than that I'm just using this
using this to set the value right and
using the setter value right and if you
guys don't know how this works
iced if you guys don't know how this
works I strongly like when you go and
have a look at the runge-kutta 4 scheme
you land us a key for scheme you
understand this nicely ok that's about
this subroutine okay and here the this
addition will take care of all the k1 k2
k3 k4 complications that will come in
the that is the 4 steps that will come
in the regular scheme simple as that
so that this way kind of makes our job
easy with writing just a single line to
do the Sun Runge Kutta scheme ok now
that being said ok that being said now
let's look at the gain no plot file I
have it's nothing it's very simple just
like the previous one I'm just using
this plot command and the file name
using oneness to do what we'll do is
that will take the first row first
column to be the x-axis the second
column to be the y axis and we'll just
plot now to make this run I have a small
shell script ok ok let me add this guy
right on the top remove dart remove star
dot exe star dot txt penny star dot Mart
there's no hard but anyway and then star
dot o most importantly ok
I know with the I'm compiling the main
file here the subroutine file here and
I'm you know collecting all the object
files under linking it to make an
executable rk4 ax dot exe here and then
running the exe file over here that's it
now I already have this Ikea for the
shell script made executable so if I
were to
all I have to do is I press f ing need
to run this make it to make it
absolutely simple okay let's see three
two one stop so if you guys notice this
this is the say this is the same
simulation we got last time and if you
guys notice this way we'll be moving
rather slowly especially constant out of
this region especially consider in this
region
you guys notice very carefully this wave
is actually moving forward and that is
kind of indicated by this wave moving
down now rk4 method is much more
consistent and precise when compared to
the compared to leader finite difference
key find a difference forward power time
backward scheme euler difference that
our fire power time backward scheme the
scheme that we used last time so what
will happen is that the values will
remain intact that remains more or less
very tight and they will be just
proceeding forward very nicely when
compared to the previous time and they
would be working just just fine then if
you guys can see this wave is again
propagating forward now in a normal in a
nice nice pattern in a nice manner and
because of this the linking thing we did
by connecting the last point of the
previous time to the first point of the
next current time step it's kind of
moving very nicely to the one and from
one and from the other now for this kind
of values again this function will go
for a long long time okay
what I will do is that I just what I'll
do is that let's close it close this
fight anyway I terminate the video but
it because it will just go on run for
another some more time so he'll I don't
want to bore you guys with it now this
guy now with this way you know if you
guys know this you guys will be having
sufficient practice with working with
executable I mean the plot files and all
gnew plot I recommend for people who
wanted to get some quick plots and
wanted to work with the seriously with
Fortran okay new plan
and even see Keeney plot is highly
recommended for me from me there
very good and it's simple and for
starters it's extremely free and it
follows very simple programming language
that can be used so have a look at it if
you guys want well well then this being
done
that's a this being done this will take
some more time for this function to
close ok this plot to close it will take
a long time
this plot is because they're like huh
roughly four 4,000 plots to be plotted
so do you take time so so I do be just
going on so I'll just
well it'll be just going on and then I
don't want to bore you guys I'll just
close the video because I had nothing
nice to tell in this video thank you
guys for watching see you guys next to
their will hope you like this one see ya
maybe in the next tutorial we look at
some real look at some other interesting
topic ok bye