Merge nodehelper lib and example codes into charm
[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 \section{Preparing Input Files}
180 The Multiblock framework reads its description of the problem
181 domain from input "block" files, which are in a Multiblock-specific format.
182 The files are named with the pattern \kw{prefix}\kw{number}.\kw{ext}, where
183 \kw{prefix} is a arbitrary string prefix you choose; \kw{number} is the 
184 number of this block (virtual processor); and \kw{ext} is either ``mblk'',
185 which contains binary data with the block coordinates, or ``bblk'', 
186 which contains ASCII data with the block's boundary conditions.
187
188 You generate these Multiblock input files using a tool called \kw{makemblock},
189 which can be found in charm/version/pgms/charm++/makemblock.
190 \kw{makemblock} can read a description of the problem domain 
191 generated by the structured meshing program Gridgen (from Pointwise)
192 in .grd and .inp format; or read a binary .msh format.
193 \kw{makemblock} divides this input domain into the number of blocks
194 you specify, then writes out .mblk and .bblk files.
195
196 For example, to divide the single binary mesh ``in1.msh'' into
197 20 pieces ``out00001.[mb]blk''..``out00020.[mb]blk'', you'd use
198 \begin{verbatim}
199         makemblock in1.msh 20 out
200 \end{verbatim}
201
202 You would then run this mesh using 20 virtual processors.
203
204
205 \section{Multiblock Framework API Reference}
206
207 The Multiblock framework is accessed from a program via a set of routines.
208 These routines are available in both C and Fortran90 versions.
209 The C versions are all functions, and always return an error code of 
210 MBLK\_SUCCESS or MBLK\_FAILURE.
211 The Fortran90 versions are all subroutines, and take an extra integer
212 parameter ``err'' which will be set to MBLK\_SUCCESS or MBLK\_FAILURE.
213
214 \subsection{Initialization}
215 All these methods should be called from the \kw{init} function by the user. 
216 The values passed to these functions are typically read from a configuration 
217 file or computed from command-line parameters.
218
219 \vspace{0.2in}
220 \function{int MBLK\_Set\_prefix(const char *prefix);}
221 \function{subroutine MBLK\_Set\_prefix(prefix,err)}
222   \args{character*, intent(in)::prefix}
223   \args{integer, intent(out)::err}
224 This function is called to set the block filename prefix. 
225 For example, if the input block files are named ``gridX00001.mblk''
226 and ``gridX00002.mblk'', the prefix is the string ``gridX''.
227
228 \vspace{0.2in}
229 \function{int MBLK\_Set\_nblocks(const int n);}
230 \function{ subroutine MBLK\_Set\_nblocks(n,err)}
231   \args{integer, intent(in)::n}
232   \args{integer, intent(out)::err}
233 This call is made to set the number of partitioned blocks to be used.
234 Each block is read from an input file and a separate \kw{driver}
235 is spawned for each.  The number of blocks determines the available
236 parallelism; so be sure to have at least as many blocks as processors.
237 We recommend using several times more blocks than processors, to ease 
238 load balancing and allow adaptive overlap of computation and communication.
239
240 Be sure to set the number of blocks equal to the number of virtual 
241 processors (+vp command-line option).
242
243 \vspace{0.2in}
244 \function{int MBLK\_Set\_dim(const int n);}
245 \function{subroutine MBLK\_Set\_dim(n, err)}
246   \args{integer, intent(in)::n}
247   \args{integer, intent(out)::err}
248 This call is made to set the number of spatial dimensions. 
249 Only three dimensional computations are currently supported.
250
251 \subsection{Utility}
252
253 \function{int MBLK\_Get\_nblocks(int* n);}
254 \function{subroutine MBLK\_Get\_nblocks(n,err)}
255   \args{integer,intent(out)::n}
256   \args{integer,intent(out)::err}
257 Get the total number of blocks in the current computation.  
258 Can only be called from the driver routine.
259
260 \vspace{0.2in}
261 \function{int MBLK\_Get\_myblock(int* m);}
262 \function{subroutine MBLK\_Get\_myblock(m,err)}
263   \args{integer,intent(out)::m }
264   \args{integer,intent(out)::err }
265 Get the id of the current block, an integer from 0 to the number of blocks minus one. 
266 Can only be called from the driver routine. 
267
268 \vspace{0.2in}
269 \function{int MBLK\_Get\_blocksize(int* dims);}
270 \function{subroutine MBLK\_Get\_blocksize(dimsm,err)}
271   \args{integer,intent(out)::dims(3)}
272   \args{integer,intent(out)::err }
273 Get the interior dimensions of the current block, in voxels. 
274 The size of the array dims should be 3, and will be filled with
275 the $i$, $j$, and $k$ dimensions of the block. 
276 Can only be called from the driver routine. 
277
278 \vspace{0.2in}
279 \function{int MBLK\_Get\_nodelocs(const int* nodedim,double *nodelocs);}
280 \function{subroutine MBLK\_Get\_blocksize(nodedim,nodelocs,err)}
281   \args{integer,intent(in)::nodedims(3)}
282   \args{double precision,intent(out)::nodedims(3,nodedims(0),nodedims(1),nodedims(2))}
283   \args{integer,intent(out)::err }
284 Get the $(x,y,z)$ locations of the nodes of the current block. 
285 The 3-array nodedim should be the number of nodes you expect,
286 which must be exactly one more than the number of interior voxels.
287
288 \begin{figure}[h]
289 \begin{center}
290 \includegraphics[width=3in]{fig/nodeloc}
291 \end{center}
292 \caption{The C node and voxel $(i,j,k)$ numbering for a 2 x 2 voxel block.
293 For the fortran numbering, add 1 to all indices.
294 Ghost voxels are omitted.}
295 \label{fig:indexing}
296 \end{figure}
297
298 You cannot obtain the locations of ghost nodes via this routine.
299 To get the locations of ghost nodes, create a node-centered
300 field containing the node locations and do an update field.
301 Can only be called from the driver routine. 
302
303 \vspace{0.2in}
304 \function{double MBLK\_Timer(void);}
305 \function{function double precision :: MBLK\_Timer()}
306
307 Return the current wall clock time, in seconds.  Resolution is
308      machine-dependent, but is at worst 10ms.
309
310 \vspace{0.2in}
311 \function{void MBLK\_Print\_block(void);}
312 \function{subroutine MBLK\_Print\_block()}
313 Print a debugging representation of the framework's information
314 about the current block.
315
316 \vspace{0.2in}
317 \function{void MBLK\_Print(const char *str);}
318 \function{subroutine MBLK\_Print(str)}
319 \args{  character*, intent(in) :: str}
320      Print the given string, prepended by the block id if called from the 
321      driver. Works on all machines; unlike \kw{printf} or
322      \kw{print *}, which may not work on all parallel machines.
323
324
325 \subsection{Internal Boundary Conditions and Block Fields}
326
327 The Multiblock framework handles the exchange of boundary values 
328 between neighboring blocks.
329 The basic mechanism to do this exchange is the {\em field}-- numeric data items 
330 associated with each cell of a block. These items must be arranged in a regular
331 3D grid; but otherwise we make no assumptions about the meaning of a
332 field.
333
334 You create a field once, with \kw{MBLK\_Create\_Field}, then pass the resulting 
335 field ID to \kw{MBLK\_Update\_Field} (which does the
336 overlapping block communication) and/or \kw{MBLK\_Reduce\_Field} (which 
337 applies a reduction over block values).
338
339
340 \vspace{0.2in}
341 \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);}
342 \function{subroutine MBLK\_Create\_Field(dimensions, isVoxel,base\_type, vec\_len, offset, dist, err)}
343   \args{integer, intent(in)  :: dimensions, isVoxel, base\_type, vec\_len, offset, dist}
344   \args{integer, intent(out) :: fid, err}
345      Creates and returns a Multiblock field ID, which can be passed to
346 \kw{MBLK\_Update\_Field} and \kw{MBLK\_Reduce\_Field.}  Can only be called from
347 \kw{driver().}  
348
349      \kw{dimensions} describes the size of the array the field is in.
350         Dimensions is itself an array of size 3, giving the $i$, $j$, and $k$ sizes.
351         The size should include the ghost regions-- i.e., pass the actual allocated
352         size of the array.
353      \kw{isVoxel} describes whether the data item is to be associated with
354         a voxel (1, a volume-centered value) or the nodes (0, a node-centered
355         value). 
356      \kw{base\_type} describes the type of each data item, one of:
357
358      \begin{itemize}
359         \item \kw{MBLK\_BYTE}-- unsigned char, INTEGER*1, or CHARACTER*1
360         \item \kw{MBLK\_INT}-- int or INTEGER*4
361         \item \kw{MBLK\_REAL}-- float or REAL*4
362         \item \kw{MBLK\_DOUBLE}-- double, DOUBLE PRECISION, or REAL*8
363      \end{itemize}
364
365      \kw{vec\_len} describes the number of data items associated with each
366      cell, an integer at least 1.
367
368      \kw{offset} is the byte offset from the start of the array to the
369      first interior cell's data items, a non-negative integer.  
370      This can be calculated using the \kw{offsetof()} function; normally with
371      \kw{offsetof(array(1,1,1),array(interiorX,interiorY,interiorZ))}.  
372      Be sure to skip over any ghost regions.
373
374      \kw{dist} is the byte offset from the first cell's data items to the
375      second, a positive integer (normally the size of the data items).
376      This can also be calculated using \kw{offsetof()}; normally with 
377         \kw{offsetof(array(1,1,1),array(2,1,1))}.
378
379      \kw{fid} is the identifier for the field that is created by the function.
380
381 \vspace{0.2in}
382 In the example below, we register a single double-precision value with
383 each voxel.  The ghost region is 2 cells deep along all sides.
384
385 \begin{alltt}
386         !In Fortran
387         double precision, allocatable :: voxData(:,:,:)
388         integer :: size(3), ni,nj,nk
389         integer :: fid, err
390
391         !Find the dimensions of the grid interior
392         MBLK_Get_blocksize(size,err);
393
394         !Add ghost region width to the interior dimensions 
395         size=size+4;  ! 4 because of the 2-deep region on both sides 
396
397         !Allocate and initialize the grid 
398         allocate(voxData(size(1),size(2),size(3)))
399         voxData=0.0
400
401         !Create a field for voxData
402         call MBLK_Create_field(&
403                &size,1, MBLK_DOUBLE,3,&
404                &offsetof(grid(1,1,1),grid(3,3,3)),&
405                &offsetof(grid(1,1,1),grid(2,1,1)),fid,err)      
406         
407
408 \end{alltt}
409      This example uses the Fortran-only helper routine \kw{offsetof}, which
410      returns the offset in bytes of memory between its two given
411      variables.  C users can use the built-in \kw{sizeof} keyword or pointer
412      arithmetic to achieve the same result.
413
414 \vspace{0.2in}
415 \function{void MBLK\_Update\_field(const int fid,int ghostwidth, void *grid);}
416 \function{subroutine MBLK\_Update\_field(fid,ghostwidth, grid,err)}
417   \args{integer, intent(in)  :: fid, ghostwidth}
418   \args{integer,intent(out) :: err}
419   \args{varies, intent(inout) :: grid}
420
421      Update the values in the ghost regions specified when the
422      field was created.  This call sends this block's interior region out,
423      and receives this block's boundary region from adjoining blocks.
424
425      Ghostwidth controls the thickness of the ghost region. To only exchange
426      one cell on the boundary, pass 1.  To exchange two cells, pass 2. 
427      To include diagonal regions, make the ghost width negative.
428      A ghost width of zero would communicate no data.
429
430 \begin{figure}[h]
431 \begin{center}
432 \includegraphics[width=2in]{fig/ghostwidth}
433 \end{center}
434 \caption{The 2D ghost cells communicated for various ghost widths.  The heavy line
435 is the block interior boundary-- this is the lower left portion of the block.}
436 \label{fig:ghostwidth}
437 \end{figure}      
438
439      \kw{MBLK\_Update\_field} can only be called from driver, and to be useful,
440      must be called from every block's driver routine.
441
442      \kw{MBLK\_Update\_field} blocks till the field has been updated.
443      After this routine returns, the given field will updated.
444      If the update was successful MBLK\_SUCCESS is returned and 
445      MBLK\_FAILURE is returned in case of error.
446
447 \vspace{0.2in}
448 \function{void MBLK\_Iupdate\_field(const int fid,int ghostwidth, void *ingrid, void* outgrid);}
449 \function{subroutine MBLK\_Iupdate\_field(fid,ghostwidth, ingrid, outgrid,err)}
450   \args{integer, intent(in)  :: fid, ghostwidth}
451   \args{integer,intent(out) :: err}
452   \args{varies,intent(in) :: ingrid}
453   \args{varies,intent(out) :: outgrid}
454      Update the values in the ghost regions which were specified when the
455      field was created. For the example above the ghost regions will be 
456      updated once for each step in the time loop.
457
458      \kw{MBLK\_Iupdate\_field} can only be called from driver, and to be useful,
459      must be called from every block's driver routine.
460
461      \kw{MBLK\_Iupdate\_field} is a non blocking call similar to MPI\_IRecv.
462      After the routine returns the update may not yet be complete; and the
463      outgrid may be in an inconsistent state.  Before using the values the 
464      status of the
465      update must be checked using \kw{MBLK\_Test\_update} or 
466      \kw{MBLK\_Wait\_update}.
467
468      There can be only one outstanding iupdate call in progress at any time.
469      
470 \vspace{0.2in}
471 \function{int MBLK\_Test\_update(int *status);}
472 \function{subroutine MBLK\_Test\_update(status,err)}
473 \args{integer, intent(out) :: status,err}
474      \kw{MBLK\_Test\_update} is a call that is used in association with 
475      \kw{MBLK\_Iupdate\_field} from the driver sub routine.  It tests whether
476       the preceeding iupdate has completed or not.
477      \kw{status} is returned as MBLK\_DONE if the update was completed or 
478       MBLK\_NOTDONE if the update is still pending.
479      Rather than looping if the update is still pending, call \kw{MBLK\_Wait\_update}
480      to relinquish the CPU.
481
482 \vspace{0.2in}
483 \function{void MBLK\_Wait\_update(void);}
484 \function{subroutine MBLK\_Wait\_update()}
485      \kw{MBLK\_Wait\_update} call is a blocking call and is used in assoication 
486      with \kw{MBLK\_Iupdate\_field} call. It blocks until the update is completed.
487
488 \vspace{0.2 in}
489 \function{void MBLK\_Reduce\_field(int fid,void *grid, void *out,int op);}
490 \function{subroutine MBLK\_Reduce\_field(fid,grid,outVal,op)}
491   \args{integer, intent(in)  :: fid,op}
492   \args{varies, intent(in) :: grid}
493   \args{varies, intent(out) :: outVal}
494      Combine a field from each block, according to op, across all blocks.
495      Only the interior values of the field will be combined; not the ghost cells.
496      After \kw{Reduce\_Field} returns, all blocks will have identical
497 values in \kw{outVal,} which must be \kw{vec\_len} copies of \kw{base\_type}.
498
499
500      May only be called from driver, and to complete, must be called
501      from every chunk's driver routine.
502
503      \kw{op} must be one of:
504
505 \begin{itemize}
506         \item \kw{MBLK\_SUM}-- each element of \kw{outVal} will be the sum 
507 of the corresponding fields of all blocks
508         \item \kw{MBLK\_MIN}-- each element of \kw{outVal} will be the 
509 smallest value among the corresponding field of all blocks
510         \item \kw{MBLK\_MAX}-- each element of \kw{outVal} will be the largest 
511 value among the corresponding field of all blocks
512 \end{itemize}
513
514 \vspace{0.2in}
515 \function{void MBLK\_Reduce(int fid,void *inVal,void *outVal,int op);}
516 \function{subroutine MBLK\_Reduce(fid,inVal,outVal,op)}
517   \args{integer, intent(in)  :: fid,op}
518   \args{varies, intent(in) :: inVal}
519   \args{varies, intent(out) :: outVal}
520      Combine a field from each block, acoording to \kw{op}, across all blocks.
521 \kw{Fid} is only used for the \kw{base\_type} and \kw{vec\_len}-- offset and
522 \kw{dist} are not used.  After this call returns, all blocks will have
523 identical values in \kw{outVal}.  Op has the same values and meaning as
524 \kw{MBLK\_Reduce\_Field}.
525      May only be called from driver, and to complete, must be called
526      from every blocks driver routine.
527
528 \vspace{0.2in}
529
530 \subsection{External Boundary Conditions}
531
532 Most problems include some sort of boundary conditions.  These conditions
533 are normally applied in the ghost cells surrounding the actual computational domain.
534 Examples of boundary conditions are imposed values, reflection walls, symmetry planes,
535 inlets, and exits.
536
537 The Multiblock framework keeps track of where boundary conditions are to be applied.
538 You register a subroutine that the framework will call to apply each type of 
539 external boundary condition.
540
541 \vspace{0.2in}
542 \function{int MBLK\_Register\_bc(const int bcnum, int ghostWidth, const MBLK\_BcFn bcfn);}
543 \function{subroutine MBLK\_Register\_bc(bcnum, ghostwidth, bcfn, err)}
544 \args{integer,intent(in) :: bcnum, ghostWidth}
545 \args{integer,intent(out) :: err}
546 \args{subroutine :: bcfn}
547 This is call is used to bind an external boundary condition subroutine,
548 written by you, to a boundary condition number.
549 \kw{MBLK\_Register\_bc} should only be called from the driver.
550
551 \begin{itemize}
552         \item \kw{bcnum} The boundry condtion number to be associated with the 
553         function.
554         \item \kw{ghostWidth} The width of the ghost cells where this boundry
555         condition is to be applied.
556         \item \kw{bcfn} The user subroutine to be called to apply this
557         boundry condition.
558 \end{itemize}
559
560 When you ask the framework to apply boundary conditions, it will call this
561 routine.  The routine should be declared like:
562
563 \begin{alltt}
564         !In Fortran
565         subroutine applyMyBC(param1,param2,start,end)
566         varies :: param1, param2
567         integer :: start(3), end(3)
568         end subroutine
569
570         /* In C */
571         void applyMyBC(void *param1,void *param2,int *start,int *end);  
572 \end{alltt}
573  
574 \kw{param1} and \kw{param2} are not used by the framework-- they are
575 passed in unmodified from \kw{MBLK\_Apply\_bc} and \kw{MBLK\_Apply\_bc\_all}.
576 \kw{param1} and \kw{param2} typically contain the block data and dimensions.  
577
578 \kw{start} and \kw{end} are 3-element arrays that give the $i$,$j$, $k$ 
579 block locations where the boundary condition
580 is to be applied.  They are both inclusive and both relative to the 
581 block interior-- you must shift them over your ghost cells.  The C versions
582 are 0-based (the first index is zero); the Fortran versions are 1-based
583 (the first index is one).  
584
585 For example, a Fortran subroutine to apply the constant value 1.0 across the
586 boundary, with a 2-deep ghost region, would be:
587
588 \begin{alltt}
589         !In Fortran
590         subroutine applyMyBC(grid,size,start,end)
591         integer :: size(3), i,j,k
592         double precision :: grid(size(1),size(2),size(3))
593         integer :: start(3), end(3)
594         start=start+2 ! Back up over ghost region
595         end=end+2
596         do i=start(1),end(1)
597         do j=start(2),end(2)
598         do k=start(3),end(3)
599             grid(i,j,k)=1.0
600         end do
601         end do
602         end do
603
604         end subroutine  
605 \end{alltt}
606
607
608 \vspace{0.2in}
609 \function{ int MBLK\_Apply\_bc(const int bcnum, void *param1,void *param2);}
610 \function{subroutine MBLK\_Apply\_bc(bcnum, param1,param2,err)}
611   \args{integer,intent(in)::bcnum}
612   \args{varies,intent(inout)::param1}
613   \args{varies,intent(inout)::param2}
614   \args{integer,intent(out)::err}
615 \kw{MBLK\_Apply\_bc} call is made to apply all boundry condition functions
616 of type \kw{bcnum} to the block.  \kw{param1} and \kw{param2} are passed unmodified
617 to the boundary condition function.
618
619 \vspace{0.2in}
620 \function{ int MBLK\_Apply\_bc\_all(void* param1, void* param2);} 
621 \function{subroutine MBLK\_Apply\_bc\_all(param1,param2, err)}
622   \args{integer,intent(out)::err}
623   \args{varies,intent(inout)::param1}
624   \args{varies,intent(inout)::param2}
625 This call is same as \kw{MBLK\_Apply\_bc} except it applies all 
626 external boundary conditions to the block.
627
628
629 \subsection{Migration}
630
631 The \charmpp\ runtime framework includes an automated, run-time load balancer,
632 which will automatically monitor the performance of your parallel program.
633 If needed, the load balancer can ``migrate'' mesh chunks from heavily-loaded
634 processors to more lightly-loaded processors, improving the load balance and
635 speeding up the program.  For this to be useful, pass the \kw{+vpN} argument
636 with a larger number of blocks \kw{N} than processors
637 Because this is somewhat involved, you may refrain from calling 
638 \kw{MBLK\_Migrate} and migration will never take place.
639
640 The runtime system can automatically move your thread stack to the new
641 processor, but you must write a PUP function to move any global or
642 heap-allocated data to the new processor (global data is declared at file scope
643 or \kw{static} in C and \kw{COMMON} in Fortran77; heap allocated data comes
644 from C \kw{malloc}, C++ \kw{new}, or Fortran90 \kw{ALLOCATE}).  A PUP
645 (Pack/UnPack) function performs both packing (converting heap data into a
646 message) and unpacking (converting a message back into heap data).  All your
647 global and heap data must be collected into a single block (\kw{struct} in C;
648 user-defined \kw{TYPE} in Fortran) so the PUP function can access it all.
649
650 Your PUP function will be passed a pointer to your heap data block and a
651 special handle called a ``pupper'', which contains the network message to be
652 sent.  Your PUP function returns a pointer to your heap data block.  In a PUP
653 function, you pass all your heap data to routines named \kw{pup\_type}, where
654 type is either a basic type (such as int, char, float, or double) or an array
655 type (as before, but with a ``s'' suffix).  Depending on the direction of
656 packing, the pupper will either read from or write to the values you pass--
657 normally, you shouldn't even know which.  The only time you need to know the
658 direction is when you are leaving a processor or just arriving.
659 Correspondingly, the pupper passed to you may be deleting (indicating that you
660 are leaving the processor, and should delete your heap storage after packing),
661 unpacking (indicating you've just arrived on a processor, and should allocate
662 your heap storage before unpacking), or neither (indicating the system is
663 merely sizing a buffer, or checkpointing your values).
664
665 PUP functions are much easier to write than explain-- a simple C heap block
666 and the corresponding PUP function is:
667
668 \begin{alltt}
669      typedef struct {
670        int n1;/*Length of first array below*/
671        int n2;/*Length of second array below*/
672        double *arr1; /*Some doubles, allocated on the heap*/
673        int *arr2; /*Some ints, allocated on the heap*/
674      } my_block;
675  
676      my_block *pup_my_block(pup_er p,my_block *m)
677      {
678        if (pup_isUnpacking(p)) m=malloc(sizeof(my_block));
679        pup_int(p,\&m->n1);
680        pup_int(p,\&m->n2);
681        if (pup_isUnpacking(p)) {
682          m->arr1=malloc(m->n1*sizeof(double));
683          m->arr2=malloc(m->n2*sizeof(int));
684        }
685        pup_doubles(p,m->arr1,m->n1);
686        pup_ints(p,m->arr2,m->n2);
687        if (pup_isDeleting(p)) {
688          free(m->arr1);
689          free(m->arr2);
690          free(m);
691        }
692        return m;
693      }
694 \end{alltt}
695
696 This single PUP function can be used to copy the \kw{my\_block} data into a
697 message buffer and free the old heap storage (deleting pupper); allocate
698 storage on the new processor and copy the message data back (unpacking pupper);
699 or save the heap data for debugging or checkpointing.
700
701 A Fortran block TYPE and corresponding PUP routine is as follows:
702
703 \begin{alltt}
704      MODULE my_block_mod
705        TYPE my_block
706          INTEGER :: n1,n2x,n2y
707          REAL*8, POINTER, DIMENSION(:) :: arr1
708          INTEGER, POINTER, DIMENSION(:,:) :: arr2
709        END TYPE
710      END MODULE
711  
712      SUBROUTINE pup_my_block(p,m)
713        IMPLICIT NONE
714        USE my_block_mod
715        USE pupmod
716        INTEGER :: p
717        TYPE(my_block) :: m
718        call pup_int(p,m%n1)
719        call pup_int(p,m%n2x)
720        call pup_int(p,m%n2y)
721        IF (pup_isUnpacking(p)) THEN
722          ALLOCATE(m%arr1(m%n1))
723          ALLOCATE(m%arr2(m%n2x,m%n2y))
724        END IF
725        call pup_doubles(p,m%arr1,m%n1)
726        call pup_ints(p,m%arr2,m%n2x*m%n2y)
727        IF (pup_isDeleting(p)) THEN
728          DEALLOCATE(m%arr1)
729          DEALLOCATE(m%arr2)
730        END IF
731      END SUBROUTINE
732 \end{alltt}
733
734 \function{int MBLK\_Register(void *block, MBLK\_PupFn pup\_ud, int* rid)}
735 \function{subroutine MBLK\_Register(block,pup\_ud, rid)}
736     \args{integer, intent(out)::rid}
737     \args{TYPE(varies), POINTER :: block}
738     \args{SUBROUTINE :: pup\_ud}
739      Associates the given data block and PUP function.  Returns a block
740      ID, which can be passed to \kw{MBLK\_Get\_registered} later.  Can only be
741      called from driver.  It returns MBLK\_SUCESS if the call was successful
742      and MBLK\_FAILURE in case of error. For the declarations above, you call
743      \kw{MBLK\_Register} as:
744
745 \begin{alltt}
746           /*C/C++ driver() function*/
747           int myId, err;
748           my_block *m=malloc(sizeof(my_block));
749           err =MBLK_Register(m,(MBLK_PupFn)pup_my_block,&rid);
750  
751           !- Fortran driver subroutine
752           use my_block_mod
753           interface
754             subroutine pup_my_block(p,m)
755               use my_block_mod
756               INTEGER :: p
757               TYPE(my_block) :: m
758             end subroutine
759           end interface
760           TYPE(my_block) :: m
761           INTEGER :: myId,err
762           MBLK_Register(m,pup_my_block,myId,err)
763 \end{alltt}
764
765      Note that Fortran blocks must be allocated on the stack in driver;
766      while C/C++ blocks may be allocated on the heap.
767 \vspace{0.2in}
768
769 \function{void MBLK\_Migrate()}
770 \function{subroutine MBLK\_Migrate()}
771      Informs the load balancing system that you are ready to be
772      migrated, if needed.  If the system decides to migrate you, the
773      PUP function passed to \kw{MBLK\_Register} will be called with a sizing
774      pupper, then a packing, deleting pupper.  Your stack (and pupped
775      data) will then be sent to the destination machine, where your PUP
776      function will be called with an unpacking pupper.  \kw{MBLK\_Migrate}
777      will then return, whereupon you should call \kw{MBLK\_Get\_registered} to
778      get your unpacked data block.  Can only be called from driver.
779
780
781
782 \function{int MBLK\_Get\_Userdata(int n, void** block)}
783      Return your unpacked userdata after migration-- that is, the
784      return value of the unpacking call to your PUP function.  Takes
785      the userdata ID returned by \kw{MBLK\_Register}.  Can be called from
786      driver at any time.
787
788      Since Fortran blocks are always allocated on the stack, the system
789      migrates them to the same location on the new processor, so no
790      \kw{Get\_registered} call is needed from Fortran.
791 \input{index}
792 \end{document}