ReplicatableVector.cpp

00001 /*
00002 
00003 Copyright (C) University of Oxford, 2005-2011
00004 
00005 University of Oxford means the Chancellor, Masters and Scholars of the
00006 University of Oxford, having an administrative office at Wellington
00007 Square, Oxford OX1 2JD, UK.
00008 
00009 This file is part of Chaste.
00010 
00011 Chaste is free software: you can redistribute it and/or modify it
00012 under the terms of the GNU Lesser General Public License as published
00013 by the Free Software Foundation, either version 2.1 of the License, or
00014 (at your option) any later version.
00015 
00016 Chaste is distributed in the hope that it will be useful, but WITHOUT
00017 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018 FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
00019 License for more details. The offer of Chaste under the terms of the
00020 License is subject to the License being interpreted in accordance with
00021 English Law and subject to any action against the University of Oxford
00022 being under the jurisdiction of the English Courts.
00023 
00024 You should have received a copy of the GNU Lesser General Public License
00025 along with Chaste. If not, see <http://www.gnu.org/licenses/>.
00026 
00027 */
00028 
00029 #include "ReplicatableVector.hpp"
00030 #include "Exception.hpp"
00031 #include "PetscTools.hpp"
00032 
00033 #include <cassert>
00034 #include <iostream>
00035 
00036 // Private methods
00037 
00038 void ReplicatableVector::RemovePetscContext()
00039 {
00040     if (mToAll != NULL)
00041     {
00042         VecScatterDestroy(mToAll);
00043         mToAll = NULL;
00044     }
00045 
00046     if (mReplicated != NULL)
00047     {
00048         VecDestroy(mReplicated);
00049         mReplicated = NULL;
00050     }
00051 
00052     if (mpData != NULL)
00053     {
00054         delete[] mpData;
00055         mpData = NULL;
00056     }
00057 }
00058 
00059 // Constructors & destructors
00060 
00061 ReplicatableVector::ReplicatableVector()
00062     : mpData(NULL),
00063       mSize(0),
00064       mToAll(NULL),
00065       mReplicated(NULL)
00066 {
00067 }
00068 
00069 ReplicatableVector::ReplicatableVector(Vec vec)
00070     : mpData(NULL),
00071       mSize(0),
00072       mToAll(NULL),
00073       mReplicated(NULL)
00074 {
00075     ReplicatePetscVector(vec);
00076 }
00077 
00078 ReplicatableVector::ReplicatableVector(unsigned size)
00079     : mpData(NULL),
00080       mSize(0),
00081       mToAll(NULL),
00082       mReplicated(NULL)
00083 {
00084     Resize(size);
00085 }
00086 
00087 ReplicatableVector::~ReplicatableVector()
00088 {
00089     RemovePetscContext();
00090 }
00091 
00092 
00093 // Vector interface methods
00094 
00095 unsigned ReplicatableVector::GetSize()
00096 {
00097     return mSize;
00098 }
00099 
00100 void ReplicatableVector::Resize(unsigned size)
00101 {
00102     // PETSc stuff will be out of date
00103     RemovePetscContext();
00104 
00105     mSize = size;
00106 
00107     try
00108     {
00109         mpData = new double[mSize];
00110     }
00111     catch(std::bad_alloc &badAlloc)
00112     {
00113 #define COVERAGE_IGNORE
00114         std::cout << "Failed to allocate a ReplicatableVector of size " << size  << std::endl;
00115         PetscTools::ReplicateException(true);
00116         throw badAlloc;
00117 #undef COVERAGE_IGNORE
00118     }
00119     PetscTools::ReplicateException(false);
00120 }
00121 
00122 double& ReplicatableVector::operator[](unsigned index)
00123 {
00124     assert(index < mSize);
00125     return mpData[index];
00126 }
00127 
00128 
00129 // The workhorse methods
00130 
00131 void ReplicatableVector::Replicate(unsigned lo, unsigned hi)
00132 {
00133     // Create a PetSC vector with the array containing the distributed data
00134     Vec distributed_vec;
00135     VecCreateMPIWithArray(PETSC_COMM_WORLD, hi-lo, this->GetSize(), &mpData[lo], &distributed_vec);
00136 
00137     // Now do the real replication
00138     ReplicatePetscVector(distributed_vec);
00139 
00140     // Clean up
00141     VecDestroy(distributed_vec);
00142 }
00143 
00144 void ReplicatableVector::ReplicatePetscVector(Vec vec)
00145 {
00146     // If the size has changed then we'll need to make a new context
00147     PetscInt isize;
00148     VecGetSize(vec, &isize);
00149     unsigned size=isize;
00150 
00151     if (this->GetSize() != size)
00152     {
00153         Resize(size);
00154     }
00155     if (mReplicated == NULL)
00156     {
00157         // This creates mToAll (the scatter context) and mReplicated (to store values)
00158         VecScatterCreateToAll(vec, &mToAll, &mReplicated);
00159     }
00160 
00161     // Replicate the data
00162 //PETSc-3.x.x or PETSc-2.3.3
00163 #if ( (PETSC_VERSION_MAJOR == 3) || (PETSC_VERSION_MAJOR == 2 && PETSC_VERSION_MINOR == 3 && PETSC_VERSION_SUBMINOR == 3)) //2.3.3 or 3.x.x
00164     VecScatterBegin(mToAll, vec, mReplicated, INSERT_VALUES, SCATTER_FORWARD);
00165     VecScatterEnd  (mToAll, vec, mReplicated, INSERT_VALUES, SCATTER_FORWARD);
00166 #else
00167 //PETSc-2.3.2 or previous
00168     VecScatterBegin(vec, mReplicated, INSERT_VALUES, SCATTER_FORWARD, mToAll);
00169     VecScatterEnd  (vec, mReplicated, INSERT_VALUES, SCATTER_FORWARD, mToAll);
00170 #endif
00171 
00172     // Information is now in mReplicated PETSc vector
00173     // Copy into mData
00174     double* p_replicated;
00175     VecGetArray(mReplicated, &p_replicated);
00176     for (unsigned i=0; i<size; i++)
00177     {
00178         mpData[i] = p_replicated[i];
00179     }
00180 }

Generated on Tue May 31 14:31:41 2011 for Chaste by  doxygen 1.5.5