
\documentstyle[12pt,lslide]{article}
\title{Case Studies on Parallelization of NASA's Semi-Lagrangian \\ 
        Transport Code : Results on Interpolation Techniques}
\organization{Manu Natarajan \\
              Northeast Parallel Architectures Center \\
              Syracuse University \\
              Syracuse, New York}
\date{}

\begin{document}
\titlepage

\section{ Introduction }

\begin{slide}
\subsection{ Terminology }

\begin{itemize}
\item circulate - message passing equivalent of ``cshift''.  Shifting of \
      thin, long slices of datum.
\item component   - resolution of a vector along the three axes, \
                     x, y, \( \phi \)
\item cshift  - synonymous to circulation
\item global index - in index range of undivided datum
\item hop  - number of nodes to be communicated with  
\item inter-nodal distance - synonymous with hop
\item lasagna - strip or thin long slice
\item local indices   - index varies from 1 to width of segment
\item nodes   - synonymous with processors
\item northern buffer - 3-element wide strip of datum added in sets of 
      points approach to an existing partition to allow a local interpolation
\item owner node  - node on which required data resides
\item partitions - disjoint sets of datum distributed among nodes 
      which form a partition of the original undivided datum 
\item segment    - contiguous parts of datum distributed on nodes which may 
      overlap with one another, i.e. not necessarily disjoint. In this 
      context it additionally refers to that slice of datum that permits a 
      local interpolation on a node.
\item sets of points - in this context refers to all points which are at a 
      specified displacement (vector) relative to a set of reference points.
\item strips - slices of data
\item wrap around     - synonymous with pole-crossover
\end{itemize}
\end{slide}
\begin{slide}
\subsection{ Bilinear Interpolation }

\begin{itemize}
 \item Standard Algorithm -  \( (1.0 - delx) * ((1.0 - dely) * w1 + \
       dely * w2) + delx * ((1.0 - dely) * w3 + dely * w4) \), where

  \begin{itemize}
    \item{w1, w2, w3, w4 are the four reference points on the grid}
    \item{delx, dely, is the distance of point at which interpolation
    is to be carried from the ``origin'' of the {\em grid-square} }
  \end{itemize}
 \item Economic  -  only 6 multiplications
 \item Therefore used where more accuracy is not required
\end{itemize}

\end{slide}

\begin{slide}
\subsection{ Bicubic Interpolation }

\begin{itemize}
\item Two-step One-dimensional Lagrange Polynomial (cubic) Interpolation.  
\item Formula : \(P(y) = \sum_{i=0}^{3}f_{i}L_{i}(y) \), where
\begin {itemize}
   \item {y is the point at which interpolation takes place.}
   \item {P(y) represents the desired value.}
   \item {\(f_{i}\) are the known values on the grid at points \(y_{i}\).}
   \item {\(L_{i}(y)\) are the  Lagrange polynomials such that 
   \(L_{i}(y) =  \Pi_{k \neq i, k = 0}^{3} (y - y_{k})/(y_{i} - y_{k}) \)}
\end{itemize}
\item {More involved computation - 32 multiplications}
\item {Used where accuracy is required}
\end{itemize}

\end{slide}

\begin{slide}
\subsection{ History - Part 1 : Data Parallel Version}

\begin{itemize}
\item Measured performance of data parallel version using {\em forall }
      construct
  \begin{itemize}
   \item Obtained poorer performance of the lasagna decomposition along E-W
     direction compared to the square patch decomposition.
   \item Initially thought to be due to lack of knowledge on part of 
     compiler about nature of data access.
   \item Also due to poor utilization of nodes.  That is run-time system did 
         not spread datum around equally to all processors due to small 
         grid size.
  \end{itemize}
\item Wanted to use the global-local paradigm of programming to perform a 
      lasagna decomposition.  Initial testing indicated that the data 
      distribution was not even (explained the reason for the poorer 
      perfomance of the lasagna decomposition).
\end {itemize}

\end{slide}

\begin{slide}
\subsection{ History - Part 2 : Message Passing Version}

\begin{itemize}

  \item  Chose the native message passing library (cmmd)
         because did not want additional overhead
         in using MPI.  Also chose lasagna decomposition due to simplicity.

  \item In message passing version, initially attempted the use of Fortran 90
        for the program that ran on the nodes.  It was not successful due to
        linker problems.
  \item Also attempted the usage of the hostless style of programming, there
        was some problem allocating memory for the program.
  \item Finally settled for the host-node style of programming using only 
        Fortran 77.
\end{itemize}

\end{slide}

\begin{slide}
\subsection{ Naive Algorithm}

\begin{itemize}
\item {Do for all (ii,jj)}
\begin{itemize}	
\item	Find the 4 points (ii,jj+2), (ii,jj+3), (ii+1,jj+2), (ii+1,jj+3)
\item   Obtain the values of the wind component at above points.
\begin{itemize}
	\item If the value is not available locally, {\em fetch} from
         the appropriate neighboring node
\end{itemize}
\item  Find the interpolation point.
\item  Compute the value at interpolation point.
\end{itemize}
\item End of do loop
\end{itemize}
\end{slide}

\begin{slide}

\subsection{ Pre-Determine Degree Of Circulation }

The algorithm for {\em the pre-determination of degree of circulation} is below
\begin{itemize}
  \item { ny = Number of nodes (processors) used }
  \item { Let 'myaddr' denote the address of the node executing the code}
  \item {Let 'ownp0' denote nodal address where value 'jj(i,j)+0' resides }
  \item {Let 'ownp3' denote nodal address where value 'jj(i,j)+3' resides }
  \item { Number of elements in E-W direction for a partition is imp3 }
  \item { Number of elements in N-S direction for a partition is range }
\end{itemize}
\end{slide}

\begin{slide}
\begin{itemize}
  \item {for j = 1, range}
  \begin {itemize}
    \item {for i = 1, imp3}
      \begin {itemize}
        \item {ownp0 = jj(i,j) / partition width }
        \item {ownp3 = (jj(i,j) + 3) / partition width }
        \item { Find the inter-nodal distance in {\em northern direction}
                from 'myaddr' to both ownp0 and ownp3}
	\begin{itemize}
          \item {distp0 = ownp0 - myaddr}
          \item {distp3 = ownp3 - myaddr}
          \item {if (distp0 .lt. 0),distp0 = ny + distp0}
          \item {if (distp3 .lt. 0),distp3 = ny + distp3}
        \end{itemize}

	\item {Find the shortest inter-nodal distance for both ownp0, ownp3, 
               and update the local maximum inter-nodal distances if required}
      \end{itemize}
    \end{itemize}

    \item {Find the global max. northern dist. (reduction operation)}
    \item {Find the global max. southern dist. (reduction operation)}
\end {itemize}

\end{slide}

\begin{slide}
\subsection{ Reconstruction of Sets of Points of Approach }
\begin{itemize}
\item {\em First level refinement} for the algorithm.

 \begin{itemize}
   \item { Number of elements in E-W direction for a partition is imp3 }
   \item { Number of elements in N-S direction for a partition is range }
   \item { Pre-determine the degree of circulation }
   \item { Construct the northern buffer. }
   \item {Select elements from segment}
   \item {Fetch segments and select elements from southern direction}
   \item {Fetch segments and select elements from northern direction}
   \item {Perform local interpolation}
  \end{itemize}
 \end{itemize}
\end{slide}

\begin{slide}
\subsection{ Reconstruction of Segment Approach}
\begin{itemize}
\item { First Level Refinement }
\begin{itemize}
   \item { Number of elements in E-W direction for a partition is imp3 }
   \item { Number of elements in N-S direction for a partition is range }
   \item { Pre-determine the degree of circulation }
   \item { Copy the partition into a large buffer called segment}
   \item { Fetch partitions from the southern direction and append to
	   segment. }
   \item { Fetch partitions from the northern direction and append to
	   segment. }
   \item { Convert global indices to local indices. }
   \item { Perform the interpolation locally. }
\end{itemize}
\end{itemize}
\end{slide}

\begin{slide}
\begin{itemize}
\item { Gathering From Southern Nodes}
\begin {itemize}
       \item {Current partition bounds = Original partition bounds }
       \item {Copy partition into temporary buffer (depart)}
       \item {For number of hops in southern direction }
       \begin {itemize} 
	 \item {Send upper \& lower bounds of current partition to northern
                node.}     
         \item {Receive upper \& lower bounds of new partition from southern 
                node.}
         \item {Send 'depart' to northern node.}
         \item {Receive new partition from southern node in buffer 'landing'.}
         \item {Update parameters needed for generation of local indices.}
	 \item {Move up the previously filled columns in the segment by width
                of the new partition (southern partitions occupy lower 
                indices in the segments [Fig. 8]). }
         \item {Append the new partition to the 'southern end' of the segment.}
         \item {Update width of segment.}
         \item {Copy 'landing' into 'depart'}
       \end {itemize}
      \end {itemize}
\end{itemize}
\end{slide}

\begin{slide}
\begin{itemize}
\item { Gathering From Northern Nodes}

\begin {itemize}
       \item {Current partition bounds = Original partition bounds }
       \item {Copy partition into temporary buffer (depart)}
       \item {For number of hops in northern direction }
       \begin {itemize} 
	 \item {Send upper \& lower bounds of current partition to southern
                node.}     
         \item {Receive upper \& lower bounds of new partition from northern 
                node.}
         \item {Send 'depart' to southern node.}
         \item {Receive new partition from northern node in buffer 'landing'.}
         \item {Update parameters needed for generation of local indices.}
         \item {Append 'landing' to the 'northern end' of the segment; 
                northern partitions occupy higher indices in the segment
                 [Fig 8].}
         \item {Update width of segment.}
         \item {Copy 'landing' into 'depart'}
       \end {itemize}
      \end {itemize}
\end{itemize}
\end{slide}

\begin{slide}

\subsection{ Conversion to Local Indices }
\begin{itemize}
  \item {co\_lower = lower bound of non-wrapped part of segment }
  \item {wr\_lower = lower bound of wrapped part of segment }
  \item {wrapd = width of wrapped part of segment.}
  \item {Let 'njj' denote the set of new (local) indices}
  \item {If there is no wrap around in either direction}
    \begin{itemize}      
      \item {for i = 1, imp3}
	\begin{itemize}
	  \item {for j = 1, range}
	    \begin{itemize}      
		\item {njj(i,j) = jj(i,j) - co\_lower + 1}
            \end{itemize}
         \end{itemize}
     \end{itemize}

  \item {If there is wrap around in southern direction}
    \begin{itemize}      
      \item {for i = 1, imp3}
	\begin {itemize}
	  \item {for j = 1, range}
	    \begin{itemize}      
                \item {if jj(i,j) .ge. wr\_lower
                       then njj(i,j) = jj(i,j) - wr\_lower + 1}
		\item {Else njj(i,j) = jj(i,j) + wrapd}
            \end{itemize}
         \end{itemize}
     \end{itemize}

  \item {If there is wrap around in northern direction}
    \begin{itemize}      
      \item {for i = 1, imp3}
	\begin {itemize}
	  \item {for j = 1, range}
	    \begin{itemize}      
                \item {if jj(i,j) refers to the wrapped part
                       then njj(i,j) = jj(i,j) + jmp3 - co\_lower + 1}
		\item {Else njj(i,j) = jj(i,j) - co\_lower + 1}
            \end{itemize}
         \end{itemize}
     \end{itemize}

\end{itemize}

\end{slide}

\end{document}


