MSA-Arrays-Ruby

__NOTITLE__

Table of Contents = Arrays (Ruby Version) =

The DRY Principle
Dan: What a difference a second-order scheme makes! Clearly, we can the same accuracy with far fewer calculations than with a first-order scheme. I wonder whether we can go to even higher order schemes, like third order or fourth order or who knows what order.

Erica: Yes, in stellar dynamics, depending on the application, various orders are used. In star cluster calculations, for example, traditionally a fourth-order scheme has been the most popular. In contrast, in planetary dynamics, people routine use even higher-order schemes, like tenth or twelfth order schemes.

But before we go any further in climbing the ladder of integration orders, I really want to write a leapfrog code. The modified Euler version that we have discovered was interesting as such, but in astrophysics at least, the leapfrog is used much more often. Presumably it has some advantages, and in any case, I'd like to see how it behaves.

Dan: What is a leapfrog?

Erica: It is probably the most popular integration algorithm, not only in astrophysics, but also in molecular dynamics, as well as in other fields. It combines simplicity with long-term stability. Shall I sketch it out?

Carol: Before we do your `before', I have an urgent wish: I want to clean up the last code we've written. If we just go on adding extra lines to produce higher-order codes, pretty soon the code becomes a bunch of spaghetti.

Look, everything we do is spelled out separately for the x coordinate, and the again for the y coordinate and then once again for the z coordinate. That's just plain silly. It violates the most basic principle in software writing, the DRY principle: Don't Repeat Yourself.

Dan: What's so wrong with repeating yourself?

Carol: Lot's of things are wrong with that! First of all, repetitions make a code unnecessarily long, and therefore harder to read. Secondly, if you want to modify a feature of a code, it is very difficult to do so correctly if that same feature is repeated elsewhere in the code,even if it is repeated in a place nearby. It is very easy to overlook the repetition, and only modify the first instance that you encounter.

Related to that is a third point: even the first time around that you write a code, if you start repeating yourself, it is quite likely that you make a mistake. ..

Erica. . . as we did in our very first code, with our typo!

Carol: Yes, indeed, I'd forgotten that already. Yes, that was a classic example of the type of penalty you can get for violating the DRY principle!

Vector Notation
Dan: When we were drawing pictures, we could look at the vectors themselves, but when we started coding, we had to go back to the components of the vectors. Are you suggesting to introduce some graphical way of coding, in which we can specify what happens directly to the vectors, as arrows, rather than to their separate x, y, and z components?

Carol: Well, in a way. Until the middle of the previous century, mathematicians often wrote vector equations in component form. But then they started more and more to use vector notation, by writing down symbols that stood for vectors as such, without any mention of components. On the level of the mathematical equations we have written down, we have used vector notation right from the beginning: we happily wrote things like r$2$=r$1$+v$1$dt on paper, but then we tediously wrote the same thing on our computer screen as:

x2 = x1 + vx1*dt y2 = y1 + vy1*dt z2 = z1 + vz1*dt

Erica: So you would like to write a computer code with lines like

r2 = r1 + v1*dt

where it would be understood that r2, etc., would be an object with the three components { x2, y2, z2 }.

Carol: Yes, exactly! But for that to work, a lot more should be understood. For example, it should also be understood that the simple + symbol is now a much more complicated addition operator. It should be clear to the computer that each of the components of r1 should be added to the equivalent component of the second expression, v1*dt. And in that last expression the * symbol should in turn be understood to be a more complicated multiplication operator. Multiplying a vector v1 with the scalar quantity dt should be understood as multiplying each of the components of the vector with the same scalar.

Dan: I like the idea of simplifying the code, and making it look more like the pen-and-paper expressions, but boy, the computer will have to understand a lot! Let me write down what you just said, to see whether I at least understand it.

Writing in the code a = b + c</tt> for three vector quantities a, b, c</tt> should be translated automatically into the following code fragment

ax = bx + cx ay = by + cy  az = bz + cz

where  ax</tt> is the first component of the vector  a</tt>,  ay</tt> is its second component, and so on. And writing in the code a = 3*b</tt> will be translated into

ax = 3 * bx ay = 3 * by  az = 3 * bz

Carol: Yes, exactly. That would be nice, wouldn't it?

Erica: Well, let's try to make that work. The first thing that comes to mind is to use arrays. If we represent a vector by an array, then each element of an array can contain a component of the vector.

Dan: That makes sense. I hope Ruby has arrays, just like Fortran?

Carol: Ruby sure does. But, as you can guess, they are far more powerful. A single array can contain objects of different types in different elements, and the length of an array can grow and shrink.

Dan: Seems like overkill to me. But who cares, let's get started.

Arrays
Erica: Before rewriting our modified Euler code, let us start with the simplest case, and rewrite our original forward Euler code, euler.rb</tt>.

Carol: Here, let me translate that code, line for line, into array notation. That way we can make sure that we perform the same calculations. Here is file euler_array_try.rb</tt>:

include Math r = [1, 0, 0] v = [0, 0.5, 0] dt = 0.01 print(r[0], " ", r[1], "  ", r[2], "  ") print(v[0], " ", v[1], "  ", v[2], "\n") 1000.times{ r2 = r[0]*r[0] + r[1]*r[1] + r[2]*r[2] r3 = r2 * sqrt(r2) a[0] = - r[0] / r3                               a[1] = - r[1] / r3   a[2] = - r[2] / r3   r[0] += v[0]*dt r[1] += v[1]*dt r[2] += v[2]*dt v[0] += a[0]*dt v[1] += a[1]*dt v[2] += a[2]*dt print(r[0], " ", r[1], "  ", r[2], "  ") print(v[0], " ", v[1], "  ", v[2], "\n") }

As you can see, I have simply replaced  x</tt> by r[0]</tt>,  y</tt> by r[1]</tt>,  z</tt> by r[2]</tt>, and similarly for the velocities and accelerations, I have replaced  vx</tt> by v[0]</tt> and <tt> ax</tt> by <tt>a[0]</tt>, and so on for the other elements.

Dan: From your example, I guess that arrays start with element zero?

Carol: Ah, yes, that's true, like in C, where the first element of an array <tt> foo</tt> is <tt>foo[0]</tt>, unlike Fortran, where the first element is <tt>foo[1]</tt>.

Erica: I noticed that right at the start of the program, you have done a bit more already than simply replacing <tt>x = 1</tt> by <tt>r[0] = 1</tt>, <tt>y = 0</tt> by <tt>r[1] = 0</tt>, and <tt>z = 0</tt> by <tt>r[2] = 0</tt>.

Instead, you have directly assigned all three values in one line by writing <tt>r = [1, 0, 0]</tt>.

Carol: That's true too. I guess I'm getting already familiar enough with Ruby that I had not noticed that I had made a shortcut. Yes, this is a nice example of the compact way in which Ruby can assign values. In fact, this line, which is the first line of the program after the <tt> include</tt> statement, defines <tt> r</tt> as having the type `array', and in fact an array of three elements. Similarly, the second line defines <tt> v</tt>, too, as an array containing three elements.

In the case of <tt> r</tt>, its three elements are integers, and in the case of <tt> v</tt>, the first and last elements are integers, and the middle element is a floating point number. As I had mentioned before, Ruby allows each array element to have its own dynamic type. In practice, though, as soon as we start calculating with these numbers, most of them will quickly become floating point numbers. Adding an integer and a floating point number, or multiplying an integer with a floating point number, is defined as giving a floating point number as a result, as you would expect.

Declaration
Dan: That all sounds reasonable. Shall we check whether we get the same result as we did before?

Carol: Good idea. I'll just print the last line, so that we can compare that one with our old result:

|gravity> ruby euler_array_try.rb | tail -1 euler_array_try.rb:13: undefined local variable or method `a' for main:Object (NameError) from euler_array_try.rb:10:in `times' from euler_array_try.rb:10 1 0  0  0  0.5  0

Dan: Well, something went wrong right at the start. It seems like there is a problem with the acceleration array <tt> a</tt>. The only output line we got was from the print statement before entering the loop.

Carol: Ah, I see. The two variables <tt> r</tt> and <tt> v</tt> are recognized as arrays, because they are defined as such, in the first two assignment lines of the program. The line:

r = [1, 0, 0]

says clearly: <tt> r</tt> is an array with three elements, and the elements are 1, 0, and 0.

In contrast, the first time that <tt> a</tt> is used occurs in the line:

<a class=codetext href=.euler_array_try.rb.html#segment-2> a[0] = - r[0] / r3                            </a>

and here we are not specifying what <tt> a</tt> is; we are not assigning anything to <tt> a</tt>. Instead, we are assigning a value to a element of <tt> a</tt>, as if <tt> a</tt> had already been defined as an object that has elements.

Erica: In C, you could just declare <tt> a</tt> to be an array. But you have told us before that in Ruby, because of dynamic typing, there was no need to declare the type of a variable. What gives?

Carol: Well, in this case we do need to give some extra information, otherwise the Ruby interpreter cannot possibly know what we mean. And yes, here we effectively need to declare <tt> a</tt> as an array.

Carol: How should we do that? We could give fake values, say <tt>a = [0, 0, 0]</tt> right at the beginning of the program.

Carol: That would be confusing, since a reader would wonder what the meaning would be of those three zeroes. In fact, in Ruby there is no need to specify how many elements an array has. All we need to say is that the variable <tt> a</tt> has the type of an array. Or more precisely, in Ruby's terminology: the class of <tt> a</tt> is <tt> Array</tt>, which is one of the built-in classes.

Classes
Dan: What is a class?

Carol: In Ruby, the word class is used to talk about the collection of possible things that have a particular type.

Dan: What exactly is a type?

Carol: The number 3 is an example of the type integer, the string "hello" is an example of the type string, and an array [3, "hello"] is an example of the type array.

In Ruby, roughly speaking, each concrete example of a type is called an object, and the collection of all possible objects of a given type is called a class. Each object of that type is called an instance of the class that corresponds to that type. In Ruby, the character string "hello" is an object, an instance of class String.

Dan: What about numbers?

Carol: In Ruby there is the class of <tt> Fixnum</tt>, which is the class of integers, numbers such as 7 or -5, and the class <tt> Float</tt>, which is the class of floating point numbers, such as 1.41421, and so on. Some classes, such as arrays, allow objects to contain objects of other classes. For example, the object <tt>[-5, 1.41421]</tt> is an array, which means it is an instance of the class <tt> Array</tt>, while the two elements are instances of the classes <tt> Fixnum</tt> and <tt> Float</tt>, respectively.

Erica: That's very different from how it is done in C++, where all the elements of an array always have the same type. In C++, you can have an array <tt>[-5, 8, 3, ...]</tt> and an array <tt>[1.41421, 3.14, ...]</tt> but you certainly cannot mix the types, by giving each element of the array a different type.

Carol: Well, in Ruby you can. This is one of the many ways in which Ruby is far more flexible than C++.

Coming back to Dan's question about classes, let's see how we can interrogate Ruby about class membership:

|gravity> irb a = [-5, 1.41421] [-5, 1.41421] a.class Array a[0].class Fixnum a[1].class Float a[1].class.class Class a.class.class Class quit

Erica: Ah, so each class, like <tt> Fixnum</tt> or <tt> Array</tt>, or whatever is a member of the class <tt> Class</tt>?

Carol: Yes, more precisely, every constant or variable in Ruby is an object, and therefore must be an instance of a class. The number 3 is an object of class <tt> Fixnum</tt>, and the class <tt> Fixnum</tt> is an object of class <tt> Class</tt>. Note the convention: in Ruby the name of a class always starts with a capital letter.

Erica: And what about <tt> Class</tt>? What is that an instance of?

Carol: Try it, ask <tt> irb</tt>.

Erica: Okay:

|gravity> irb 3 3 3.class Fixnum 3.class.class Class 3.class.class.class Class quit

So <tt> Class</tt> is a member of class <tt> Class</tt>. That sounds circular!

Carol: But it isn't. This is an example of the strength of Ruby. Like Lisp and other seemingly circular languages, Ruby has the power to invoke itself, without arbitrary boundaries between metalevels. This is one of my favorite topics. Shall I explain how it works?

Dan: Not today. I got enough of an idea of what a class could be. Let's keep writing code.