File: cartrank.c

package info (click to toggle)
lam 7.1.2-1.4
  • links: PTS
  • area: main
  • in suites: lenny
  • size: 54,608 kB
  • ctags: 17,034
  • sloc: ansic: 156,264; sh: 9,976; cpp: 7,699; makefile: 5,589; perl: 476; fortran: 260; asm: 83
file content (117 lines) | stat: -rw-r--r-- 2,491 bytes parent folder | download | duplicates (10)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
/*
 * Copyright (c) 2001-2002 The Trustees of Indiana University.  
 *                         All rights reserved.
 * Copyright (c) 1998-2001 University of Notre Dame. 
 *                         All rights reserved.
 * Copyright (c) 1994-1998 The Ohio State University.  
 *                         All rights reserved.
 * 
 * This file is part of the LAM/MPI software package.  For license
 * information, see the LICENSE file in the top level directory of the
 * LAM/MPI source distribution.
 * 
 * $HEADER$
 *
 * $Id: cartrank.c,v 6.5 2002/12/11 19:15:09 jsquyres Exp $
 *
 *	Function:	- translate coordinates to rank
 *	Accepts:	- communicator
 *			- coordinates array
 *			- ptr rank (returned value)
 *	Returns:	- MPI_SUCCESS or error code
 */

#include <blktype.h>
#include <mpi.h>
#include <mpisys.h>

/*@

MPI_Cart_rank - Determines process rank in communicator given Cartesian
                location

Input Parameters:
+ comm - communicator with cartesian structure (handle) 
- coords - integer array (of size  'ndims') specifying the cartesian 
  coordinates of a process 

Output Parameter:
. prank - rank of specified process (integer) 

.N fortran

.N Errors
.N MPI_SUCCESS
.N MPI_ERR_COMM
.N MPI_ERR_TOPOLOGY
.N MPI_ERR_ARG

.N ACK
@*/
int MPI_Cart_rank(MPI_Comm comm, int *coords, int *prank)
{
	int		rank;
	int		dim;
	int		ord;
	int		factor;
	int		i;
	int		*d;
	int		*c;

	lam_initerr();
	lam_setfunc(BLKMPICARTRANK);
/*
 * Check the arguments.
 */
	if (comm == MPI_COMM_NULL) {
		return(lam_errfunc(MPI_COMM_WORLD,
			BLKMPICARTRANK, lam_mkerr(MPI_ERR_COMM, EINVAL)));
	}

	if (LAM_IS_INTER(comm)) {
		return(lam_errfunc(comm,
			BLKMPICARTRANK, lam_mkerr(MPI_ERR_COMM, EINVAL)));
	}

	if (!LAM_IS_CART(comm)) {
		return(lam_errfunc(comm, BLKMPICARTRANK,
				lam_mkerr(MPI_ERR_TOPOLOGY, EINVAL)));
	}

	if ((coords == 0) || (prank == 0)) {
		return(lam_errfunc(comm,
			BLKMPICARTRANK, lam_mkerr(MPI_ERR_ARG, EINVAL)));
	}
/*
 * Loop over coordinates computing the rank.
 */
	factor = 1;
	rank = 0;
	i = comm->c_topo_ndims - 1;
	d = comm->c_topo_dims + i;
	c = coords + i;

	for (; i >= 0; --i, --c, --d) {
		dim = (*d > 0) ? *d : -(*d);
		ord = *c;
		if ((ord < 0) || (ord >= dim)) {
			if (*d > 0) {
				return(lam_errfunc(comm, BLKMPICARTRANK,
					lam_mkerr(MPI_ERR_ARG, EINVAL)));
			}

			ord %= dim;
			if (ord < 0) {
				ord += dim;
			}
		}

		rank += factor * ord;
		factor *= dim;
	}

	*prank = rank;

	lam_resetfunc(BLKMPICARTRANK);
	return(MPI_SUCCESS);
}