66 if (std::abs(Value) <
chopVal_)
return 0;
74 const int m = 2147483647;
80 int test = a * lo - r * hi;
91 const double Modulus = 2147483647.0;
92 const double DbleOne = 1.0;
93 const double DbleTwo = 2.0;
96 randdouble = DbleTwo * (randdouble / Modulus) - DbleOne;
115 int NumDoubleCompanions,
double ** DoubleCompanions,
116 int NumIntCompanions,
int ** IntCompanions,
117 int NumLongLongCompanions,
long long ** LongLongCompanions)
122 T *
const list = Keys;
127 for (
int j=0; j<max; j++)
129 for (
int k=j; k>=0; k-=m)
131 if ((SortAscending && list[k+m] >= list[k]) ||
132 ( !SortAscending && list[k+m] <= list[k]))
137 for (i=0; i<NumDoubleCompanions; i++) {
138 double dtemp = DoubleCompanions[i][k+m];
139 DoubleCompanions[i][k+m] = DoubleCompanions[i][k];
140 DoubleCompanions[i][k] = dtemp;
142 for (i=0; i<NumIntCompanions; i++) {
143 int itemp = IntCompanions[i][k+m];
144 IntCompanions[i][k+m] = IntCompanions[i][k];
145 IntCompanions[i][k] = itemp;
147 for (i=0; i<NumLongLongCompanions; i++) {
148 long long LLtemp = LongLongCompanions[i][k+m];
149 LongLongCompanions[i][k+m] = LongLongCompanions[i][k];
150 LongLongCompanions[i][k] = LLtemp;
159 int NumDoubleCompanions,
double ** DoubleCompanions,
160 int NumIntCompanions,
int ** IntCompanions,
161 int NumLongLongCompanions,
long long ** LongLongCompanions)
163 Sort<int>(SortAscending, NumKeys, Keys,
164 NumDoubleCompanions, DoubleCompanions,
165 NumIntCompanions, IntCompanions,
166 NumLongLongCompanions, LongLongCompanions);
170 int NumDoubleCompanions,
double ** DoubleCompanions,
171 int NumIntCompanions,
int ** IntCompanions,
172 int NumLongLongCompanions,
long long ** LongLongCompanions)
174 Sort<long long>(SortAscending, NumKeys, Keys,
175 NumDoubleCompanions, DoubleCompanions,
176 NumIntCompanions, IntCompanions,
177 NumLongLongCompanions, LongLongCompanions);
181 int NumDoubleCompanions,
double ** DoubleCompanions,
182 int NumIntCompanions,
int ** IntCompanions,
183 int NumLongLongCompanions,
long long ** LongLongCompanions)
185 Sort<double>(SortAscending, NumKeys, Keys,
186 NumDoubleCompanions, DoubleCompanions,
187 NumIntCompanions, IntCompanions,
188 NumLongLongCompanions, LongLongCompanions);
193 int NumDoubleCompanions,
double ** DoubleCompanions,
194 int NumIntCompanions,
int ** IntCompanions)
199 int *
const list = Keys;
204 for (
int j=0; j<max; j++)
206 for (
int k=j; k>=0; k-=m)
208 if ((SortAscending && list[k+m] >= list[k]) ||
209 ( !SortAscending && list[k+m] <= list[k]))
211 int temp = list[k+m];
214 for (i=0; i<NumDoubleCompanions; i++) {
215 double dtemp = DoubleCompanions[i][k+m];
216 DoubleCompanions[i][k+m] = DoubleCompanions[i][k];
217 DoubleCompanions[i][k] = dtemp;
219 for (i=0; i<NumIntCompanions; i++) {
220 int itemp = IntCompanions[i][k+m];
221 IntCompanions[i][k+m] = IntCompanions[i][k];
222 IntCompanions[i][k] = itemp;
231 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES // FIXME 235 bool high_rank_proc_owns_shared)
249 int* owner_procs =
new int[numMyElems];
252 0, 0, high_rank_proc_owns_shared);
256 int* myOwnedElems =
new int[numMyElems];
257 int numMyOwnedElems = 0;
259 for(
int i=0; i<numMyElems; ++i) {
260 int GID = myElems[i];
261 int owner = owner_procs[i];
263 if (myPID == owner) {
264 myOwnedElems[numMyOwnedElems++] = GID;
268 Epetra_Map one_to_one_map(-1, numMyOwnedElems, myOwnedElems,
271 delete [] myOwnedElems;
272 delete [] owner_procs;
275 return(one_to_one_map);
277 #endif // EPETRA_NO_32BIT_GLOBAL_INDICES 280 template<
typename int_type>
291 bool isRoot = usermap.
Comm().
MyPID()==root;
295 int globalquickreturn = 0;
303 usermap.
Comm().
MinAll(&quickreturn, &globalquickreturn, 1);
305 if (globalquickreturn==1) {
312 int numMyElements = 0;
313 if(usermap.
MaxAllGID64()+1 > std::numeric_limits<int>::max())
314 throw "Epetra_Util::Create_Root_Map: cannot fit all gids in int";
315 if (isRoot) numMyElements = (int)(usermap.
MaxAllGID64()+1);
321 throw usermap.
ReportError(
"usermap must have unique GIDs",-1);
327 Epetra_Map allGidsMap((int_type) -1, numMyElements, (int_type) 0, comm);
329 for (
int i=0; i<numMyElements; i++) allGids[i] = (int_type) usermap.
GID64(i);
331 if(usermap.
MaxAllGID64() > std::numeric_limits<int>::max())
332 throw "Epetra_Util::Create_Root_Map: cannot fit all gids in int";
335 int n1 = 0;
if (isRoot) n1 = numGlobalElements;
336 Epetra_Map allGidsOnRootMap((int_type) -1, n1, (int_type) 0, comm);
339 allGidsOnRoot.Import(allGids, importer,
Insert);
341 Epetra_Map rootMap((int_type)-1, allGidsOnRoot.MyLength(), allGidsOnRoot.Values(), (int_type)usermap.
IndexBase64(), comm);
345 int n1 = numGlobalElements;
349 allGidsOnRoot.Import(allGids, importer,
Insert);
351 Epetra_Map rootMap((int_type) -1, allGidsOnRoot.MyLength(), allGidsOnRoot.Values(), (int_type)usermap.
IndexBase64(), comm);
360 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 362 return TCreate_Root_Map<int>(usermap, root);
366 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES 368 return TCreate_Root_Map<long long>(usermap, root);
372 throw "Epetra_Util::Create_Root_Map: GlobalIndices type unknown";
376 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES // FIXME 380 bool high_rank_proc_owns_shared)
396 int* owner_procs =
new int[numMyElems*2];
397 int* sizes = owner_procs+numMyElems;
400 0, sizes, high_rank_proc_owns_shared);
404 int* myOwnedElems =
new int[numMyElems*2];
405 int* ownedSizes = myOwnedElems+numMyElems;
406 int numMyOwnedElems = 0;
408 for(
int i=0; i<numMyElems; ++i) {
409 int GID = myElems[i];
410 int owner = owner_procs[i];
412 if (myPID == owner) {
413 ownedSizes[numMyOwnedElems] = sizes[i];
414 myOwnedElems[numMyOwnedElems++] = GID;
421 delete [] myOwnedElems;
422 delete [] owner_procs;
425 return(one_to_one_map);
427 #endif // EPETRA_NO_32BIT_GLOBAL_INDICES 435 int nnz = CRS_rowptr[NumRows];
437 for(
int i = 0; i < NumRows; i++){
438 int start=CRS_rowptr[i];
439 if(start >= nnz)
continue;
441 double* locValues = &CRS_vals[start];
442 int NumEntries = CRS_rowptr[i+1] - start;
443 int* locIndices = &CRS_colind[start];
450 for(
int j = 0; j < max; j++) {
451 for(
int k = j; k >= 0; k-=m) {
452 if(locIndices[k+m] >= locIndices[k])
454 double dtemp = locValues[k+m];
455 locValues[k+m] = locValues[k];
456 locValues[k] = dtemp;
457 int itemp = locIndices[k+m];
458 locIndices[k+m] = locIndices[k];
459 locIndices[k] = itemp;
474 int nnz = CRS_rowptr[NumRows];
475 int new_curr=CRS_rowptr[0], old_curr=CRS_rowptr[0];
477 for(
int i = 0; i < NumRows; i++){
478 int start=CRS_rowptr[i];
479 if(start >= nnz)
continue;
481 double* locValues = &CRS_vals[start];
482 int NumEntries = CRS_rowptr[i+1] - start;
483 int* locIndices = &CRS_colind[start];
491 for(
int j = 0; j < max; j++) {
492 for(
int k = j; k >= 0; k-=m) {
493 if(locIndices[k+m] >= locIndices[k])
495 double dtemp = locValues[k+m];
496 locValues[k+m] = locValues[k];
497 locValues[k] = dtemp;
498 int itemp = locIndices[k+m];
499 locIndices[k+m] = locIndices[k];
500 locIndices[k] = itemp;
507 for(
int j=CRS_rowptr[i]; j < CRS_rowptr[i+1]; j++) {
508 if(j > CRS_rowptr[i] && CRS_colind[j]==CRS_colind[new_curr-1]) {
509 CRS_vals[new_curr-1] += CRS_vals[j];
511 else if(new_curr==j) {
515 CRS_colind[new_curr] = CRS_colind[j];
516 CRS_vals[new_curr] = CRS_vals[j];
521 CRS_rowptr[i] = old_curr;
525 CRS_rowptr[NumRows] = new_curr;
531 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 544 const int *RemoteLIDs = Importer.
RemoteLIDs();
555 if(use_minus_one_for_local)
556 for(i=0;i <N; i++) gpids[i]=std::make_pair(-1,Importer.
TargetMap().
GID(i));
558 for(i=0;i <N; i++) gpids[i]=std::make_pair(mypid,Importer.
TargetMap().
GID(i));
562 for(i=0,j=0;i<NumReceives;i++){
563 int pid=ProcsFrom[i];
564 for(k=0;k<LengthsFrom[i];k++){
565 if(pid!=mypid) gpids[RemoteLIDs[j]].first=pid;
581 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES 594 const int *RemoteLIDs = Importer.
RemoteLIDs();
605 if(use_minus_one_for_local)
606 for(i=0;i <N; i++) gpids[i]=std::make_pair(-1,Importer.
TargetMap().
GID64(i));
608 for(i=0;i <N; i++) gpids[i]=std::make_pair(mypid,Importer.
TargetMap().
GID64(i));
612 for(i=0,j=0;i<NumReceives;i++){
613 int pid=ProcsFrom[i];
614 for(k=0;k<LengthsFrom[i];k++){
615 if(pid!=mypid) gpids[RemoteLIDs[j]].first=pid;
642 const int *RemoteLIDs = Importer.
RemoteLIDs();
653 if(use_minus_one_for_local)
654 for(i=0; i<N; i++) pids[i]=-1;
656 for(i=0; i<N; i++) pids[i]=mypid;
660 for(i=0,j=0;i<NumReceives;i++){
661 int pid=ProcsFrom[i];
662 for(k=0;k<LengthsFrom[i];k++){
663 if(pid!=mypid) pids[RemoteLIDs[j]]=pid;
682 RemotePIDs.resize(0);
698 for(i=0,j=0;i<NumReceives;i++){
699 int pid=ProcsFrom[i];
700 for(k=0;k<LengthsFrom[i];k++){
707 RemotePIDs.resize(0);
724 unsigned start = 0, end = len - 1;
726 while(end - start > 1) {
727 unsigned mid = (start + end) >> 1;
728 if (list[mid] < item) start = mid;
732 if (list[start] == item)
return(start);
733 if (list[end] == item)
return(end);
735 if (list[end] < item) {
740 if (list[start] < item) insertPoint = end;
741 else insertPoint = start;
751 return Epetra_Util_binary_search<int>(item, list, len, insertPoint);
756 const long long* list,
760 return Epetra_Util_binary_search<long long>(item, list, len, insertPoint);
776 unsigned start = 0, end = len - 1;
778 while(end - start > 1) {
779 unsigned mid = (start + end) >> 1;
780 if (aux_list[list[mid]] < item) start = mid;
784 if (aux_list[list[start]] == item)
return(start);
785 if (aux_list[list[end]] == item)
return(end);
787 if (aux_list[list[end]] < item) {
792 if (aux_list[list[start]] < item) insertPoint = end;
793 else insertPoint = start;
805 return Epetra_Util_binary_search_aux<int>(item, list, aux_list, len, insertPoint);
811 const long long* aux_list,
815 return Epetra_Util_binary_search_aux<long long>(item, list, aux_list, len, insertPoint);
823 int & M,
int & N,
int & nz,
int * & ptr,
824 int * & ind,
double * & val,
int & Nrhs,
825 double * & rhs,
int & ldrhs,
826 double * & lhs,
int & ldlhs) {
866 for (
int i=0; i<M; i++) ptr[i+1] = ptr[i] + Graph.
NumMyIndices(i);
874 template void Epetra_Util::Sort<int> (bool, int,
int *, int,
double **, int,
int **, int,
long long **);
875 template void Epetra_Util::Sort<long long>(bool, int,
long long *, int,
double **, int,
int **, int,
long long **);
Epetra_MultiVector: A class for constructing and using dense multi-vectors, vectors and matrices in p...
int MyGlobalElements(int *MyGlobalElementList) const
Puts list of global elements on this processor into the user-provided array.
Epetra_Map: A class for partitioning vectors and matrices.
bool UniqueGIDs() const
Returns true if map GIDs are 1-to-1.
bool ConstantStride() const
Returns true if this multi-vector has constant stride between vectors.
const Epetra_CrsGraph & Graph() const
Returns a reference to the Epetra_CrsGraph object associated with this matrix.
static Epetra_Map Create_OneToOne_Map(const Epetra_Map &usermap, bool high_rank_proc_owns_shared=false)
Epetra_Util Create_OneToOne_Map function.
int NumMyRows() const
Returns the number of matrix rows owned by the calling processor.
int NumMyCols() const
Returns the number of entries in the set of column-indices that appear on this processor.
virtual Epetra_Directory * CreateDirectory(const Epetra_BlockMap &Map) const =0
Create a directory object for the given Epetra_BlockMap.
int NumRemoteIDs() const
Returns the number of elements that are not on the calling processor.
static const double chopVal_
#define EPETRA_CHK_ERR(a)
int Stride() const
Returns the stride between vectors in the multi-vector (only meaningful if ConstantStride() is true)...
int IndexBase() const
Index base for this map.
const Epetra_BlockMap & TargetMap() const
Returns the TargetMap used to construct this importer.
unsigned int Seed() const
Get seed from Random function.
MPI implementation of Epetra_Distributor.
virtual int ReportError(const std::string Message, int ErrorCode) const
Error reporting method.
virtual int MinAll(double *PartialMins, double *GlobalMins, int Count) const =0
Epetra_Comm Global Min function.
long long NumGlobalElements64() const
virtual int GetDirectoryEntries(const Epetra_BlockMap &Map, const int NumEntries, const int *GlobalEntries, int *Procs, int *LocalEntries, int *EntrySizes, bool high_rank_sharing_procs=false) const =0
GetDirectoryEntries : Returns proc and local id info for non-local map entries.
Epetra_Import: This class builds an import object for efficient importing of off-processor elements...
int * RemoteLIDs() const
List of elements in the target map that are coming from other processors.
static int SortAndMergeCrsEntries(int NumRows, int *CRS_rowptr, int *CRS_colind, double *CRS_vals)
Epetra_Util SortAndMergeCrsEntries function.
static int GetPids(const Epetra_Import &Importer, std::vector< int > &pids, bool use_minus_one_for_local)
Epetra_Util GetPids function.
virtual int MyPID() const =0
Return my process ID.
int MakeDataContiguous()
Eliminates memory that is used for construction. Make consecutive row index sections contiguous...
unsigned int RandomInt()
Returns a random integer on the interval (0, 2^31-1)
long long MaxAllGID64() const
bool IndicesAreContiguous() const
If matrix indices are packed into single array (done in OptimizeStorage()) return true...
Epetra_Directory: This class is a pure virtual class whose interface allows Epetra_Map and Epetr_Bloc...
static int GetPidGidPairs(const Epetra_Import &Importer, std::vector< std::pair< int, int > > &gpids, bool use_minus_one_for_local)
Epetra_Util GetPidGidPairs function.
static int SortCrsEntries(int NumRows, const int *CRS_rowptr, int *CRS_colind, double *CRS_vals)
Epetra_Util SortCrsEntries function.
int SetSeed(unsigned int Seed_in)
Set seed for Random function.
static Epetra_BlockMap Create_OneToOne_BlockMap(const Epetra_BlockMap &usermap, bool high_rank_proc_owns_shared=false)
Epetra_Util Create_OneToOne_Map function.
Epetra_Comm: The Epetra Communication Abstract Base Class.
long long GID64(int LID) const
static Epetra_Map Create_Root_Map(const Epetra_Map &usermap, int root=0)
Epetra_Util Create_Root_Map function.
int NumVectors() const
Returns the number of vectors in the multi-vector.
bool GlobalIndicesInt() const
Returns true if map create with int NumGlobalElements.
int NumMyElements() const
Number of elements on the calling processor.
static Epetra_Map TCreate_Root_Map(const Epetra_Map &usermap, int root)
Epetra_GIDTypeVector: A class for constructing and using dense "int" and "long long" vectors on a par...
Epetra_BlockMap: A class for partitioning block element vectors and matrices.
const Epetra_Comm & Comm() const
Access function for Epetra_Comm communicator.
int Epetra_Util_binary_search_aux(T item, const int *list, const T *aux_list, int len, int &insertPoint)
Utility function to perform a binary-search on a list of data.
Epetra_Distributor & Distributor() const
int Epetra_Util_binary_search(T item, const T *list, int len, int &insertPoint)
Utility function to perform a binary-search on a list of data.
int NumMyNonzeros() const
Returns the number of nonzero entries in the calling processor's portion of the matrix.
long long IndexBase64() const
static int GetRemotePIDs(const Epetra_Import &Importer, std::vector< int > &RemotePIDs)
Epetra_Util GetRemotePIDs.
virtual int NumProc() const =0
Returns total number of processes.
const int * LengthsFrom() const
Number of values we're receiving from each proc.
double RandomDouble()
Returns a random double on the interval (-1.0,1.0)
Epetra_CrsMatrix: A class for constructing and using real-valued double-precision sparse compressed r...
int NumMyIndices(int Row) const
Returns the current number of nonzero entries in specified local row on this processor.
static double Chop(const double &Value)
Epetra_Util Chop method. Return zero if input Value is less than ChopValue.
static void Sort(bool SortAscending, int NumKeys, T *Keys, int NumDoubleCompanions, double **DoubleCompanions, int NumIntCompanions, int **IntCompanions, int NumLongLongCompanions, long long **LongLongCompanions)
Epetra_Util Sort Routine (Shell sort)
int NumReceives() const
The number of procs from which we will receive data.
bool LinearMap() const
Returns true if the global ID space is contiguously divided (but not necessarily uniformly) across al...
const int * ProcsFrom() const
A list of procs sending values to this proc.
bool GlobalIndicesLongLong() const
Returns true if map create with long long NumGlobalElements.
int GID(int LID) const
Returns global ID of local ID, return IndexBase-1 if not found on this processor. ...
Epetra_CrsGraph: A class for constructing and using sparse compressed row graphs. ...
Epetra_LocalMap: A class for replicating vectors and matrices across multiple processors.
int Epetra_Util_ExtractHbData(Epetra_CrsMatrix *A, Epetra_MultiVector *LHS, Epetra_MultiVector *RHS, int &M, int &N, int &nz, int *&ptr, int *&ind, double *&val, int &Nrhs, double *&rhs, int &ldrhs, double *&lhs, int &ldlhs)
Harwell-Boeing data extraction routine.