FEI Version of the Day
HexBeamCR.cpp
00001 /*
00002 // @HEADER
00003 // ************************************************************************
00004 //             FEI: Finite Element Interface to Linear Solvers
00005 //                  Copyright (2005) Sandia Corporation.
00006 //
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation, the
00008 // U.S. Government retains certain rights in this software.
00009 //
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Alan Williams (william@sandia.gov) 
00038 //
00039 // ************************************************************************
00040 // @HEADER
00041 */
00042 
00043 
00044 #include <fei_macros.hpp>
00045 
00046 #include <test_utils/HexBeamCR.hpp>
00047 
00048 HexBeamCR::HexBeamCR(int W, int D, int DofPerNode,
00049      int decomp, int numProcs, int localProc)
00050   : HexBeam(W, D, DofPerNode, decomp, numProcs, localProc)
00051 {
00052   totalNumElems_ = W*W*D;
00053   totalNumNodes_ = (W+1)*(W+1)*(D+1) + (W+1)*(W+1)*(2*numProcs_-1);
00054 
00055   numGlobalDOF_ = totalNumNodes_*dofPerNode_;
00056 
00057   numLocalSlices_ = D/numProcs;
00058   int remainder = D%numProcs;
00059   firstLocalSlice_ = localProc_*numLocalSlices_;
00060 
00061   switch(decomp) {
00062   case HexBeamCR::OneD:
00063     if (D < numProcs) {
00064       fei::console_out() << "HexBeamCR: too many processors." << FEI_ENDL;
00065       inErrorState_ = true;
00066       break;
00067     }
00068     if (localProc < remainder) {
00069       ++numLocalSlices_;
00070       ++firstLocalSlice_;
00071     }
00072 
00073     localNumNodes_ = numNodesPerSlice_*(numLocalSlices_+2);
00074 
00075     localNumElems_ = numElemsPerSlice_*numLocalSlices_;
00076     numLocalDOF_ = localNumNodes_*dofPerNode_;
00077 
00078     if (localProc > 0) {
00079       firstLocalElem_ = localProc*numLocalSlices_*numElemsPerSlice_;
00080 
00081       firstLocalNode_ = localProc*localNumNodes_;
00082 
00083       if (remainder <= localProc && remainder > 0) {
00084   firstLocalElem_ += remainder*numElemsPerSlice_;
00085   firstLocalNode_ += remainder*numNodesPerSlice_;
00086       }
00087     }
00088 
00089     break;
00090 
00091   case HexBeamCR::TwoD:
00092   case HexBeamCR::ThreeD:
00093   default:
00094     fei::console_out() << "HexBeamCR: invalid decomp option: " << decomp
00095    <<" aborting." << FEI_ENDL;
00096     std::abort();
00097   }
00098 
00099   localCRslice_ = firstLocalSlice_ + numLocalSlices_/2;
00100   numLocalCRs_ = numNodesPerSlice_;
00101   if (localProc_ < numProcs_-1) {
00102     numLocalCRs_ += numNodesPerSlice_;
00103   }
00104 
00105   numNodesPerCR_ = 2;
00106 }
00107 
00108 HexBeamCR::~HexBeamCR()
00109 {
00110 }
00111 
00112 int HexBeamCR::getElemConnectivity(int elemID, int* nodeIDs)
00113 {
00114   if (elemID < firstLocalElem_ || elemID > firstLocalElem_+localNumElems_) {
00115     return(-1);
00116   }
00117 
00118   int whichGlobalSlice = elemID/numElemsPerSlice_;
00119   int elemX = elemID%W_;
00120   int elemY = (elemID%(W_*W_))/W_;
00121   //FEI_COUT << "whichGlobalSlice: " << whichGlobalSlice << FEI_ENDL;
00122   int firstElemNode = (whichGlobalSlice + localProc_*2)*numNodesPerSlice_
00123                      + elemY*(W_+1) + elemX;
00124 
00125   if (whichGlobalSlice >= localCRslice_) {
00126     firstElemNode += numNodesPerSlice_;
00127   }
00128 
00129   nodeIDs[0] = firstElemNode;
00130   nodeIDs[1] = firstElemNode+1;
00131   nodeIDs[2] = firstElemNode+W_+1;
00132   nodeIDs[3] = nodeIDs[2]+1;
00133 
00134   nodeIDs[4] = nodeIDs[0]+numNodesPerSlice_;
00135   nodeIDs[5] = nodeIDs[1]+numNodesPerSlice_;
00136   nodeIDs[6] = nodeIDs[2]+numNodesPerSlice_;
00137   nodeIDs[7] = nodeIDs[3]+numNodesPerSlice_;
00138 
00139   return(0);
00140 }
00141 
00142 int HexBeamCR::getCRNodes(int** nodeIDs)
00143 {
00144   int offset = 0;
00145   int firstCRnode = (localCRslice_+localProc_*2)*numNodesPerSlice_;
00146   int i;
00147   for(i=0; i<numNodesPerSlice_; ++i) {
00148     nodeIDs[offset][0] = firstCRnode+i;
00149     nodeIDs[offset++][1] = firstCRnode+i+numNodesPerSlice_;
00150   }
00151 
00152   if (localProc_ >= numProcs_-1) return(0);
00153 
00154   int nextCRnode = firstLocalNode_ + localNumNodes_ - numNodesPerSlice_;
00155   for(i=0; i<numNodesPerSlice_; ++i) {
00156     nodeIDs[offset][0] = nextCRnode+i;
00157     nodeIDs[offset++][1] = nextCRnode+i+numNodesPerSlice_;    
00158   }
00159 
00160   return(0);
00161 }
00162 
00163 int HexBeamCR::getElemStiffnessMatrix(int elemID, double* elemMat)
00164 {
00165   if (elemID < firstLocalElem_ || elemID > firstLocalElem_+localNumElems_) {
00166     return(-1);
00167   }
00168 
00169   int i, len = nodesPerElem_*dofPerNode_*nodesPerElem_*dofPerNode_;
00170 
00171   for(i=0; i<len; ++i) {
00172     elemMat[i] = 0.0;
00173   }
00174 
00175   //Should set up some semi-realistic stiffness-matrix coefficients here...
00176   //For now just use arbitrary numbers and set it up so the matrix won't be
00177   //too ill-conditioned. (This is intended for an assembly test more than
00178   //a solver test.)
00179 
00180   //Now set the diagonal to 4.0
00181   len = nodesPerElem_*dofPerNode_;
00182   for(i=0; i<len; ++i) {
00183     int offset = i*len+i;
00184     elemMat[offset] = 4.0;
00185   }
00186 
00187   //Now set some off-diagonals
00188   for(i=0; i<len; ++i) {
00189     int offset = i*len+i;
00190     if (i>1) {
00191       elemMat[offset-2] = -0.5;
00192     }
00193 
00194     if (i<len-2) {
00195       elemMat[offset+2] = -0.5;
00196     }
00197 
00198     if (i>3) {
00199       elemMat[offset-4] = -0.1;
00200     }
00201     if (i<len-4) {
00202       elemMat[offset+4] = -0.1;
00203     }
00204   }
00205 
00206   return(0);
00207 }
00208 
00209 int HexBeamCR::getElemLoadVector(int elemID, double* elemVec)
00210 {
00211   if (elemID < firstLocalElem_ || elemID > firstLocalElem_+localNumElems_) {
00212     return(-1);
00213   }
00214 
00215   int i, len = nodesPerElem_*dofPerNode_;
00216   for(i=0; i<len; ++i) {
00217     elemVec[i] = 1.0;
00218   }
00219 
00220   return(0);
00221 }
00222 
00223 int HexBeamCR::getNumBCNodes()
00224 {
00225   int numBCNodes = (numLocalSlices_+1)*(W_+1);
00226   return( numBCNodes );
00227 }
00228 
00229 int HexBeamCR::getBCNodes(int numNodes, int* nodeIDs)
00230 {
00231   if (numNodes != getNumBCNodes()) {
00232     return(-1);
00233   }
00234 
00235   int firstBCNode = firstLocalNode_ + W_;
00236 
00237   for(int i=0; i<numNodes; ++i) {
00238     nodeIDs[i] = firstBCNode + W_+1;
00239   }
00240 
00241   return(0);
00242 }
00243 
00244 int HexBeamCR::getBCGammaValues(int numBCDofs, double* gamma)
00245 {
00246   if (numBCDofs != getNumBCNodes()*dofPerNode_) {
00247     return(-1);
00248   }
00249 
00250   for(int i=0; i<numBCDofs; ++i) {
00251     gamma[i] = 2.0;
00252   }
00253 
00254   return(0);
00255 }
00256 
00257 int HexBeamCR::getNumSharedNodes()
00258 {
00259   if (numProcs_ < 2) return(0);
00260 
00261   int numSharedNodes = numNodesPerSlice_;
00262   if (localProc_ > 0 && localProc_ < numProcs_-1) {
00263     numSharedNodes += numNodesPerSlice_;
00264   }
00265 
00266   return(numSharedNodes);
00267 }
00268 
00269 int HexBeamCR::getSharedNodes(int numSharedNodes,
00270           int*& sharedNodes,
00271           int*& numSharingProcsPerNode,
00272           int**& sharingProcs)
00273 {
00274   if (numProcs_ < 2) return(0);
00275 
00276   if (numSharedNodes != getNumSharedNodes()) {
00277     return(-1);
00278   }
00279 
00280   sharedNodes = new int[numSharedNodes];
00281   numSharingProcsPerNode = new int[numSharedNodes];
00282   sharingProcs = new int*[numSharedNodes];
00283   int* sharingProcVals = new int[numSharedNodes];
00284   if (sharedNodes == NULL || numSharingProcsPerNode == NULL ||
00285       sharingProcs == NULL || sharingProcVals == NULL) {
00286     return(-1);
00287   }
00288 
00289   int i;
00290   for(i=0; i<numSharedNodes; ++i) {
00291     numSharingProcsPerNode[i] = 1;
00292     sharingProcs[i] = &(sharingProcVals[i]);
00293   }
00294 
00295   int firstSharedNode = firstLocalNode_+numNodesPerSlice_*(numLocalSlices_+2);
00296   int offset = 0;
00297   //FEI_COUT << localProc_ << ": firstSharedNode: " << firstSharedNode << FEI_ENDL;
00298   if (localProc_ < numProcs_ - 1) {
00299     for(i=0; i<numNodesPerSlice_; ++i) {
00300       sharedNodes[offset] = firstSharedNode+i;
00301       sharingProcs[offset++][0] = localProc_+1;
00302     }
00303   }
00304 
00305   firstSharedNode = firstLocalNode_;
00306   //FEI_COUT << localProc_ << ":+1 firstSharedNode: " << firstSharedNode << FEI_ENDL;
00307   if (localProc_ > 0) {
00308     for(i=0; i<numNodesPerSlice_; ++i) {
00309       sharedNodes[offset] = firstSharedNode+i;
00310       sharingProcs[offset++][0] = localProc_-1;
00311     }
00312   }
00313 
00314   return(0);
00315 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends