%\section{Sound Simulation and Visualization}
%\label{chapter:sound}

The goal of this chapter is to explore the techniques used when accelerating
an application of the GPU with CUDA. The goal of this application was to serve
as a proof of concept for a larger application. For this prototype, a
visualization of the of the sound waves in the room was desired so that phenomena
like standing waves could be quickly identified. In order to produce the
visualization quickly, some form of acceleration would be required. The GPU was
chosen as the accelerator using CUDA to produce the data for the simulation.

\section{Background}

As discussed in Chapter~\ref{chapter:background}, a factor in the rise of the use
in GPGPU has been the adoption of GPGPU as an accelerator for many different
simulations. Of particular interest to this application are simulations which model the movement
of waves through a medium. CUDA has been shown to greatly accelerate the
computation of solutions for Finite-Difference Time-Domain (FDTD) methods
\cite{livesey2012development,zunoubi2010cuda}. FDTD methods work
by breaking the computational domain, commonly a physical region being modeled,
into a grid. Then the solutions for the chosen equations are
calculated for each point in the grid for each time step. The benefits of the
method are: values are calculated for each point in the grid which makes
visualizations easier to create, a large range of frequencies can be simulated
with proper input allowing phenomena to present themselves and the data
parallelism inherent in the methods aids in accelerating the computation of
results.

\subsection{Computational Acoustics}

Computational Acoustics is largely based on the physical modeling of the
propagation of sound waves. This propagation is governed by the linearized Euler
equations\cite{tam1993dispersion}. There have been many methods developed to
model the propagation of sound waves; these methods fall into one of two categories:
geometric and numeric. 

The geometric methods tend to be based on a technique called ray tracing, more
commonly used to create computer graphics. These ray tracing methods are able to
quickly model the acoustic properties of a room by assuming that high frequency
sound waves behave similarly to a ray of light. However, these methods have a critical
flaw. Because the wavelengths of audible sounds are on the order of
commonly found objects, sound waves exhibit diffraction (bending) where light
would not. Because of this recently these geometric methods have fallen out of
favor as more accurate models have become accessible.

The numeric methods attempt to directly solve the Wave Equation. The benefit of
this approach is that phenomena such as diffraction are accounted for. However,
numerically approximating the wave equation can become expensive, especially for FDTD
methods. This is due to the requirements for memory and computation scaling up
as the sampling rate desired increases. With the the time step T defined as
$T = 1/f_{s}$ and the grid spacing g defined as $g = \frac{344 * \sqrt{3}}{f_{s}}$
it can be seen that as the sampling rate increases not only does the time step of
the simulation get smaller, but the size of the grid also increase. These factors
combined have meant that until recently FDTD solutions for the acoustic wave
equations have been prohibitively expensive.

Recently solutions using CUDA have shown promising results in
accelerating 3D FDTD solutions to the wave equation\cite{sheaffer2010fdtd,webb2011computing}.


\subsection{CUDA Streams}

As discussed in Chapter~\ref{chapter:background}, CUDA expose the computational
power of the GPU through a C programming model. It additionally provides an API
for scheduling multiple streams of execution on the GPU. This allows the hardware
scheduler present on CUDA enabled GPUs to more efficiently use all of the compute
resources available. When using streams, the scheduler is able to concurrently
execute independent tasks. In order to fully use streams for all  including
memory transfers, the CUDA driver must allocate all memory on the host that will
be used for CUDA API calls. This ensures that the memory is pinned (page-locked)
so that Direct Memory Access (DMA) can be enabled.

However, using streams is the only way to get the full performance of the GPU using CUDA,
and it comes with some concerns. The first concern is that pinned memory can
not be paged out and can therefore impact the amount of paging for any other
memory that is virtually allocated. This means that if too much pinned memory is
allocated, other components of an application may see a performance loss. In
addition care must be taken to order CUDA memory transfers and kernel launches
in such a way that the scheduler is able to properly schedule each of the actions
concurrently\cite{nvidiawebinar}. Additionally some advanced memory features like
shared memory or texture memory can become restricted when using streams.

\subsection{OpenGL}

OpenGL is an API for rendering 2D and 3D computer graphics. It was originally
developed by Silicon Graphics Inc. It provides a portable API for creating
graphics which is independent of the underlying hardware. It is widely used for
a variety of applications ranging from Computer Aided Design (CAD) to scientific
visualizations.

At its core, the OpenGL API controls a state machine. This state machine maintains
the hardware and is in charge of ensuring the proper resources are loaded at the
proper time. It is important to note that for performance reasons this state
machine was not designed to be used with multiple threads. The OpenGL driver can
be made to run in a multithreaded way, but the driver does nothing to protect the
programmer from race conditions. The API has calls for uploading geometry data,
texture data and more to the GPU. In addition, it also exposes a compiler for the
GLSL shading language.

As discussed in Chapter~\ref{chapter:background} the GLSL shading language
(shaders) gives the application programmer more control over the functions of
OpenGL. The programmer can use these shaders to accomplish advanced rendering
techniques such as lighting effects, producing rolling ocean waves or
programmatically generating textures or geometry. 

One technique of interest to the application discussed in this chapter is the
rendering of volumetric data. There are many techniques for accomplishing this
task\cite{ikits2004volume}. A popular method for volume rendering uses the
texture mapping hardware of the GPU and alpha blending to render volume metric
data. The technique involves rendering many semi-transparent slices of the 
volumetric data as shown in Figure~\ref{fig:sound:volume}. 

\begin{figure}
    \centering
    \includegraphics[height=\textheight,width=\textwidth,keepaspectratio]{figures/sound/vol_texMap3Drot.png}
    \caption{Texture Mapped Volume Rendering: a technique which uses alpha blending to quickly render volumetric data.\cite{evolveVolume}}
    \label{fig:sound:volume}
\end{figure}


\section{Overview}

The application presented in this chapter uses a 3D FDTD method for simulating
sound in a room. A couple of assumptions are made. A simplified model of the sound
absorbence is used in lieu of a more computationally expensive one and that all
sources of sound are omnidirectional.

The application is structured as a chain of producers and consumers. Once the
simulation has begun, the simulation manager begins to produce simulation data.
This raw form of the data is both unsuitable for visualization and is also located
in a section of GPU memory which the OpenGL driver is unable to access directly.
The simulation manager therefore publishes this raw data to a queue for processing.
The memory manager takes this raw data, copies it to the CPU and puts it into a
form suitable for visualization using OpenGL. These processed frames are
published to a queue until they can be uploaded to the OpenGL driver. The
renderer is responsible for all the calls to the OpenGL API. After uploading the
current frame as a texture, the frame is drawn using texture mapped volume
rendering. An overview of this structure can be seen in
Figure~\ref{fig:sound:arch}

\begin{figure}
    \centering
    \includegraphics[height=\textheight,width=\textwidth,keepaspectratio]{figures/sound/SoundArchitecture.pdf}
    \caption{An overview of the structure of the application}
    \label{fig:sound:arch}
\end{figure}

\section{Implementation}

The application is structure into three major components: the renderer, the
simulation manager and the memory manager. Each of these components is run on its
own thread to allow as many concurrent actions to occur as possible. An overview
of each component's execution is shown in Figure~\ref{fig:sound:threads}. Care was
taken to design each component in such a way that a component could be redesigned
and replaced easily. For instance if the data set would not fit on a single GPU
then the simulation manager could be rewritten to accommodate this and the rest of
the application could remain unchanged.

\begin{figure}
    \centering
    \includegraphics[height=\textheight,width=\textwidth,keepaspectratio]{figures/sound/soundThreads.pdf}
    \caption{An overview of the work done by each of the three threads in the application.}
    \label{fig:sound:threads}
\end{figure}

\subsection{Renderer/Main Thread}

The first action taken by the main thread is to initialize the state of the OpenGL
driver, create a window and initialize a OpenGL context within that window. These
actions are all required to begin making any other calls using the OpenGL API.
Assuming that there are no errors, the next thing done is the loading of the model
information and test sound sample. For this prototype, a test room was made. For
the simulation it is important that the model contains both geometry and absorbence
information for each point in the grid. The geometry data consists of a 3D array
of values encoding whether or that point in the grid is connected its neighbors.
An example of a possible encoding for the 2D case is shown in 
Figure~\ref{fig:sound:encoding}. Once both the room and sound sample are loaded
from disk, the data is used to construct the initial state of the simulation.
After the simulation is initialized, both the simulation manager and memory
manager are started and the rendering tasks begin.

\begin{figure}
    \centering
    \includegraphics[height=\textheight,width=\textwidth,keepaspectratio]{figures/sound/SoundEncoding.pdf}
    \caption{An example encoding using 4 bits, the 2 most significant bits represent whether a neighbor exists on the y axis and the 2 least significant bits represent whether a neighbor exists on the x axis.}
    \label{fig:sound:encoding}
\end{figure}

The first of the rendering tasks is to initialize the state of OpenGL, this
includes loading the shaders required to perform texture mapped volume rendering
as well as preparing OpenGL geometry data required for the technique. Once this
initialization has been done, the render thread falls into the render loop.

The render loop primarily consists of two actions: uploading frames to OpenGL memory and
performing the volume rendering. Uploading frames to OpenGL memory is the first
task done. The renderer checks to see if any processed frames have arrived in its
queue, shown in Figure~\ref{fig:sound:data}. If a processed frame has been
published then it is uploaded to OpenGL memory and the now blank frame is
published back to the memory manager. The renderer then volume renders the
current frame that is in OpenGL. This loop continues until a user inputs the exit
command or the simulation reaches its completion.

\subsection{Simulation Manager}

The simulation manager's first action is to allocate all the memory required for
the simulation (room information, blank simulation state arrays and input/output
arrays). It then transfers the encoded grid representing the room, the absorbence
values for each point in the grid and the input for the simulation. Once all of the CUDA
memory has been allocated and all of the initial simulation state has been
transferred onto the GPU, the simulation manager waits to be told to start.

Once the signal to start is received, the simulation acquires three simulation
blanks, representing the t$_{n-1}$, t$_{n}$ and t$_{n+1}$ simulation states, from 
the simulation blanks queue. With that last bit of preparation done, the simulation
manager drops into the main simulation loop.

The first step of the simulation is to update all of the input into the simulation.
Each input source has a location associated with it that is updated. After that
the simulation kernel is run. The kernel uses an equation in the form of
Equation~\ref{eq:cases}, where p$_n$ is the acoustic pressure for t$_n$, to numerically
approximate the propagation of the wave.
\begin{multline} \label{eq:cases}
p_{n+1}(x,y,z) = (1/3 [ p_{n}(x+1,y,z) + p_{n}(x-1,y,z) + p_{n}(x,y+1,z) + \\
p_{n}(x,y-1,z) + p_{n}(x,y,z+1) + p_{n}(x,y,z-1) ]) - p_{n-1}(x,y,z)
\end{multline} 
The form of the equation depends on the
value of the encoded geomotry for that point on the grid. Care is
taken to ensure that memory accesses are sequential for the warps assigned to the
the kernel. The kernel is run in its own stream, which allows the hardware scheduler
to schedule both the kernel and any copies that the memory manager is scheduling.
If the kernel was not run in its own stream, then the hardware scheduler would not
have enough information to schedule the two tasks concurrently. Once the kernel
has finished, any listeners present in the model room are updated. The simulation
manager then decides whether it's time to publish a raw frame to the unprocessed
frame queue as shown in Figure~\ref{fig:sound:data}.


\subsection{Memory Manager}

The memory manager is the simplest of the three components of the application.
As seen in Figure~\ref{fig:sound:arch} the memory manager acts a a bridge
between the simulation and the visualization. The memory manage takes the raw
frames published by the simulation and normalizes them before publishing them
to the renderer. The only major concern for the memory manager is that memory
transfers off of the GPU must be run in a separate stream of execution than the
simulation kernel. If this is not done, then any memory transfers will block the
execution of the simulation kernel.

\subsection{Communication}

\begin{figure}
    \centering
    \includegraphics[height=\textwidth,width=\textheight,keepaspectratio,angle=90]{figures/sound/SoundDataFlow.pdf}
    \caption{The communication between the three components (shown in yellow) using the various communication queues (shown in blue).}
    \label{fig:sound:data}
\end{figure}

The application uses three threads to run the three components concurrently.
These threads communicate using thread safe queues. The reason that queues were
chosen as the data structure for the message passing between the thread is that
they both naturally preserve ordering and allow the simulation and visualization
to not outpace each other. The application is essentially a producer-consumer
problem where both the producer and consumer are fighting for the same
resource. If there was not some way to limit which component controls most of the
GPU time, then the contention for the GPU would cause problems in either the 
simulation or the visualization. Additionally, the use of queue as the inter-thread
communication medium makes any future attempts to use multiple machines or devices
easier by allowing the queue to hide the origin of any information.

\section{Results}

The prototype was tested on a computer with an Intel Core2 Quad core CPU Q9450 and 
a NVidia GeForce GTX 480 GPU which contains 15 groupings of 32 compute cores for
a total of 480 compute cores and has 1.5 GB of on board memory. Two different
test signals were used, one a music sample and the other a 8kHz sin wave. The
test room modeled was 12m x 12m x 3m with a pillar in the center of the room.
When using the sin wave the simulation modeled standing waves where standing waves
would be expected to form. Figure~\ref{fig:sound:initial} and
Figure~\ref{fig:sound:running} show the visualization during a test run using the
music sample.

\begin{figure}
    \centering
    \includegraphics[height=\textheight,width=\textwidth,keepaspectratio]{figures/sound/InitialSound.png}
    \caption{The initial state before simulation initialization.}
    \label{fig:sound:initial}
\end{figure}

\begin{figure}
    \centering
    \includegraphics[height=\textheight,width=\textwidth,keepaspectratio]{figures/sound/SoundRunning.png}
    \caption{A rendering from the simulation while running.}
    \label{fig:sound:running}
\end{figure}


\subsection{Future Work}

This chapter presented an application that benefited greatly from the use of
GPU acceleration. The problem presented exhibited data parallelism which assisted
in the implementation of the kernel. This acceleration allowed the simulation to
run at a quick enough pace to facilitate the creation of a visualization alongside
the simulation. 

Despite this, there are issues that should be addressed in future work. Currently,
the simulation is limited in size and frequency range due to memory concerns.
This could be remedied by replacing the current single GPU simulator with a
simulator that uses multiple GPUs if present on the computer or even a clustered
simulator. Additionally, the visualization requires that the user has some
understanding of the mechanics of wave propagation to be useful. If a user wanted
to use the simulator and visualizer to place speakers in a room, it might be more
helpful for the visualizer to analyze the frames coming from the simulation to
programmatically find good speaker placements.
