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