My Project Code is Working
After one of the most stressful work-related weekends I have had yet, I finally have my project code bug free (or so I hope), and I am in the business of producing pretty pictures. In fact, one such picture is below (click on the picture to see a larger version).
I spent about four solid days debugging my code, at first trying to get the code to compile, then trying to get the code past the seemingly endless number of arithmetic exceptions it would generate, and finally trying to arrive at a plausible answer. It turns out my solution agrees very well with three other people in the class, so I am reasonably confident I have the bugs worked out. The total project code is approximately 1,000 lines long, and as a result of all of the debugging attempts, I could practically write it again from memory.
What was the final culprit, you ask? Well, the bug that took tens of hours to find was an extra factor of "dx" (the spacing between vertices on the computational grid) in one of my linear operators. How could I miss such an easy bug? Well, I got a little too crafty for my own good. I had debugged the module that performed all of my matrix multiplications by comparing my output with the output from a Mathematica script written by my advisor. I had chosen a nice test vector and nice values of "dx" and "dy" so as not to have two errors cancel themselves. Everything worked beautifully.
In one of my more desperate moments, I had thought maybe the garbage I was seeing (which turned out to be the result of another bug I subsequently fixed) was roundoff error accumulating somewhere I did not expect it. In the module that performed all the matrix multiplications, I had left off the "dx" and "dy" factors until the end of the subroutine, and then I would scale my output by the appropriate factor (e.g. I would calculate some output vector "output", and then at the end of the subroutine do "output=(1d0/dx)*output" or whatever was appropriate). I decided to go ahead and comment that line out and add the "(1d0/dx)" term to each line of the subroutine, hoping that would somehow solve my problem. Once I had compared a subroutine against the Mathematica output and saw that it worked, I would copy that subroutine back into the module for my main program.
To make a long story not quite so long, I forgot to comment out the last line in one of the subroutines I had already copied to the module for my main program, so that I had an erroneous "output=(1d0/dx)*output", which was wrecking everything. I had assumed the matrix multiplication module was completely bug free (because the output had agreed with the Mathematica script and I had commented out the lines correctly in the tested version), so I did not put a lot of energy into combing through the copied version in my actual project code during the debug process. Lesson learned.
13 Comments:
Good job on the debugging, and nice picture. However, you've left out a few things:
1) Title for pretty picture
2) Labels on your axes
3) Explanation of what the pretty picture is representing (fluid flow? heat flow?)
4) What language you're programming in (C/C++ ?)
Oooh... pretty.
5) An absolute path to the picture so it loads right both on the main page and on the post's individual page. :)
I must say, your story of debugging matches well with my brief exposure to computational simulation work during my undergraduate time.
Based on that, I can safely say that I have the utmost respect for people such as you, who can survive the frustration of things just not working out. I'd break something. Expensive.
Yes, the pretty picture was mainly just for looks. I generated it pretty quickly after the code compiled, so I just cropped it and posted it to the web without editing it the way I will to include it in my project report.
1) I don't have a good title yet.
2) The graph is non-dimensional, so the labels are "non-dimensional y" and "non-dimensional x".
3) The pretty picture shows temperature contours (flooded) and streamfunction contours (lines) for air in a box that is heated from below and cooled on the top and sides.
4) Fortran.
5) What do you mean? I was worried I was somehow making it load differently based on platform. It seems to load fine for me (Firefox running on Linu). What are you seeing on your end? How do I add an absolute path?
Linu=Linux
Okay, I think I fixed the image path issue. Let me know if that works better.
yay! I love pretty pictures. Good job dude!
Yikes! FORTRAN! Double-Yikes! It's rather amazing that you can probably do the exact same thing in 10-15 lines of Matlab that you did with 1000 lines of Fortran. Are you using Fortran90/95 or *shudder* Fortran77? On the scale of programming language attractiveness, between 0 and butt-ugly, Fortran is only a little bit more attractive than Perl (butt-ugly) and Assembly Language/Machine Code (retina destroyingly ugly).
Let me try to put this as politely as I can. Don't pass judgement until you know exactly what's being computed. Trust me, you can not do the same thing in 10-15 lines of Matlab that I did in 1,000 lines of Fortran code. People in the class have tried to code this project in Matlab; it doesn't work well. Plus, one point of the project is to adhere to an order N matrix multiplication scheme, and it is reasonably difficult to force Matlab to strictly follow this requirement. And, granted, if I sacrificed readability and tried to optimize the code, I could probably cut it down to 600 or 700 lines or so, but it is much easier to look at and debug the way it is currently written. Fortran is also orders of magnitude faster than Matlab. I am using Fortran 95, not 77, only because of the allocatable arrays feature. Fortran is faster than C(++) for most of what we do computationally, so there are plenty of reasons to prefer it. I don't see why it's garnered the reputation it seems to have, but I have had no problem with it to date.
I don't see why we couldn't adhere to the O(NM) requirement, using Matlab, after all it does just what you ask it to do. I think its just prejudice in our friend T.
For a project like this, a matlab code wouldn't have be a single line shorter. Matlab will speed you up (as measured in #lines) when you can use it's intrinsic functions, such as the FFT. We'd still have to loop over every single node and its mother.
Unless there's a Matlab routine called SolveMyCFD_Project? Kidding.
The difference between doing exactly the same N*M ops in Matlab and Fortran, is that the latter eats Matlabs' lunch, timewise. Just try it, e.g. with a nested loop. The main difference stems from the fact that you're running (potentially optimized) compiled code when running your Fortran program.
I think most of Fortrans 'punch card' reputation stems from the many archaic features of F77, such as the non-free form source code. F90 is very nice, and this is coming from someone who was intent on using C for computation.
BTW I love Matlab.
CFD in MATLAB is challenging... It is not impossible but is too inefficient. I must make myself clear that I love MATLAB more than anything else but it is simply not suited for large scale computation.
General public is easily tricked that MATLAB is capable of doing anything. However, even some of the canned routines proved to me that they are useless.
This is obvious if you search for the functions that have been removed in previous versions of matlab, such as f---- (...you can check yourself).
If you are not familiar with numerical analysis, it is very easy to accept "trash" as "the" solution, i.e., try inverting a hilbert matrix in MATLAB. In order to avoid these pitfalls, it is still necessary to write our own codes in Fortran/C/C++. I personally see nothing wrong with Fortran being used. It is a great programming language. I have not heard of C/C++ being significantly better than Fortran for scientific computation.
Previously, I have written an "efficient" (may be questionable... but) MATLAB code that fully utilizes the sparse nature of the operators (to stick to O(MN) operations) and also using "mcc" function. But there is "NO" way MATLAB can achieve the speed of Fortran.
I personally believe that one needs to code one's own code (functions in the case of MATLAB) to control the accuracy and have accountability in the solution. If you want to report a simulation based on 15 lines in matlab, you better be darn prepared to justify your solution (unless 1% error or whatever level is satisfactory for you). I personally would not trust a 15-line-program from MATLAB. And I think anyone in the scientific community with a reasonable background would agree.
Although I passionately love MATLAB (when it is the appropriate time to use it), it is so sad that some of the general public thinks 15 lines of MATLAB can solve anything in the entire universe...
I'm glad you guys jumped in, because both of you know much more about Matlab than I do (and thus know more about how it might handle a project of this magnitude). I think the only reason it would be hard to get Matlab to strictly adhere to the O(N) matrix multiplication scheme is because of laziness. Would I have spent tens of hours debugging my matrix multiplications in Matlab? Probably not. I likely would have cheated and just had it do the matrix multiplications I wasn't doing correctly and move on. Other than that aspect, you're right, it will do exactly what Fortran will do, only eons and eons slower.
Plus, one of the intentions of doing this particular project was so we knew how to write a code that was scalable to 3D and much higher spatial and temporal resolutions. Matlab wouldn't stand a chance in that setting.
My comments about Fortran vs. Matlab were meant (partially) as a joke. Yes, Fortran will run faster, but the current version of Matlab has everything you could want, plus the kitchen sink, and some other stuff (including code readability). At Rice all of my applied math courses used Matlab, including the optimizations courses. Additionally, one of my profs, Dr. Danny Sorensen (of ARPACK fame and practically a god when it comes to large scale numerical linear algebra), uses Matlab for at least a portion of his work.
Anyway, I guess what I'm getting at is that Matlab will do everything that you can do in FORTRAN, but your code will look better and it might run a little bit slower (but not *that* much slower if you do it right, and it helps to run from the command line rather than using the memory-hungry GUI). I have nothing but great respect for people who are willing (and able) to use Fortran. My experiences with it, however, have not been great (both in terms of quantity and quality).
Oh yeah, as for my 1000 lines -> 10 lines comment, it was exaggerated. But I *have* had C codes that were a few hundred lines (maybe 500-700) condense down into less than 100 lines with insignificant slowdown.
Post a Comment
<< Home