Setting the initial water height values

Before iterating the solution forward in time, we need to specify the initial state of the water height. In other words, we need to set the initial values of the array h. A common choice in toy models like ours is a Gaussian (bell) shape (figure 2.6).

Figure 2.6 Initial values of water height

This is a well-defined, smooth function that we can calculate using the exponential, as shown in the following listing.

Listing 2.11 Initializing the water height with a Gaussian shape

integer, parameter :: icenter = 25         ❶
real, parameter :: decay = 0.02            ❶
 
do i = 1, grid_size                        ❷
  h(i) = exp(-decay * (i - icenter)**2)    ❸
end do

❶ Central index and decay factor of the shape

❷ Loops over all elements, from index 1 to grid_size

❸ Calculates the value and assigns it to each element of the array

Here, we have the first practical use of the following:

• A do-loop to iterate over array elements. Since we’ve declared h as an array of size grid_size, this loop will iterate over all elements of the array.
• Arithmetic operators -, *, and ** (power).
• An intrinsic mathematical function exp (exponential). This and other intrinsic math functions are readily available to use in Fortran programs and don’t need to be imported in any special way.
• Arithmetic assignment (=) of the result of the expression on the right side to the variable on the left side. The value of the left side is updated only after the whole expression on the right side has been evaluated.

Parameters icenter and decay control the position and width of the water height perturbation, respectively. The integer icenter is the array index at which the perturbation is centered. For example, when i equals icenter, the exponent argument reduces to zero, and h(i) equals 1. The real parameter decay determines the rate of exponential decay. Larger values yield a thinner perturbation.

Can our array assignment be done in parallel?

Recall our discussion of embarrassingly parallel problems in the previous chapter. We said that a problem is embarrassingly parallel if there’s no data dependency between individual iterations. Take a look at our expression for the initial value of h. Could we distribute this workload among multiple processors?

Begin the practice of asking that question for every computational problem, formula, or equation that you encounter. Over time, you’ll find more opportunity to distribute the computation, or at least mark sections of the code that are safe for the compiler to vectorize. Fortran offers a special do loop for this purpose, called do concurrent. It guarantees to the compiler that there’s no dependency between individual iterations and that they can be executed out of order, as we’ll see in the next subsection.