FEI Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
fei_EqnComm.cpp
Go to the documentation of this file.
1/*--------------------------------------------------------------------*/
2/* Copyright 2007 Sandia Corporation. */
3/* Under the terms of Contract DE-AC04-94AL85000, there is a */
4/* non-exclusive license for use of this work by or on behalf */
5/* of the U.S. Government. Export of this program may require */
6/* a license from the United States Government. */
7/*--------------------------------------------------------------------*/
8
9#include "fei_EqnComm.hpp"
10#include "fei_sstream.hpp"
11
12namespace fei {
13
14EqnComm::EqnComm(MPI_Comm comm, int numLocalEqns)
15 : comm_(comm),
16 globalOffsets_(2)
17{
18#ifndef FEI_SER
19
21 MPI_Comm_size(comm, &numProcs);
22 MPI_Comm_rank(comm, &localProc);
23
24 std::vector<int> local(numProcs*2, 0);
25 int* global = &local[0] + numProcs;
26
27 if (numLocalEqns < 0) {
28 throw std::runtime_error("fei::EqnComm ERROR, negative numLocalEqns not allowed.");
29 }
30
31 local[localProc] = numLocalEqns;
32
33 MPI_Allreduce(&local[0], global, numProcs, MPI_INT, MPI_MAX, comm_);
34
35 globalOffsets_.resize(numProcs+1);
36
37 int offset = 0;
38 for(int i=0; i<numProcs; ++i) {
39 globalOffsets_[i] = offset;
40 offset += global[i];
41 }
42 globalOffsets_[numProcs] = offset;
43
44#else
45 globalOffsets_[0] = 0;
46 globalOffsets_[1] = numLocalEqns;
47#endif
48}
49
50EqnComm::EqnComm(MPI_Comm comm, int numLocalEqns, const std::vector<int>& globalOffsets)
51 : comm_(comm),
52 globalOffsets_(globalOffsets)
53{
54}
55
57{
58}
59
60const std::vector<int>&
62{
63 return(globalOffsets_);
64}
65
66int
68{
69// std::vector<int>::const_iterator
70// iter = std::lower_bound(globalOffsets_.begin(), globalOffsets_.end(),
71// eqn);
72// int proc = iter - globalOffsets_.begin() - 1;
73// if (*iter==eqn) ++proc;
74 int proc = 0;
75 for(unsigned p=1; p<globalOffsets_.size(); ++p) {
76 if (eqn < globalOffsets_[p]) {
77 proc = p-1;
78 break;
79 }
80 }
81
82#ifndef NDEBUG
83 int numProcs = globalOffsets_.size()-1;
84 if (proc >= numProcs) {
86 osstr << "fei::EqnComm::getOwnerProc: input eqn="<<eqn<<", proc="<<proc
87 << ", ERROR, proc should be in [0.."<<numProcs-1<<"].";
88 throw std::runtime_error(std::string(osstr.str().c_str()));
89 }
90#endif
91
92 return((int)proc);
93}
94
95}//namespace fei
96
const std::vector< int > & getGlobalOffsets() const
Definition: fei_EqnComm.cpp:61
std::vector< int > globalOffsets_
Definition: fei_EqnComm.hpp:32
MPI_Comm comm_
Definition: fei_EqnComm.hpp:31
EqnComm(MPI_Comm comm, int numLocalEqns)
Definition: fei_EqnComm.cpp:14
virtual ~EqnComm()
Definition: fei_EqnComm.cpp:56
int getOwnerProc(int eqn) const
Definition: fei_EqnComm.cpp:67
#define MPI_Comm
Definition: fei_mpi.h:56
#define FEI_OSTRINGSTREAM
Definition: fei_sstream.hpp:32
int localProc(MPI_Comm comm)
int numProcs(MPI_Comm comm)