Use more sensible relative path to get to doc/ directory.
[charm.git] / doc / mblock / manual.tex
1 \documentclass[10pt]{article}
2 \usepackage{../pplmanual}
3 \input{../pplmanual}
4
5 \makeindex
6
7 \title{\charmpp\\ Multiblock Framework\\ Manual}
8 \version{1.0}
9 \credits{
10 This version of \charmpp{} Multiblock Framework was developed
11 by Orion Lawlor and Milind Bhandarkar.
12 }
13
14 \begin{document}
15
16 \maketitle
17
18 \section{Motivation}
19 A large class of problems can be solved by first decomposing the
20 problem domain into a set of structured grids.  For simplicity,
21 each structured grid is often made rectangular, when it is called a {\em block}.
22 These blocks may face one another or various parts of the outside world,
23 and taken together comprise a {\em multiblock computation}.
24
25 There are two main types of multiblock computations-- implicit and explicit.
26 In an implicit computation a global matrix, which represents the entire
27 problem domain, is formed and solved.  Implicit computations require
28 a fast sparse matrix solver, and are typically used for steady-state
29 problems.  In an explicit computation, the solution proceeds locally,
30 computing new values based on the values of nearby points. Explict
31 computations often have stability criteria, and are typically used for
32 time-dependent problems.
33
34 The \charmpp{} multiblock framework allows you to write a parallel 
35 explicit multiblock program,
36 in C or Fortran 90, by concentrating on what happens to a single block
37 of the domain.  Boundary condition housekeeping and ``ghost cell'' exchange
38 are all handled transparently by the framework.
39 Using the multiblock framework also allows you to take advantage of all the
40 features of \charmpp, including adaptive computation and communication 
41 overlap, run-time load balancing,  performance
42 monitoring and visualization, and checkpoint/restart, with no additional
43 effort.
44
45
46 \section{Introduction/Terminology}
47 A {\em block} is a distorted rectangular grid that represents a 
48 portion of the problem domain.  A volumetric cell in the grid
49 is called a {\em voxel}.  Each exterior side of a block
50 is called a {\em face}. Each face may consist of several 
51 rectangular {\em patches}, which all abut the same block 
52 and experience the same boundary conditions.
53
54 \begin{figure}[h]
55 \begin{center}
56 \includegraphics[width=3in]{fig/terminology}
57 \end{center}
58 \caption{Terminology used by the framework.}
59 \label{fig:terminology}
60 \end{figure}
61
62 For example, Figure~\ref{fig:terminology} shows a 3D
63 4x8x7-voxel block, with a face and 6x3 patch indicated.
64
65 The computational domain is tiled with such blocks, which 
66 are required to be conformal-- the voxels must match exactly.
67 The blocks need not be the same size or orientation, however,
68 as illustrated in the 2D domain of Figure~\ref{fig:decompose}.
69
70 \begin{figure}[h]
71 \begin{center}
72 \includegraphics[width=4in]{fig/decompose}
73 \end{center}
74 \caption{A 2D domain decomposed into three blocks: A (5x3), B (3x6), 
75 and C (5x4). Also shows the computation as seen from block A.}
76 \label{fig:decompose}
77 \end{figure}
78
79 Figure~\ref{fig:decompose} also shows the computation 
80 from the point of view of block A, which has two external
81 boundary conditions (on the left and top sides) and two 
82 ``internal'' boundary conditions (on the right and bottom sides).
83 During the computation, the external boundary conditions can
84 be imposed independent of any other blocks; while the internal 
85 boundary conditions must be obtained from the other blocks.
86
87 To simplify the computation on the interior, these boundary conditions
88 are typically written into special extra ``ghost'' (or dummy) cells
89 around the outside of the real interior cells.  The array indexing
90 for these ghost cells is illustrated in Figure~\ref{fig:indexing}.
91
92 \begin{figure}[h]
93 \begin{center}
94 \includegraphics[width=2in]{fig/indexing}
95 \end{center}
96 \caption{The ghost cells around a 5x3-voxel 2D block}
97 \label{fig:indexing}
98 \end{figure}
99
100
101 The Multiblock framework manages all the boundary conditions--
102 both internal and external.  Internal boundary conditions are 
103 sent across processors, and require you to register the data 
104 ``fields'' you wish exchanged.  External boundary conditions are
105 not communicated, but require you to register a function to
106 apply that boundary condition to your data.  Either type of boundary 
107 condition can have arbitrary thickness.
108
109 Finally, the Multiblock framework manages nothing {\em but} boundary conditions.
110 The rest of the computation, such as deciding on and implementing
111 timestepping, stencils, numerics, and interpolation schemes are all 
112 left up to the user.  
113
114 \section{Input Files}
115 The Multiblock framework reads, in parallel, a partitioned set of 
116 blocks from block input files.  Each block consists of a file
117 with extension ``.mblk'' for the interior data (grid coordinates
118 and initial conditions) and ``.bblk'' for the boundary condition data
119 (patches where boundaries should be applied).
120
121 These block files are generated with a separate, offline tool
122 called ``makemblock'', which is documented elsewhere.
123
124 \section{Structure of a Multiblock Framework Program}
125
126 A Multiblock framework program consists of several subroutines: \kw{init}, \kw{driver},\kw{finalize}, and external boundary condition subroutines.  
127
128 \kw{init} and \kw{finalize} are called by the Multiblock framework only on the first processor -- these routines typically do specialized I/O, startup and shutdown tasks.  
129
130 A separate \kw{driver} subroutine runs for each block, and does the main work of the program.  Because there may be several blocks per processor, several \kw{driver} routines may be executing as threads simultaniously.  
131
132 The boundary condition subroutines are called by the framework after a 
133 request from \kw{driver}.
134
135
136 \begin{alltt}
137      subroutine init
138           read configuration data
139      end subroutine
140
141      subroutine bc1
142           apply first type of boundary condition
143      end subroutine bc1
144
145      subroutine bc2
146           apply second type of boundary condition
147      end subroutine bc2
148
149      subroutine driver
150           allocate and initialize the grid
151           register boundary condition subroutines bc1 and bc2
152           time loop
153                apply external boundary conditions
154                apply internal boundary conditions
155                perform serial internal computation
156           end time loop
157      end subroutine
158
159      subroutine finalize
160            write results
161      end subroutine
162 \end{alltt}
163
164
165
166 \section{Compilation and Execution}
167
168 A Multiblock framework program is a \charmpp\ program, so you must begin by
169 downloading the latest source version of \charmpp\ from
170 {\tt http://charm.cs.uiuc.edu/}.  Build the source with 
171 {\tt ./build MBLOCK version} or {\tt cd} into the build directory, 
172 {\tt version/tmp}, and type {\tt make MBLOCK}.
173 To compile a MULTIBLOCK program, pass the {\tt -language mblock} (for C) or 
174 {\tt -language mblockf} (for Fortran) option to {\tt charmc}.
175
176 In a charm installation, see charm/version/pgms/charm++/mblock/
177 for example and test programs.
178
179
180 \section{Multiblock Framework API Reference}
181
182 The Multiblock framework is accessed from a program via a set of routines.
183 These routines are available in both C and Fortran90 versions.
184 The C versions are all functions, and always return an error code of 
185 MBLK\_SUCCESS or MBLK\_FAILURE.
186 The Fortran90 versions are all subroutines, and take an extra integer
187 parameter ``err'' which will be set to MBLK\_SUCCESS or MBLK\_FAILURE.
188
189 \subsection{Initialization}
190 All these methods should be called from the \kw{init} function by the user. 
191 The values passed to these functions are typically read from a configuration 
192 file or computed from command-line parameters.
193
194 \vspace{0.2in}
195 \function{int MBLK\_Set\_prefix(const char *prefix);}
196 \function{subroutine MBLK\_Set\_prefix(prefix,err)}
197   \args{character*, intent(in)::prefix}
198   \args{integer, intent(out)::err}
199 This function is called to set the block filename prefix. 
200 For example, if the input block files are named ``gridX000001.mblk''
201 and ``gridX000002.mblk'', the prefix is the string ``gridX''.
202
203 \vspace{0.2in}
204 \function{int MBLK\_Set\_nblocks(const int n);}
205 \function{ subroutine MBLK\_Set\_nblocks(n,err)}
206   \args{integer, intent(in)::n}
207   \args{integer, intent(out)::err}
208 This call is made to set the number of partitioned blocks to be used.
209 Each block is read from an input file and a separate \kw{driver}
210 is spawned for each.  The number of blocks determines the available
211 parallelism; so be sure to have at least as many blocks as processors.
212 We recommend using several times more blocks than processors, to ease 
213 load balancing and allow adaptive overlap of computation and communication.
214
215 \vspace{0.2in}
216 \function{int MBLK\_Set\_dim(const int n);}
217 \function{subroutine MBLK\_Set\_dim(n, err)}
218   \args{integer, intent(in)::n}
219   \args{integer, intent(out)::err}
220 This call is made to set the number of spatial dimensions. 
221 Only three dimensional computations are currently supported.
222
223 \subsection{Utility}
224
225 \function{int MBLK\_Get\_nblocks(int* n);}
226 \function{subroutine MBLK\_Get\_nblocks(n,err)}
227   \args{integer,intent(out)::n}
228   \args{integer,intent(out)::err}
229 Get the total number of blocks in the current computation.  
230 Can only be called from the driver routine.
231
232 \vspace{0.2in}
233 \function{int MBLK\_Get\_myblock(int* m);}
234 \function{subroutine MBLK\_Get\_myblock(m,err)}
235   \args{integer,intent(out)::m }
236   \args{integer,intent(out)::err }
237 Get the id of the current block, an integer from 0 to the number of blocks minus one. 
238 Can only be called from the driver routine. 
239
240 \vspace{0.2in}
241 \function{int MBLK\_Get\_blocksize(int* dims);}
242 \function{subroutine MBLK\_Get\_blocksize(dimsm,err)}
243   \args{integer,intent(out)::dims(3)}
244   \args{integer,intent(out)::err }
245 Get the interior dimensions of the current block, in voxels. 
246 The size of the array dims should be 3, and will be filled with
247 the $i$, $j$, and $k$ dimensions of the block. 
248 Can only be called from the driver routine. 
249
250 \vspace{0.2in}
251 \function{double MBLK\_Timer(void);}
252 \function{function double precision :: MBLK\_Timer()}
253
254 Return the current wall clock time, in seconds.  Resolution is
255      machine-dependent, but is at worst 10ms.
256
257 \vspace{0.2in}
258 \function{void MBLK\_Print\_block(void);}
259 \function{subroutine MBLK\_Print\_block()}
260 Print a debugging representation of the framework's information
261 about the current block.
262
263 \vspace{0.2in}
264 \function{void MBLK\_Print(const char *str);}
265 \function{subroutine MBLK\_Print(str)}
266 \args{  character*, intent(in) :: str}
267      Print the given string, prepended by the block id if called from the 
268      driver. Works on all machines; unlike \kw{printf} or
269      \kw{print *}, which may not work on all parallel machines.
270
271
272 \subsection{Internal Boundary Conditions and Block Fields}
273
274 The Multiblock framework handles the exchange of boundary values 
275 between neighboring blocks.
276 The basic mechanism to do this exchange is the {\em field}-- numeric data items 
277 associated with each cell of a block. These items must be arranged in a regular
278 3D grid; but otherwise we make no assumptions about the meaning of a
279 field.
280
281 You create a field once, with \kw{MBLK\_Create\_Field}, then pass the resulting 
282 field ID to \kw{MBLK\_Update\_Field} (which does the
283 overlapping block communication) and/or \kw{MBLK\_Reduce\_Field} (which 
284 applies a reduction over block values).
285
286
287 \vspace{0.2in}
288 \function{int MBLK\_Create\_Field(int *dimensions,int isVoxel,const int base\_type,const int vec\_len,const int offset,const int dist, int *fid);}
289 \function{subroutine MBLK\_Create\_Field(dimensions, isVoxel,base\_type, vec\_len, offset, dist, err)}
290   \args{integer, intent(in)  :: dimensions, isVoxel, base\_type, vec\_len, offset, dist}
291   \args{integer, intent(out) :: fid, err}
292      Creates and returns a Multiblock field ID, which can be passed to
293 \kw{MBLK\_Update\_Field} and \kw{MBLK\_Reduce\_Field.}  Can only be called from
294 \kw{driver().}  
295
296      \kw{dimensions} describes the size of the array the field is in.
297         Dimensions is itself an array of size 3, giving the $i$, $j$, and $k$ sizes.
298         The size should include the ghost regions-- i.e., pass the actual allocated
299         size of the array.
300      \kw{isVoxel} describes whether the data item is to be associated with
301         a voxel (1, a volume-centered value) or the grid corners (0, a corner-centered
302         value). 
303      \kw{base\_type} describes the type of each data item, one of:
304
305      \begin{itemize}
306         \item \kw{MBLK\_BYTE}-- unsigned char, INTEGER*1, or CHARACTER*1
307         \item \kw{MBLK\_INT}-- int or INTEGER*4
308         \item \kw{MBLK\_REAL}-- float or REAL*4
309         \item \kw{MBLK\_DOUBLE}-- double, DOUBLE PRECISION, or REAL*8
310      \end{itemize}
311
312      \kw{vec\_len} describes the number of data items associated with each
313      cell, an integer at least 1.
314
315      \kw{offset} is the byte offset from the start of the array to the
316      first interior cell's data items, a non-negative integer.  
317      This can be calculated using the \kw{offsetof()} function; normally with
318      \kw{offsetof(array(1,1,1),array(interiorX,interiorY,interiorZ))}.  
319      Be sure to skip over any ghost regions.
320
321      \kw{dist} is the byte offset from the first cell's data items to the
322      second, a positive integer (normally the size of the data items).
323      This can also be calculated using \kw{offsetof()}; normally with 
324         \kw{offsetof(array(1,1,1),array(2,1,1))}.
325
326      \kw{fid} is the identifier for the field that is created by the function.
327
328 \vspace{0.2in}
329 In the example below, we register a single double-precision value with
330 each voxel.  The ghost region is 2 cells deep along all sides.
331
332 \begin{alltt}
333         !In Fortran
334         double precision, allocatable :: voxData(:,:,:)
335         integer :: size(3), ni,nj,nk
336         integer :: fid, err
337
338         !Find the dimensions of the grid interior
339         MBLK_Get_blocksize(size,err);
340
341         !Add ghost region width to the interior dimensions 
342         size=size+4;  ! 4 because of the 2-deep region on both sides 
343
344         !Allocate and initialize the grid 
345         allocate(voxData(size(1),size(2),size(3)))
346         voxData=0.0
347
348         !Create a field for voxData
349         call MBLK_Create_field(&
350                &size,1, MBLK_DOUBLE,3,&
351                &offsetof(grid(1,1,1),grid(3,3,3)),&
352                &offsetof(grid(1,1,1),grid(2,1,1)),fid,err)      
353         
354
355 \end{alltt}
356      This example uses the Fortran-only helper routine \kw{offsetof}, which
357      returns the offset in bytes of memory between its two given
358      variables.  C users can use the built-in \kw{sizeof} keyword or pointer
359      arithmetic to achieve the same result.
360
361 \vspace{0.2in}
362 \function{void MBLK\_Update\_field(const int fid,int ghostwidth, void *grid);}
363 \function{subroutine MBLK\_Update\_field(fid,ghostwidth, grid,err)}
364   \args{integer, intent(in)  :: fid, ghostwidth}
365   \args{integer,intent(out) :: err}
366   \args{varies, intent(inout) :: grid}
367
368      Update the values in the ghost regions specified when the
369      field was created.  This call sends this block's interior region out,
370      and receives this block's boundary region from adjoining blocks.
371
372      Ghostwidth controls the thickness of the ghost region. To only exchange
373      one cell on the boundary, pass 1.  To exchange two cells, pass 2. 
374      To include diagonal regions, make the ghost width negative.
375      A ghost width of zero would communicate no data.
376
377 \begin{figure}[h]
378 \begin{center}
379 \includegraphics[width=2in]{fig/ghostwidth}
380 \end{center}
381 \caption{The 2D ghost cells communicated for various ghost widths.  The heavy line
382 is the block interior boundary-- this is the lower left portion of the block.}
383 \label{fig:ghostwidth}
384 \end{figure}      
385
386      \kw{MBLK\_Update\_field} can only be called from driver, and to be useful,
387      must be called from every block's driver routine.
388
389      \kw{MBLK\_Update\_field} blocks till the field has been updated.
390      After this routine returns, the given field will updated.
391      If the update was successful MBLK\_SUCCESS is returned and 
392      MBLK\_FAILURE is returned in case of error.
393
394 \vspace{0.2in}
395 \function{void MBLK\_Iupdate\_field(const int fid,int ghostwidth, void *ingrid, void* outgrid);}
396 \function{subroutine MBLK\_Iupdate\_field(fid,ghostwidth, ingrid, outgrid,err)}
397   \args{integer, intent(in)  :: fid, ghostwidth}
398   \args{integer,intent(out) :: err}
399   \args{varies,intent(in) :: ingrid}
400   \args{varies,intent(out) :: outgrid}
401      Update the values in the ghost regions which were specified when the
402      field was created. For the example above the ghost regions will be 
403      updated once for each step in the time loop.
404
405      \kw{MBLK\_Iupdate\_field} can only be called from driver, and to be useful,
406      must be called from every block's driver routine.
407
408      \kw{MBLK\_Iupdate\_field} is a non blocking call similar to MPI\_IRecv.
409      After the routine returns the update may not yet be complete; and the
410      outgrid may be in an inconsistent state.  Before using the values the 
411      status of the
412      update must be checked using \kw{MBLK\_Test\_update} or 
413      \kw{MBLK\_Wait\_update}.
414
415      There can be only one outstanding iupdate call in progress at any time.
416      
417 \vspace{0.2in}
418 \function{int MBLK\_Test\_update(int *status);}
419 \function{subroutine MBLK\_Test\_update(status,err)}
420 \args{integer, intent(out) :: status,err}
421      \kw{MBLK\_Test\_update} is a call that is used in association with 
422      \kw{MBLK\_Iupdate\_field} from the driver sub routine.  It tests whether
423       the preceeding iupdate has completed or not.
424      \kw{status} is returned as MBLK\_DONE if the update was completed or 
425       MBLK\_NOTDONE if the update is still pending.
426      Rather than looping if the update is still pending, call \kw{MBLK\_Wait\_update}
427      to relinquish the CPU.
428
429 \vspace{0.2in}
430 \function{void MBLK\_Wait\_update(void);}
431 \function{subroutine MBLK\_Wait\_update()}
432      \kw{MBLK\_Wait\_update} call is a blocking call and is used in assoication 
433      with \kw{MBLK\_Iupdate\_field} call. It blocks until the update is completed.
434
435 \vspace{0.2 in}
436 \function{void MBLK\_Reduce\_field(int fid,void *grid, void *out,int op);}
437 \function{subroutine MBLK\_Reduce\_field(fid,grid,outVal,op)}
438   \args{integer, intent(in)  :: fid,op}
439   \args{varies, intent(in) :: grid}
440   \args{varies, intent(out) :: outVal}
441      Combine a field from each block, according to op, across all blocks.
442      Only the interior values of the field will be combined; not the ghost cells.
443      After \kw{Reduce\_Field} returns, all blocks will have identical
444 values in \kw{outVal,} which must be \kw{vec\_len} copies of \kw{base\_type}.
445
446
447      May only be called from driver, and to complete, must be called
448      from every chunk's driver routine.
449
450      \kw{op} must be one of:
451
452 \begin{itemize}
453         \item \kw{MBLK\_SUM}-- each element of \kw{outVal} will be the sum 
454 of the corresponding fields of all blocks
455         \item \kw{MBLK\_MIN}-- each element of \kw{outVal} will be the 
456 smallest value among the corresponding field of all blocks
457         \item \kw{MBLK\_MAX}-- each element of \kw{outVal} will be the largest 
458 value among the corresponding field of all blocks
459 \end{itemize}
460
461 \vspace{0.2in}
462 \function{void MBLK\_Reduce(int fid,void *inVal,void *outVal,int op);}
463 \function{subroutine MBLK\_Reduce(fid,inVal,outVal,op)}
464   \args{integer, intent(in)  :: fid,op}
465   \args{varies, intent(in) :: inVal}
466   \args{varies, intent(out) :: outVal}
467      Combine a field from each block, acoording to \kw{op}, across all blocks.
468 \kw{Fid} is only used for the \kw{base\_type} and \kw{vec\_len}-- offset and
469 \kw{dist} are not used.  After this call returns, all blocks will have
470 identical values in \kw{outVal}.  Op has the same values and meaning as
471 \kw{MBLK\_Reduce\_Field}.
472      May only be called from driver, and to complete, must be called
473      from every blocks driver routine.
474
475 \vspace{0.2in}
476
477 \subsection{External Boundary Conditions}
478
479 Most problems include some sort of boundary conditions.  These conditions
480 are normally applied in the ghost cells surrounding the actual computational domain.
481 Examples of boundary conditions are imposed values, reflection walls, symmetry planes,
482 inlets, and exits.
483
484 The Multiblock framework keeps track of where boundary conditions are to be applied.
485 You register a subroutine that the framework will call to apply each type of 
486 external boundary condition.
487
488 \vspace{0.2in}
489 \function{int MBLK\_Register\_bc(const int bcnum, int ghostWidth, const MBLK\_BcFn bcfn);}
490 \function{subroutine MBLK\_Register\_bc(bcnum, ghostwidth, bcfn, err)}
491 \args{integer,intent(in) :: bcnum, ghostWidth}
492 \args{integer,intent(out) :: err}
493 \args{subroutine :: bcfn}
494 This is call is used to bind an external boundary condition subroutine,
495 written by you, to a boundary condition number.
496 \kw{MBLK\_Register\_bc} should only be called from the driver.
497
498 \begin{itemize}
499         \item \kw{bcnum} The boundry condtion number to be associated with the 
500         function.
501         \item \kw{ghostWidth} The width of the ghost cells where this boundry
502         condition is to be applied.
503         \item \kw{bcfn} The user subroutine to be called to apply this
504         boundry condition.
505 \end{itemize}
506
507 When you ask the framework to apply boundary conditions, it will call this
508 routine.  The routine should be declared like:
509
510 \begin{alltt}
511         !In Fortran
512         subroutine applyMyBC(param1,param2,start,end)
513         varies :: param1, param2
514         integer :: start(3), end(3)
515         end subroutine
516
517         /* In C */
518         void applyMyBC(void *param1,void *param2,int *start,int *end);  
519 \end{alltt}
520  
521 \kw{param1} and \kw{param2} are not used by the framework-- they are
522 passed in unmodified from \kw{MBLK\_Apply\_bc} and \kw{MBLK\_Apply\_bc\_all}.
523 \kw{param1} and \kw{param2} typically contain the block data and dimensions.  
524
525 \kw{start} and \kw{end} are 3-element arrays that give the $i$,$j$, $k$ 
526 block locations where the boundary condition
527 is to be applied.  They are both inclusive and both relative to the 
528 block interior-- you must shift them over your ghost cells.  The C versions
529 are 0-based (the first index is zero); the Fortran versions are 1-based
530 (the first index is one).  
531
532 For example, a Fortran subroutine to apply the constant value 1.0 across the
533 boundary, with a 2-deep ghost region, would be:
534
535 \begin{alltt}
536         !In Fortran
537         subroutine applyMyBC(grid,size,start,end)
538         integer :: size(3), i,j,k
539         double precision :: grid(size(1),size(2),size(3))
540         integer :: start(3), end(3)
541         start=start+2 ! Back up over ghost region
542         end=end+2
543         do i=start(1),end(1)
544         do j=start(2),end(2)
545         do k=start(3),end(3)
546             grid(i,j,k)=1.0
547         end do
548         end do
549         end do
550
551         end subroutine  
552 \end{alltt}
553
554
555 \vspace{0.2in}
556 \function{ int MBLK\_Apply\_bc(const int bcnum, void *param1,void *param2);}
557 \function{subroutine MBLK\_Apply\_bc(bcnum, param1,param2,err)}
558   \args{integer,intent(in)::bcnum}
559   \args{varies,intent(inout)::param1}
560   \args{varies,intent(inout)::param2}
561   \args{integer,intent(out)::err}
562 \kw{MBLK\_Apply\_bc} call is made to apply all boundry condition functions
563 of type \kw{bcnum} to the block.  \kw{param1} and \kw{param2} are passed unmodified
564 to the boundary condition function.
565
566 \vspace{0.2in}
567 \function{ int MBLK\_Apply\_bc\_all(void* param1, void* param2);} 
568 \function{subroutine MBLK\_Apply\_bc\_all(param1,param2, err)}
569   \args{integer,intent(out)::err}
570   \args{varies,intent(inout)::param1}
571   \args{varies,intent(inout)::param2}
572 This call is same as \kw{MBLK\_Apply\_bc} except it applies all 
573 external boundary conditions to the block.
574
575
576 \subsection{Migration}
577
578 The \charmpp\ runtime framework includes an automated, run-time load balancer,
579 which will automatically monitor the performance of your parallel program.
580 If needed, the load balancer can ``migrate'' mesh chunks from heavily-loaded
581 processors to more lightly-loaded processors, improving the load balance and
582 speeding up the program.  For this to be useful, pass the \kw{+vpN} argument
583 with a larger number of blocks \kw{N} than processors
584 Because this is somewhat involved, you may refrain from calling 
585 \kw{MBLK\_Migrate} and migration will never take place.
586
587 The runtime system can automatically move your thread stack to the new
588 processor, but you must write a PUP function to move any global or
589 heap-allocated data to the new processor (global data is declared at file scope
590 or \kw{static} in C and \kw{COMMON} in Fortran77; heap allocated data comes
591 from C \kw{malloc}, C++ \kw{new}, or Fortran90 \kw{ALLOCATE}).  A PUP
592 (Pack/UnPack) function performs both packing (converting heap data into a
593 message) and unpacking (converting a message back into heap data).  All your
594 global and heap data must be collected into a single block (\kw{struct} in C;
595 user-defined \kw{TYPE} in Fortran) so the PUP function can access it all.
596
597 Your PUP function will be passed a pointer to your heap data block and a
598 special handle called a ``pupper'', which contains the network message to be
599 sent.  Your PUP function returns a pointer to your heap data block.  In a PUP
600 function, you pass all your heap data to routines named \kw{pup\_type}, where
601 type is either a basic type (such as int, char, float, or double) or an array
602 type (as before, but with a ``s'' suffix).  Depending on the direction of
603 packing, the pupper will either read from or write to the values you pass--
604 normally, you shouldn't even know which.  The only time you need to know the
605 direction is when you are leaving a processor or just arriving.
606 Correspondingly, the pupper passed to you may be deleting (indicating that you
607 are leaving the processor, and should delete your heap storage after packing),
608 unpacking (indicating you've just arrived on a processor, and should allocate
609 your heap storage before unpacking), or neither (indicating the system is
610 merely sizing a buffer, or checkpointing your values).
611
612 PUP functions are much easier to write than explain-- a simple C heap block
613 and the corresponding PUP function is:
614
615 \begin{alltt}
616      typedef struct {
617        int n1;/*Length of first array below*/
618        int n2;/*Length of second array below*/
619        double *arr1; /*Some doubles, allocated on the heap*/
620        int *arr2; /*Some ints, allocated on the heap*/
621      } my_block;
622  
623      my_block *pup_my_block(pup_er p,my_block *m)
624      {
625        if (pup_isUnpacking(p)) m=malloc(sizeof(my_block));
626        pup_int(p,\&m->n1);
627        pup_int(p,\&m->n2);
628        if (pup_isUnpacking(p)) {
629          m->arr1=malloc(m->n1*sizeof(double));
630          m->arr2=malloc(m->n2*sizeof(int));
631        }
632        pup_doubles(p,m->arr1,m->n1);
633        pup_ints(p,m->arr2,m->n2);
634        if (pup_isDeleting(p)) {
635          free(m->arr1);
636          free(m->arr2);
637          free(m);
638        }
639        return m;
640      }
641 \end{alltt}
642
643 This single PUP function can be used to copy the \kw{my\_block} data into a
644 message buffer and free the old heap storage (deleting pupper); allocate
645 storage on the new processor and copy the message data back (unpacking pupper);
646 or save the heap data for debugging or checkpointing.
647
648 A Fortran block TYPE and corresponding PUP routine is as follows:
649
650 \begin{alltt}
651      MODULE my_block_mod
652        TYPE my_block
653          INTEGER :: n1,n2x,n2y
654          REAL*8, POINTER, DIMENSION(:) :: arr1
655          INTEGER, POINTER, DIMENSION(:,:) :: arr2
656        END TYPE
657      END MODULE
658  
659      SUBROUTINE pup_my_block(p,m)
660        IMPLICIT NONE
661        USE my_block_mod
662        USE pupmod
663        INTEGER :: p
664        TYPE(my_block) :: m
665        call pup_int(p,m%n1)
666        call pup_int(p,m%n2x)
667        call pup_int(p,m%n2y)
668        IF (pup_isUnpacking(p)) THEN
669          ALLOCATE(m%arr1(m%n1))
670          ALLOCATE(m%arr2(m%n2x,m%n2y))
671        END IF
672        call pup_doubles(p,m%arr1,m%n1)
673        call pup_ints(p,m%arr2,m%n2x*m%n2y)
674        IF (pup_isDeleting(p)) THEN
675          DEALLOCATE(m%arr1)
676          DEALLOCATE(m%arr2)
677        END IF
678      END SUBROUTINE
679 \end{alltt}
680
681 \function{int MBLK\_Register(void *block, MBLK\_PupFn pup\_ud, int* rid)}
682 \function{subroutine MBLK\_Register(block,pup\_ud, rid)}
683     \args{integer, intent(out)::rid}
684     \args{TYPE(varies), POINTER :: block}
685     \args{SUBROUTINE :: pup\_ud}
686      Associates the given data block and PUP function.  Returns a block
687      ID, which can be passed to \kw{MBLK\_Get\_registered} later.  Can only be
688      called from driver.  It returns MBLK\_SUCESS if the call was successful
689      and MBLK\_FAILURE in case of error. For the declarations above, you call
690      \kw{MBLK\_Register} as:
691
692 \begin{alltt}
693           /*C/C++ driver() function*/
694           int myId, err;
695           my_block *m=malloc(sizeof(my_block));
696           err =MBLK_Register(m,(MBLK_PupFn)pup_my_block,&rid);
697  
698           !- Fortran driver subroutine
699           use my_block_mod
700           interface
701             subroutine pup_my_block(p,m)
702               use my_block_mod
703               INTEGER :: p
704               TYPE(my_block) :: m
705             end subroutine
706           end interface
707           TYPE(my_block) :: m
708           INTEGER :: myId,err
709           MBLK_Register(m,pup_my_block,myId,err)
710 \end{alltt}
711
712      Note that Fortran blocks must be allocated on the stack in driver;
713      while C/C++ blocks may be allocated on the heap.
714 \vspace{0.2in}
715
716 \function{void MBLK\_Migrate()}
717 \function{subroutine MBLK\_Migrate()}
718      Informs the load balancing system that you are ready to be
719      migrated, if needed.  If the system decides to migrate you, the
720      PUP function passed to \kw{MBLK\_Register} will be called with a sizing
721      pupper, then a packing, deleting pupper.  Your stack (and pupped
722      data) will then be sent to the destination machine, where your PUP
723      function will be called with an unpacking pupper.  \kw{MBLK\_Migrate}
724      will then return, whereupon you should call \kw{MBLK\_Get\_registered} to
725      get your unpacked data block.  Can only be called from driver.
726
727
728
729 \function{int MBLK\_Get\_Userdata(int n, void** block)}
730      Return your unpacked userdata after migration-- that is, the
731      return value of the unpacking call to your PUP function.  Takes
732      the userdata ID returned by \kw{MBLK\_Register}.  Can be called from
733      driver at any time.
734
735      Since Fortran blocks are always allocated on the stack, the system
736      migrates them to the same location on the new processor, so no
737      \kw{Get\_registered} call is needed from Fortran.
738 \input{index}
739 \end{document}