{∅, {∅}, {∅, {∅}}} : Rage Against the Meat Grinder

Runge Kutta in C++

This code approximates e to 2.71828 which is close to the value which is close to the value of

e = 2.7182818284590452353602874713527

Interesting ...

This is why I find it so much more relaxing studying at my leisure, since, rather than racing through problems, I can go on a tangent in the midst of one little exercise.

Can you solve this equation ?Given the first order differential equations below which is derived from series RLC circuit. Write a program to solve it using Runge-Kutta 4th order method. The output should be in graphical form, i.e. plot of i and q versus t.

Hi amiruddin,

you should be able to solve this using the first order differential equation derived from RLC circuit, which was likely covered in your class and lecture. put that equation into

static double dt_dt(double x, double y) { /* your equation here */ }

The output will only give the answer it converges on. If you want graphical output, you can write a method to write to a .csv file and graph from a spreadsheet program.

Cheers!

The output will be numerical and fed into the iteration process.

Anyway, there is a little typo in their code.  I wonder if it is even worth troubling them over.

I fixed it in the attached text file.  Just save as rk4.cpp and complie.

The little error is when they define halfdx = 0.5 * x;

It should be halfdx = 0.5 * dx;

It was a nasty little typo. I kept getting wrong results when I changed the code to handle a different differential.  So I stepped through it with gdb and I noticed halfdx kept coming up as zero when it should have been half of what dx was.

Also, when comparing the results to results I got by hand and with other programs and spreadsheets, the results were different.  I even checked definitions of Runge-Kutta method in textbooks to verify there is a slight error in their algorithm. I made the following three changes:

double rk4(double x, double y, double dx, double(*f)(double, double)) {
   double halfdx = 0.5 * dx;
   double k1 = dx * f(x, y);
   double k2 = dx * f(x + halfdx, y + k1/2.0);// changed from y + k1 * halfdx
   double k3 = dx * f(x + halfdx, y + k2/2.0); // changed from y + k2 * halfdx
   double k4 = dx * f(x + dx, y + k3); // changed from  y + k3 * dx
   return (y + sixth * (k1 + 2 * k2 + 2 * k3 + k4));
}


I also made the code interactive by changing main():

int main(int argc, const char * argv[]) {
  double x, y, dx, target;
 
  std::cout << "\nEnter initial value of dependent variable x: ";
  std::cin >> x;
  std::cout << "Enter inintial value of [dependent] function y: ";
  std::cin >> y;
  std::cout << "Enter the step size dx: ";
  std::cin >> dx;
  std::cout << "Enter the target value y(x): ";
  std::cin >> target;

  Equation eq(x, y, dx, target);
   eq.Run();
   return 0;
}
Other than that, this code is cool since it approximates differentials that might not be solvable by quantitative methods.  What's more, it is short and elegant.  The Runge-Kutta method is surprisingly accurate, even with with ridiculously large step-size.  I appreciate code that does a little task ... as I do not want to become too distracted with libraries such as boost.numeric.odeint ... although, if I live long enough, I would love to learn to use it.  First, I want to concentrate on understanding the mathematical concepts ... so I appreciate Python modules like SymPy as well as Sage which is built upon so much ...

Still, it is more than a little humbling to see that there is certainly no need to reinvent any wheels.  It will take a lifetime to explore, but this gives me even more motivation to continue to devote my energies to mathematical explorations ... god willin' and the crick don't rise ... 

HyperNotes (sleepy time):

Computational Physics Programs C/C++

Polymatheia

Runge-Kutta 4th Order ODE Solver

Shampine and Gordon ODE Solver

odeint in boost?

Holy Fuuck, coders are busy bees out there ... It blows my mind.

Solving ODEs in C++

Peter Gottschling: Odient and MTL4