The main goal of this chapter is to apply derived types toward refactoring our parallel tsunami simulator. Specifically, we’ll use them to abstract away the low-level code to allocate and manipulate data arrays, as well as to black-box much of the tedious code used for setting up the data structures. Derived types will help us write more expressive and concise code, which will be ever more important as we work toward making it more general and user-friendly. They’ll also allow us to expand the solver from one to two dimensions without compromising the simplicity of the main program.
At the end of the previous chapter, we left off with the core of the tsunami simulator as shown in the following listing.
Listing 8.1 The main time loop of the parallel tsunami simulator
time_loop: do n = 1, num_time_steps
h(ime)[left] = h(ils) ❶
h(ims)[right] = h(ile) ❶
sync all ❷
u = u - (u * diff(u) + g * diff(h)) / dx * dt ❸
sync all ❹
u(ime)[left] = u(ils) ❺
u(ims)[right] = u(ile) ❺
sync all ❻
h = h - diff(u * (hmean + h)) / dx * dt ❼
gather(is:ie)[1] = h(ils:ile) ❽
sync all ❽
if (this_image() == 1) write(unit=output_unit, fmt=*) n, gather ❽
end do time_loop❶ Updates halo cells for water height
❷ Waits for all images before proceeding
❸ Updates the solution for water velocity
❹ Waits for all images before proceeding
❺ Updates halo cells for velocity
❻ Waits for all images before proceeding
❼ Updates the solution for water height
❽ Gathers the water height on image 1 and prints it to the screen
In each iteration of the time loop, this code computed new values for water velocity, u, and height, h; synchronized the data with the neighboring parallel images; and wrote the water height data to the screen. This solver produced a fairly realistic-looking simulation of a propagating water wave (figure 8.1).
Although working directly with coarrays has served us well so far, and helped us learn how they work, this approach may not be ideal in the long run. That’s the case because each piece that we needed for the solver to work led to quite a bit of boilerplate code, specifically
- Calculating start and end indices on each parallel image
- Explicitly allocating coarrays
- Updating the halo points between neighboring images
- Synchronizing images to prevent race conditions
- Gathering and writing data into a file
This is bound to get worse as we extend the solver from one to two dimensions. We’ll go from two to three equations (two for each of the x– and y-velocity components and one for water height), and each of the equations will have additional terms in the y dimension. Since most of these operations are the same for each of the variables, we can express most of this code as derived type components and methods once, and reuse it whenever we need it.
Here are the main advantages of formulating our fields as derived type instances instead of bare multidimensional coarrays:
- We can use type-bound procedures for operations that are common to all fields.
- We can abstract away a lot of the low-level code, such as partitioning the computational domain, calculating start and end indices, setting initial conditions, and synchronizing the tiles in parallel mode.
- If we decide to implement more features in the future, we can write them directly in the definition of the type and its methods, keeping the main program concise and clean.
Figure 8.1 Simulation of water height initialized as a bell-shape 1 m high and about 20 m wide
At the end of this chapter, the tsunami simulator will produce a solution of a water wave radiating outward from the center of the domain, just as if we threw a pebble in a pond (figure 8.2).
Figure 8.2 Two-dimensional solution of a ripple wave radiating away from the center. The color is scaled to the minimum and maximum values of water height in each panel.
In the next section, we’ll go over the basic syntax and rules for defining, declaring, and initializing derived types. At each step, we’ll tie the new knowledge into the tsunami application and gradually build a derived type to model a physical field, such as water height and velocity.