First version of the multiblock manual
[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
20
21
22 \section{Introduction/Terminology}
23
24
25
26 \section{Structure of a MultiBlock Framework Program}
27
28 A MultiBlock framework program consists of three subroutines: \kw{init}, \kw{driver}, and \kw{finalize}.  \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.  \kw{driver} is called for every chunk on every processor, and does the main work of the program.
29
30 \begin{alltt}
31      subroutine init
32           read the configuration data llike block prefix, 
33           number of blocks and dimension 
34      end subroutine
35
36      subroutine driver
37           allocate and initialize the grid
38           register boundary condition functions
39           time loop
40                MultiBlock computations
41                update shared node fields
42                more MultiBlock computations
43           end time loop
44      end subroutine
45
46      subroutine finalize
47            write results
48      end subroutine
49 \end{alltt}
50
51 \section{Compilation and Execution}
52
53 A Multiblock framework program is a \charmpp\ program, so you must begin by
54 downloading the latest source version of \charmpp\ from
55 {\tt http://charm.cs.uiuc.edu/}.  Build the source with 
56 {\tt ./build MBLOCK version} or {\tt cd} into the build directory, 
57 {\tt version/tmp}, and type {\tt make MBLOCK}.
58 To compile a MULTIBLOCK program, pass the {\tt -language mblock} (for C) or 
59 {\tt -language mblockf} (for Fortran) option to {\tt charmc}.
60
61 In a charm installation, see charm/version/pgms/charm++/mblock/
62 for example and test programs.
63
64
65 \section{MultiBlock Framework API Reference}
66
67 \subsection{Initialization}
68 All these methods should be called from the \kw{init} function by the user.The values
69 passed to these functions are typically read from a file.
70 \function{int MBLK\_Set\_prefix(const char *prefix);}
71 \function{subroutine MBLK\_Set\_prefix(prefix,err)}
72 \args{ Character, intent(in)::prefix}
73 \args{integer, intent(out)::err}
74
75 This function is called to set the prefix. It returns MBLK\_SUCCESS in case of 
76 success else it returns MBLK\_FAILURE.
77
78 \vspace{0.2in}
79
80 \function{int MBLK\_Set\_nblocks(const int n);}
81 \function{ subroutine MBLK\_Set\_nblocks(n,err)}
82 \args{integer, intent(in)::n}
83 \args{integer, intent(out)::err}
84 This call is made to set the number of blocks to be used in the application.It 
85 returns MBLK\_SUCCESS in case of success else it returns MBLK\_FAILURE.
86 \vspace{0.2in}
87
88 \function{int MBLK\_Set\_dim(const int n);}
89 \function{subroutine MBLK\_Set\_dim(n, err)}
90 \args{integer, intent(in)::n}
91 \args{integer, intent(out):: err}
92 This call is made to set the dimension. It returns MBLK\_SUCCESS in case of 
93 success else it returns MBLK\_FAILURE.
94
95 \subsection{Utility}
96
97 \function{int MBLK\_Get\_nblocks(int* n);}
98 \function{subroutine MBLK\_Get\_nblocks(n,err)}
99 \args{  integer,intent(out)::n
100         integer,intent(out)::err}
101      Get the number of blocks in the current computation.  Can
102      only be called from the driver routine. Returns MBLK\_SUCCESS in case of
103      success and MBLK\_FAILURE in case of error.
104 \vspace{0.2in}
105
106 \function{int MBLK\_Get\_myblock(int* m);}
107 \function{subroutine MBLK\_Get\_myblock(m,err)}
108 \args{  integer,intent(out)::m
109         integer,intent(out)::err }
110
111      Get the id of the current block. Can only be called from the driver 
112      routine. Returns MBLK\_SUCCESS in case of success and MBLK\_FAILURE in 
113      case of error.
114 \vspace{0.2in}
115
116 \function{int MBLK\_Get\_blocksize(int* dims);}
117 \function{subroutine MBLK\_Get\_blocksize(dimsm,err)}
118 \args{  integer,intent(out)::dims(3)
119         integer,intent(out)::err }
120      Get the Interior dimensions, in voxels. The size of the array dims should
121      be 3. Can only be called from the driver routine. Returns MBLK\_SUCCESS 
122      in case of success and MBLK\_FAILURE in case of error.
123 \vspace{0.2in}
124
125 \function{double MBLK\_Timer(void);}
126 \function{function double precision :: MBLK\_Timer()}
127
128      Return the current wall clock time, in seconds.  Resolution is
129      machine-dependent, but is at worst 10ms.
130 \vspace{0.2in}
131
132 \function{void MBLK\_Print\_block(void);}
133 \function{subroutine MBLK\_Print\_block()}
134
135      Print a debugging representation of the current block.
136      Prints the entire blocks array, and data associated with
137      each block.
138 \vspace{0.2in}
139
140 \function{void MBLK\_Print(const char *str);}
141 \function{subroutine MBLK\_Print(str)}
142 \args{  character*, intent(in) :: str}
143
144      Print the given string, prepended by the block id if called from the 
145      driver. Works on all machines; unlike \kw{printf} or
146      \kw{print *}, which may not work on all parallel machines.
147
148
149 \subsection{Block Fields}
150
151
152 The MultiBlock framework handles the updating of the values of blocks.
153 The basic mechanism to do this update is the field-- numeric data items associated with each block. We make no assumptions about the meaning of the block data, and allow various data types and non-communicated data associated with each
154 block.  To do this, the framework must be able to find the data items
155 associated with each block in memory.
156
157 Each field represents a (set of) block data items stored in a contiguous array,
158 indexed by block number.  You create a field once, with \kw{MBLK\_Create\_Field}, then pass the resulting field ID to \kw{MBLK\_Update\_Field} (which does the
159 overlapping block communication) and/or \kw{MBLK\_Reduce\_Field} (which applies a reduction over block values).
160 \vspace{0.2in}
161
162 \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);}
163 \function{subroutine MBLK\_Create\_Field(dimensions, isVoxel,base\_type, vec\_len, offset, dist, err)}
164   \args{integer, intent(in)  :: dimensions, isVoxel, base\_type, vec\_len, offset, dist
165         integer, intent(out) :: fid, err}
166
167      Creates and returns a MultiBlock field ID, which can be passed to
168 \kw{MBLK\_Update\_Field} and \kw{MBLK\_Reduce\_Field.}  Can only be called from
169 \kw{driver().}  A field is a range of values associated with each local block--
170 the Multiblock framework uses the information you pass to find the values associated with overlapping blocks (for \kw{MBLk\_Update\_Field}) and primary blocks (for \kw{MBLK\_Reduce\_Field}). Returns MBLK\_SUCCESS in case of success and MBLK\_FAILURE in case of error.
171
172      \kw{dimensions} describes the number of dimensions should be in array of size 3
173      \kw{isVoxel} describes whether the dimenstions passed in are in voxel(=11) or not(= 0).
174      \kw{base\_type} describes the kind of data item associated with each
175      node, one of:
176
177      \begin{itemize}
178         \item \kw{MBLk\_BYTE}-- unsigned char, INTEGER*1, or CHARACTER*1
179         \item \kw{MBLK\_INT}-- int or INTEGER*4
180         \item \kw{MBLK\_REAL}-- float or REAL*4
181         \item \kw{MBLK\_DOUBLE}-- double, DOUBLE PRECISION, or REAL*8
182      \end{itemize}
183
184      \kw{vec\_len} describes the number of data items associated with each
185      node, an integer at least 1.
186
187      \kw{offset} is the byte offset from the start of the nodes array to the
188      data items, a non-negative integer.
189
190      \kw{dist} is the byte offset from the first node's data item to the
191      second, a positive integer.
192      \kw{fid} is the identifier for the field that is created by the function.
193      \kw{err} is returned as  MBLK\_SUCCESS in case of success otherwise MBLK\_FAILURE is returned.
194 \vspace{0.2in}
195
196      For example, if each block has a 3D grid, over which we are performing Successive over relaxation.You can register the ghost regions for update with:
197 \begin{alltt}
198         !In Fortran
199         
200         integer :: size(3), ni,nj,nk
201         integer :: si,sj,sk
202         integer :: fid, err
203         integer, parameter::ghostwidth=1;       
204         !Find the dimensions of the grid, from the framework to allocateit
205         MBLk_Get_blocksize(size,err);
206
207         !Add ghost region width to the dimensions obtained from the framework
208         ni=size(1)+2*ghostWidth; 
209         nj=size(2)+2*ghostWidth; 
210         nk=size(3)+2*ghostWidth;
211         si=1+ghostWidth; sj=1+ghostWidth; sk=1+ghostWidth;
212
213         ...allocate and initialize the grid 
214
215         !Create the field that needs to be updated
216         size(1)=ni; size(2)=nj; size(3)=nk;
217         call MBLK_Create_field(&
218                &size,1, MBLK_DOUBLE,1,&
219                &offsetof(grid(1,1,1),grid(si,sj,sk)),&
220                &offsetof(grid(1,1,1),grid(2,1,1)),fid,err)      
221         
222
223 \end{alltt}
224      This example uses the Fortran-only helper routine \kw{offsetof}, which
225      returns the offset in bytes of memory between its two given
226      variables.  The C version uses pointer arithmetic to achieve the
227      same result.
228 \vspace{0.2in}
229      
230
231 \function{void MBLK\_Update\_field(const int fid,int ghostwidth, void *grid);}
232 \function{subroutine MBLK\_Update\_field(fid,ghostwidth, grid,err)}
233   \args{integer, intent(in)  :: fid, ghostwidth}
234   \args{integer,intent(out) :: err}
235   \args{varies, intent(inout) :: grid}
236
237      Update the values in the ghost regions which specified when the
238      field was created. For the example above the ghost regions will be 
239      updated once for each step in the time loop.
240
241      \kw{MBLK\_Update\_field} can only be called from driver, and to be useful,
242      must be called from every block's driver routine.
243
244      \kw{MBLK\_Update\_field} blocks till the field has been updated.
245      After this routine returns, the given field will updated.
246      If the update was successful MBLK\_SUCCESS is returned and 
247      MBLk\_FAILURE is returned in case of error.
248 \vspace{0.2in}
249
250 \function{void MBLK\_Iupdate\_field(const int fid,int ghostwidth, void *ingrid, void* outgrid);}
251 \function{subroutine MBLK\_Iupdate\_field(fid,ghostwidth, ingrid, outgrid,err)}
252   \args{integer, intent(in)  :: fid, ghostwidth}
253   \args{integer,intent(out) :: err}
254   \args{varies,intent(in) :: ingrid}
255   \args{varies,intent(out) :: outgrid}
256
257
258      Update the values in the ghost regions which were specified when the
259      field was created. For the example above the ghost regions will be 
260      updated once for each step in the time loop.
261
262      \kw{MBLK\_Iupdate\_field} can only be called from driver, and to be useful,
263      must be called from every block's driver routine.
264
265      \kw{MBLK\_Iupdate\_field} is a non blocking call like MPI\_IRecv
266      After the routine returns the update is not complete but is guranteed
267      to be complete in future. So before using the values the status of the
268      update should be checked using \kw{MBLk\_Test\_update} or should wait
269      for the completion of the call using \kw{MBLK\_Wait\_update}
270      
271      If the\kw{MBLK\_Iupdate\_field} call was successful MBLK\_SUCCESS is 
272      returned and MBLK\_FAILURE is returned in case of error.
273 \vspace{0.2in}
274
275 \function{int MBLK\_Test\_update(int *status);}
276 \function{subroutine MBLK\_Test\_update(status,err)}
277 \args{integer, intent(out) :: status,err}
278
279      \kw{MBLK\_Test\_update} is a call that is used in assosiation with 
280      \kw{MBLk\_Iupdate\_field} from the driver sub routine.It tests whether
281       the field has been updated or not.
282      \kw{status} is returned as MBLK\_DONE if the update was completed or 
283       MBLK\_NOTDONE if the update is still pending.
284      \kw{err} is returned as MBLK\_SUCCESS if the call was successful or 
285      MBLK\_FAILURE if there was an error.
286
287
288 ORION: I think this function at this time is not updating the status I think it is returning MBLK\_DONE or MBLK\_NOTDONE or MBLK\_FAILURE in error instead of using the 
289 status. Please have a look at that
290 \vspace{0.2in}
291
292 \function{void MBLk\_Wait\_update(void);}
293 \function{subroutine MBLk\_Wait\_update()}
294
295      \kw{MBLk\_Wait\_update} call is a blocking call and is used in assoisation 
296      with \kw{MBLK\_Iupdate\_field} call. It blocks till the update is completed.
297 \vspace{0.2 in}
298
299
300 \function{void MBLK\_Reduce\_field(int fid,void *grid, void *out,int op);}
301 \function{subroutine MBLK\_Reduce\_field(fid,grid,outVal,op)}
302   \args{integer, intent(in)  :: fid,op}
303   \args{varies, intent(in) :: grid}
304   \args{varies, intent(out) :: outVal}
305
306      Combine a field from each block, according to op, across all blocks.
307      After \kw{Reduce\_Field} returns, all blocks will have identical
308 values in \kw{outVal,} which must be \kw{vec\_len} copies of \kw{base\_type}.
309
310      May only be called from driver, and to complete, must be called
311      from every chunk's driver routine.
312
313      \kw{op} must be one of:
314
315 \begin{itemize}
316         \item \kw{MBLk\_SUM}-- each element of \kw{outVal} will be the sum 
317 of the corresponding fields of all blocks
318         \item \kw{MBLK\_MIN}-- each element of \kw{outVal} will be the 
319 smallest value among the corresponding field of all blocks
320         \item \kw{MBLK\_MAX}-- each element of \kw{outVal} will be the largest 
321 value among the corresponding field of all blocks
322 \end{itemize}
323 \vspace{0.2in}
324
325 \function{void MBLK\_Reduce(int fid,void *inVal,void *outVal,int op);}
326 \function{subroutine MBLK\_Reduce(fid,inVal,outVal,op)}
327   \args{integer, intent(in)  :: fid,op}
328   \args{varies, intent(in) :: inVal}
329   \args{varies, intent(out) :: outVal}
330
331      Combine a field from each block, acoording to \kw{op}, across all blocks.
332 \kw{Fid} is only used for the \kw{base\_type} and \kw{vec\_len}-- offset and
333 \kw{dist} are not used.  After this call returns, all blocks will have
334 identical values in \kw{outVal}.  Op has the same values and meaning as
335 \kw{MBLK\_Reduce\_Field}.
336
337      May only be called from driver, and to complete, must be called
338      from every blocks driver routine.
339
340 \vspace{0.2in}
341
342 \subsection{Boundary Conditions}
343 In most of the applications using the MultiBlock Framework the blocks have 
344 boundary conditions for eaach block depending on the geometery. Various calls 
345 are provided in the framework to register and apply the boundry conditions.
346
347 \function{int MBLK\_Register\_bc(const int bcnum, int ghostWidth, const MBLK\_BcFn bcfn);}
348 \function{subroutine MBLK\_Register\_bc(bcnum, ghostwidth, bcfn, err)}
349 \args{integer,intent(in) :: bcnum, ghostWidth}
350 \args{integer,intent(out) :: err}
351 \args{subroutine :: bcfn}
352
353 This is call is used to register the boundry condition function, for a block,
354 with the framework.
355 \begin{itemize}
356         \item \kw{bcnum} The boundry condtion number to be associated with the 
357         function.
358         \item \kw{ghostWidth} The width of the ghostcells where this boundry
359         condition is going to be applied.
360         \item \kw{bcfn} The user function that is to be used for applying the
361         the boundry conditions.
362 \end{itemize}
363 \kw{MBLK\_Register\_bc} should only be called from the driver.
364 This call returns MBLK\_SUCCESS in case of success or else returns MBLK\_FAILURE
365 in case of an error.
366 \vspace{0.2in}
367
368 \function{ int MBLK\_Apply\_bc(const int bcnum, void *grid,void *size);}
369 \function{subroutine MBLK\_Apply\_bc(bcnum, grid, size,err)}
370 \args{ integer,intent(in)::bcnum}
371 \args{integer,intent(out)::err}
372 \args{varies,intent(inout)::grid}
373 \args{varies, intent(in)::size}
374
375 \kw{MBLK\_Apply\_bc} call is made to apply the boundry condition function
376 associated to \kw{bcnum} to the block.The grid specifies the place where 
377 the boundary condition are to be applied adn sizes array gives the dimensions
378 of the grid.
379 It returns MBLK\_SUCCESS if the call is successful else it returns
380 MBLK\_FAILURE in case of error. 
381
382 \function{ int MBLK\_Apply\_bc\_all(void* grid, void* size);} 
383 \function{subroutine MBLK\_Apply\_bc\_all(grid, size, err)}
384 \args{integer,intent(out)::err}
385 \args{varies,intent(inout)::grid}
386 \args{varies, intent(in)::size}
387 This call is same as \kw{MBLK\_Apply\_bc} except it applies all the boundary 
388 functions to the block.
389
390 \subsection{Migration}
391
392 The \charmpp\ runtime framework includes an automated, run-time load balancer,
393 which will automatically monitor the performance of your parallel program.
394 If needed, the load balancer can ``migrate'' mesh chunks from heavily-loaded
395 processors to more lightly-loaded processors, improving the load balance and
396 speeding up the program.  For this to be useful, pass the \kw{+vpN} argument
397 with a larger number of blocks \kw{N} than processors
398 Because this is somewhat involved, you may refrain from calling 
399 \kw{MBLK\_Migrate} and migration will never take place.
400
401 The runtime system can automatically move your thread stack to the new
402 processor, but you must write a PUP function to move any global or
403 heap-allocated data to the new processor (global data is declared at file scope
404 or \kw{static} in C and \kw{COMMON} in Fortran77; heap allocated data comes
405 from C \kw{malloc}, C++ \kw{new}, or Fortran90 \kw{ALLOCATE}).  A PUP
406 (Pack/UnPack) function performs both packing (converting heap data into a
407 message) and unpacking (converting a message back into heap data).  All your
408 global and heap data must be collected into a single block (\kw{struct} in C;
409 user-defined \kw{TYPE} in Fortran) so the PUP function can access it all.
410
411 Your PUP function will be passed a pointer to your heap data block and a
412 special handle called a ``pupper'', which contains the network message to be
413 sent.  Your PUP function returns a pointer to your heap data block.  In a PUP
414 function, you pass all your heap data to routines named \kw{pup\_type}, where
415 type is either a basic type (such as int, char, float, or double) or an array
416 type (as before, but with a ``s'' suffix).  Depending on the direction of
417 packing, the pupper will either read from or write to the values you pass--
418 normally, you shouldn't even know which.  The only time you need to know the
419 direction is when you are leaving a processor or just arriving.
420 Correspondingly, the pupper passed to you may be deleting (indicating that you
421 are leaving the processor, and should delete your heap storage after packing),
422 unpacking (indicating you've just arrived on a processor, and should allocate
423 your heap storage before unpacking), or neither (indicating the system is
424 merely sizing a buffer, or checkpointing your values).
425
426 PUP functions are much easier to write than explain-- a simple C heap block
427 and the corresponding PUP function is:
428
429 \begin{alltt}
430      typedef struct {
431        int n1;/*Length of first array below*/
432        int n2;/*Length of second array below*/
433        double *arr1; /*Some doubles, allocated on the heap*/
434        int *arr2; /*Some ints, allocated on the heap*/
435      } my_block;
436  
437      my_block *pup_my_block(pup_er p,my_block *m)
438      {
439        if (pup_isUnpacking(p)) m=malloc(sizeof(my_block));
440        pup_int(p,\&m->n1);
441        pup_int(p,\&m->n2);
442        if (pup_isUnpacking(p)) {
443          m->arr1=malloc(m->n1*sizeof(double));
444          m->arr2=malloc(m->n2*sizeof(int));
445        }
446        pup_doubles(p,m->arr1,m->n1);
447        pup_ints(p,m->arr2,m->n2);
448        if (pup_isDeleting(p)) {
449          free(m->arr1);
450          free(m->arr2);
451          free(m);
452        }
453        return m;
454      }
455 \end{alltt}
456
457 This single PUP function can be used to copy the \kw{my\_block} data into a
458 message buffer and free the old heap storage (deleting pupper); allocate
459 storage on the new processor and copy the message data back (unpacking pupper);
460 or save the heap data for debugging or checkpointing.
461
462 A Fortran block TYPE and corresponding PUP routine is as follows:
463
464 \begin{alltt}
465      MODULE my_block_mod
466        TYPE my_block
467          INTEGER :: n1,n2x,n2y
468          REAL*8, POINTER, DIMENSION(:) :: arr1
469          INTEGER, POINTER, DIMENSION(:,:) :: arr2
470        END TYPE
471      END MODULE
472  
473      SUBROUTINE pup_my_block(p,m)
474        IMPLICIT NONE
475        USE my_block_mod
476        USE pupmod
477        INTEGER :: p
478        TYPE(my_block) :: m
479        call pup_int(p,m%n1)
480        call pup_int(p,m%n2x)
481        call pup_int(p,m%n2y)
482        IF (pup_isUnpacking(p)) THEN
483          ALLOCATE(m%arr1(m%n1))
484          ALLOCATE(m%arr2(m%n2x,m%n2y))
485        END IF
486        call pup_doubles(p,m%arr1,m%n1)
487        call pup_ints(p,m%arr2,m%n2x*m%n2y)
488        IF (pup_isDeleting(p)) THEN
489          DEALLOCATE(m%arr1)
490          DEALLOCATE(m%arr2)
491        END IF
492      END SUBROUTINE
493 \end{alltt}
494
495 \function{int MBLK\_Register(void *block, MBLk\_PupFn pup\_ud, int* rid)}
496 \function{subroutine MBLK\_Register(block,pup\_ud, rid)}
497     \args{integer, intent(out)::rid}
498     \args{TYPE(varies), POINTER :: block}
499     \args{SUBROUTINE :: pup\_ud}
500
501      Associates the given data block and PUP function.  Returns a block
502      ID, which can be passed to \kw{MBLK\_Get\_registered} later.  Can only be
503      called from driver.  It returns MBLK\_SUCESS if the call was successful
504      and MBLK\_FAILURE in case of error. For the declarations above, you call
505      \kw{MBLK\_Register} as:
506
507 \begin{alltt}
508           /*C/C++ driver() function*/
509           int myId, err;
510           my_block *m=malloc(sizeof(my_block));
511           err =MBLK_Register(m,(MBLK_PupFn)pup_my_block,&rid);
512  
513           !- Fortran driver subroutine
514           use my_block_mod
515           interface
516             subroutine pup_my_block(p,m)
517               use my_block_mod
518               INTEGER :: p
519               TYPE(my_block) :: m
520             end subroutine
521           end interface
522           TYPE(my_block) :: m
523           INTEGER :: myId,err
524           MBLK_Register(m,pup_my_block,myId,err)
525 \end{alltt}
526
527      Note that Fortran blocks must be allocated on the stack in driver;
528      while C/C++ blocks may be allocated on the heap.
529 \vspace{0.2in}
530
531 \function{void MBLK\_Migrate()}
532 \function{subroutine MBLK\_Migrate()}
533
534      Informs the load balancing system that you are ready to be
535      migrated, if needed.  If the system decides to migrate you, the
536      PUP function passed to \kw{MBLK\_Register} will be called with a sizing
537      pupper, then a packing, deleting pupper.  Your stack (and pupped
538      data) will then be sent to the destination machine, where your PUP
539      function will be called with an unpacking pupper.  \kw{MBLK\_Migrate}
540      will then return, whereupon you should call \kw{MBLK\_Get\_registered} to
541      get your unpacked data block.  Can only be called from driver.
542
543
544
545 \function{int MBLK\_Get\_Userdata(int n, void** block)}
546
547      Return your unpacked userdata after migration-- that is, the
548      return value of the unpacking call to your PUP function.  Takes
549      the userdata ID returned by \kw{MBLK\_Register}.  Can be called from
550      driver at any time.
551
552      Since Fortran blocks are always allocated on the stack, the system
553      migrates them to the same location on the new processor, so no
554      \kw{Get\_registered} call is needed from Fortran.
555 \input{index}
556 \end{document}