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