2
* ***************************************************************************
3
* MALOC = < Minimal Abstraction Layer for Object-oriented C >
4
* Copyright (C) 1994--2006 Michael Holst
6
* This program is free software; you can redistribute it and/or modify it
7
* under the terms of the GNU General Public License as published by the
8
* Free Software Foundation; either version 2 of the License, or (at your
9
* option) any later version.
11
* This program is distributed in the hope that it will be useful,
12
* but WITHOUT ANY WARRANTY; without even the implied warranty of
13
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
14
* See the GNU General Public License for more details.
16
* You should have received a copy of the GNU General Public License along
17
* with this program; if not, write to the Free Software Foundation, Inc.,
18
* 675 Mass Ave, Cambridge, MA 02139, USA.
20
* rcsid="$Id: vmpi.c,v 1.20 2006/06/03 07:22:29 mholst Exp $"
21
* ***************************************************************************
25
* ***************************************************************************
28
* Purpose: Class Vmpi: methods.
30
* Notes: Class Vmpi is a thin object-oriented Clean C layer on top of the
31
* MPI communication library. Vmpi provides access to the minimal
32
* set of ten MPI primitives required to implement the Bank-Holst
33
* parallel adaptive algorithm, using either the Bank-Holst Oracle
34
* library, or directly.
36
* Vmpi is also just about the smallest possible communication
37
* library that still contains just enough functionality to solve
38
* an elliptic equation in parallel using a typical parallel
39
* algorithm such as domain decomposition or parallel multigrid.
41
* This minimal functionality is provided by Vmpi in ten core
42
* methods: startup, shutdown, size, rank, send, receive,
43
* broadcast, reduce, barrier, and non-blocking send (sometimes
44
* very useful). That is it, and it all fits in about 90 lines
45
* of object-oriented Clean C code (not counting comments...).
46
* To fit into such a small amount of code, Vmpi assumes that the
47
* user does all of his own datapacking into pure character arrays.
48
* The best way to make use of Vmpi is to pack every message into
49
* a character buffer in XDR format (RPC's answer to
50
* machine-independent binary format). You can accomplish this
51
* e.g. using the Vio library that sits at the bottom of MC.
53
* The 6 primary MPI routines (according to Ian Foster):
54
* -----------------------------------------------------
56
* MPI_Init : Initiate an MPI computation
57
* MPI_Finalize : Terminate a computation
58
* MPI_Comm_size : Determine number of processes
59
* MPI_Comm_rank : Determine my process identifier
60
* MPI_Send : Send a message (blocking)
61
* MPI_Recv : Receive a message (blocking)
63
* The 4 additional useful MPI routines (according to me):
64
* -------------------------------------------------------
66
* MPI_Barrier : Synchronizes all processes
67
* MPI_Bcast : Sends data from one process to all processes
68
* MPI_Reduce : Sums, etc., distributed data (result in root)
69
* MPI_Isend : Send a message (non-blocking)
71
* Author: Michael Holst
72
* ***************************************************************************
77
#if defined(HAVE_MPI_H)
81
VEMBED(rcsid="$Id: vmpi.c,v 1.20 2006/06/03 07:22:29 mholst Exp $")
84
* ***************************************************************************
85
* Class Vmpi: Inlineable methods
86
* ***************************************************************************
88
#if !defined(VINLINE_MALOC)
90
#endif /* if !defined(VINLINE_MALOC) */
92
* ***************************************************************************
93
* Class Vmpi: Non-inlineable methods
94
* ***************************************************************************
98
* ***************************************************************************
101
* Purpose: The Vmp initializer.
103
* Author: Michael Holst
104
* ***************************************************************************
106
VPUBLIC int Vmpi_init(int *argc, char ***argv)
108
#if defined(HAVE_MPI_H)
109
return (MPI_SUCCESS == MPI_Init(argc,argv));
116
* ***************************************************************************
117
* Routine: Vmpi_finalize
119
* Purpose: The Vmp finalizerr.
121
* Author: Michael Holst
122
* ***************************************************************************
124
VPUBLIC int Vmpi_finalize(void)
126
#if defined(HAVE_MPI_H)
127
return (MPI_SUCCESS == MPI_Finalize());
134
* ***************************************************************************
137
* Purpose: The Vmpi constructor.
139
* Author: Michael Holst
140
* ***************************************************************************
142
VPUBLIC Vmpi* Vmpi_ctor(void)
144
#if defined(HAVE_MPI_H)
148
VDEBUGIO("Vmpi_ctor: CREATING object..");
149
thee = Vmem_malloc( VNULL, 1, sizeof(Vmpi) );
152
#if defined(HAVE_MPI_H)
153
VASSERT( MPI_SUCCESS == MPI_Initialized(&dummy) );
154
VASSERT( MPI_SUCCESS == MPI_Comm_rank(MPI_COMM_WORLD, &(thee->mpi_rank)) );
155
VASSERT( MPI_SUCCESS == MPI_Comm_size(MPI_COMM_WORLD, &(thee->mpi_size)) );
156
Vnm_setIoTag(thee->mpi_rank, thee->mpi_size);
157
Vnm_print(2,"Vmpi_ctor: process %d of %d is ALIVE!\n",
158
thee->mpi_rank, thee->mpi_size);
160
VDEBUGIO("..done.\n");
165
* ***************************************************************************
168
* Purpose: The Vmpi destructor.
170
* Author: Michael Holst
171
* ***************************************************************************
173
VPUBLIC void Vmpi_dtor(Vmpi **thee)
175
VDEBUGIO("Vmpi_dtor: DESTROYING object..");
176
Vmem_free( VNULL, 1, sizeof(Vmpi), (void**)thee );
177
#if defined(HAVE_MPI_H)
179
VASSERT( MPI_SUCCESS == MPI_Finalize() );
182
VDEBUGIO("..done.\n");
186
* ***************************************************************************
189
* Purpose: Return my processor ID.
191
* Author: Michael Holst
192
* ***************************************************************************
194
VPUBLIC int Vmpi_rank(Vmpi *thee)
196
return thee->mpi_rank;
200
* ***************************************************************************
203
* Purpose: Return the number of processors involved.
205
* Author: Michael Holst
206
* ***************************************************************************
208
VPUBLIC int Vmpi_size(Vmpi *thee)
210
return thee->mpi_size;
214
* ***************************************************************************
217
* Purpose: An MPI barrier.
219
* Author: Michael Holst
220
* ***************************************************************************
222
VPUBLIC int Vmpi_barr(Vmpi *thee)
225
#if defined(HAVE_MPI_H)
226
rc = (MPI_SUCCESS == MPI_Barrier(MPI_COMM_WORLD));
232
* ***************************************************************************
235
* Purpose: An MPI blocking send.
237
* Author: Michael Holst
238
* ***************************************************************************
240
VPUBLIC int Vmpi_send(Vmpi *thee, int des, char *buf, int bufsize)
243
#if defined(HAVE_MPI_H)
245
rc = (MPI_SUCCESS == MPI_Send(buf,bufsize,MPI_CHAR,des,tag,MPI_COMM_WORLD));
251
* ***************************************************************************
254
* Purpose: An MPI blocking receive.
256
* Author: Michael Holst
257
* ***************************************************************************
259
VPUBLIC int Vmpi_recv(Vmpi *thee, int src, char *buf, int bufsize)
262
#if defined(HAVE_MPI_H)
267
rsrc = MPI_ANY_SOURCE;
271
rc = (MPI_SUCCESS == MPI_Recv(buf,bufsize,MPI_CHAR,rsrc,tag,
272
MPI_COMM_WORLD,&stat));
278
* ***************************************************************************
279
* Routine: Vmpi_bcast
281
* Purpose: An MPI broadcast.
283
* Author: Michael Holst
284
* ***************************************************************************
286
VPUBLIC int Vmpi_bcast(Vmpi *thee, char *buf, int bufsize)
289
#if defined(HAVE_MPI_H)
290
rc = (MPI_SUCCESS == MPI_Bcast(buf,bufsize,MPI_CHAR,0,MPI_COMM_WORLD));
296
* ***************************************************************************
297
* Routine: Vmpi_reduce
299
* Purpose: An MPI reduce.
301
* Author: Michael Holst
302
* ***************************************************************************
304
VPUBLIC int Vmpi_reduce(Vmpi *thee, char *sbuf, char *rbuf, int bufsize)
307
#if defined(HAVE_MPI_H)
309
MPI_Reduce(sbuf,rbuf,bufsize,MPI_CHAR,MPI_SUM,0,MPI_COMM_WORLD));
315
* ***************************************************************************
316
* Routine: Vmpi_isend
318
* Purpose: An MPI non-blocking send.
320
* Author: Michael Holst
321
* ***************************************************************************
323
VPUBLIC int Vmpi_isend(Vmpi *thee, int des, char *buf, int bufsize)
326
#if defined(HAVE_MPI_H)
329
rc = (MPI_SUCCESS == MPI_Isend(buf,bufsize,MPI_CHAR,des,tag,MPI_COMM_WORLD,