FEI Version of the Day
Poisson_Elem.cpp
00001 /*--------------------------------------------------------------------*/
00002 /*    Copyright 2005 Sandia Corporation.                              */
00003 /*    Under the terms of Contract DE-AC04-94AL85000, there is a       */
00004 /*    non-exclusive license for use of this work by or on behalf      */
00005 /*    of the U.S. Government.  Export of this program may require     */
00006 /*    a license from the United States Government.                    */
00007 /*--------------------------------------------------------------------*/
00008 
00009 #include <fei_iostream.hpp>
00010 #include <fei_defs.h>
00011 #include <test_utils/Poisson_Elem.hpp>
00012 
00013 #include <cmath>
00014 #include <cstdlib>
00015 #include <assert.h>
00016 
00017 //==============================================================================
00018 Poisson_Elem::Poisson_Elem() {
00019 
00020     globalElemID_ = (GlobalID)(-1);
00021     ID_IsSet_ = false;
00022     numElemNodes_ = 4;
00023     numElemRows_ = 4;
00024 
00025     nodeList_ = NULL;
00026     nodalX_ = NULL;
00027     nodalY_ = NULL;
00028 
00029     elemStiff_ = NULL;
00030     elemLoad_ = NULL;
00031     internalsAllocated_ = false;
00032 
00033     elemLength_ = 0.0;
00034     elemLengthIsSet_ = false;
00035     totalLength_ = 0.0;
00036     totalLengthIsSet_ = false;
00037 
00038     loadAllocated_ = false;
00039     stiffAllocated_ = false;
00040 }
00041 
00042 //==============================================================================
00043 Poisson_Elem::~Poisson_Elem() {
00044    deleteMemory();
00045 }
00046 
00047 //==============================================================================
00048 void Poisson_Elem::deleteMemory()
00049 {
00050     if (!internalsAllocated_) return;
00051 
00052     delete [] nodeList_;
00053     delete [] nodalX_;
00054     delete [] nodalY_;
00055 
00056     internalsAllocated_ = false;
00057 
00058     if (loadAllocated_) {
00059        delete [] elemLoad_;
00060        loadAllocated_ = false;
00061     }
00062 
00063     if (stiffAllocated_) {
00064        delete [] elemStiff_[0];
00065        delete [] elemStiff_;
00066        stiffAllocated_ = false;
00067    }
00068 }
00069 
00070 //==============================================================================
00071 int Poisson_Elem::allocateInternals(int DOF) {
00072 
00073     assert(DOF==1);
00074 
00075     nodeList_ = new GlobalID[numElemNodes_];
00076     if (!nodeList_) return(1);
00077 
00078     nodalX_ = new double[numElemNodes_];
00079     if (!nodalX_) return(1);
00080 
00081     nodalY_ = new double[numElemNodes_];
00082     if (!nodalY_) return(1);
00083 
00084     for(int i = 0; i < numElemNodes_; i++) {
00085         nodeList_[i] = (GlobalID)0;
00086         nodalX_[i] = 0.0;
00087         nodalY_[i] = 0.0;
00088     }
00089 
00090     internalsAllocated_ = true;
00091     return(0);
00092 }
00093 
00094 //==============================================================================
00095 int Poisson_Elem::allocateLoad(int DOF)
00096 {
00097     assert(DOF==1);
00098 
00099     elemLoad_ = new double[numElemNodes_];
00100     if (!elemLoad_) return(1);
00101 
00102     loadAllocated_ = true;
00103     return(0);
00104 }
00105 
00106 //==============================================================================
00107 int Poisson_Elem::allocateStiffness(int DOF)
00108 {
00109     assert(DOF==1);
00110 
00111     elemStiff_ = new double*[numElemNodes_];
00112     if (!elemStiff_) return(1);
00113 
00114     elemStiff_[0] = NULL;
00115     elemStiff_[0] = new double[numElemNodes_*numElemNodes_];
00116     if (!elemStiff_[0]) return(1);
00117 
00118     int i;
00119     for(i=0; i<numElemNodes_*numElemNodes_; i++) elemStiff_[0][i] = 0.0;
00120 
00121     for(i=1; i<numElemNodes_; i++) {
00122         elemStiff_[i] = elemStiff_[i-1] + numElemNodes_;
00123     }
00124 
00125     stiffAllocated_ = true;
00126     return(0);
00127 }
00128 
00129 //==============================================================================
00130 GlobalID* Poisson_Elem::getElemConnPtr(int& size) {
00131 
00132     if (internalsAllocated_) {
00133         size = numElemNodes_;
00134         return(nodeList_);
00135     }
00136     else {
00137         size = 0;
00138         return(NULL);
00139     }
00140 }
00141 
00142 //==============================================================================
00143 double* Poisson_Elem::getElemLoad(int& size) {
00144         
00145     if (internalsAllocated_) {
00146         size = numElemNodes_;
00147         return(elemLoad_);
00148     }   
00149     else {
00150         size = 0;
00151         return(NULL);
00152     }
00153 }
00154 
00155 //==============================================================================
00156 double** Poisson_Elem::getElemStiff(int& size) {
00157         
00158     if (internalsAllocated_) {
00159         size = numElemNodes_*numElemNodes_;
00160         return(elemStiff_);
00161     }   
00162     else {
00163         size = 0;
00164         return(NULL);
00165     }
00166 }
00167 
00168 //==============================================================================
00169 void Poisson_Elem::calculateCoords() {
00170 //
00171 //This function calculates nodal x- and y-coordinates for this element.
00172 //NOTE: element IDs are assumed to be 1-based.
00173 //
00174     if (!internalsAllocated_)
00175         messageAbort("calculateCoords: internals not allocated.");
00176     if (!elemLengthIsSet_)
00177         messageAbort("calculateCoords: elemLength not set.");
00178     if (!totalLengthIsSet_)
00179         messageAbort("calculateCoords: totalLength not set.");
00180     if (!ID_IsSet_)
00181         messageAbort("calculateCoords: elemID not set.");
00182     if (std::abs(elemLength_) < 1.e-49)
00183         messageAbort("calculateCoords: elemLength == 0.");
00184 
00185     int lowLeft = 0;
00186     int lowRight = 1;
00187     int upperRight = 2;
00188     int upperLeft = 3;
00189 
00190     int elemsPerSide = (int)std::ceil(totalLength_/elemLength_);
00191 
00192     int elemX = (int)globalElemID_%elemsPerSide;
00193     if (elemX==0) elemX = elemsPerSide;
00194 
00195     int elemY = ((int)globalElemID_ - elemX)/elemsPerSide + 1;
00196 
00197     //elemX and elemY are 1-based coordinates of this element in
00198     //the global square of elements. The origin is position (1,1),
00199     //which is at the bottom left of the square.
00200 
00201     nodalX_[lowLeft] = (elemX-1)*elemLength_;
00202     nodalX_[upperLeft] = nodalX_[lowLeft];
00203     nodalX_[lowRight] = elemX*elemLength_;
00204     nodalX_[upperRight] = nodalX_[lowRight];
00205 
00206     nodalY_[lowLeft] = (elemY-1)*elemLength_;
00207     nodalY_[lowRight] = nodalY_[lowLeft];
00208     nodalY_[upperLeft] = elemY*elemLength_;
00209     nodalY_[upperRight] = nodalY_[upperLeft];
00210 }
00211 
00212 //==============================================================================
00213 void Poisson_Elem::messageAbort(const char* str) {
00214     fei::console_out() << "Poisson_Elem: ERROR: " << str << " Aborting." << FEI_ENDL;
00215     std::abort();
00216 }
00217 
00218 //==============================================================================
00219 void Poisson_Elem::calculateLoad() {
00220 
00221     assert (numElemRows_ == 4);
00222     
00223     elemLoad_[0] = -2.0*elemLength_*elemLength_;
00224     elemLoad_[1] = -2.0*elemLength_*elemLength_;
00225     elemLoad_[2] = -2.0*elemLength_*elemLength_;
00226     elemLoad_[3] = -2.0*elemLength_*elemLength_;
00227 }
00228 
00229 //==============================================================================
00230 void Poisson_Elem::calculateStiffness() {
00231 
00232     assert (numElemRows_ == 4);
00233     
00234     elemStiff_[0][0] = 2.0;
00235     elemStiff_[0][1] = -1.0;
00236     elemStiff_[0][2] = 0.0;
00237     elemStiff_[0][3] = -1.0;
00238 
00239     elemStiff_[1][0] = -1.0;
00240     elemStiff_[1][1] = 2.0;
00241     elemStiff_[1][2] = -1.0;
00242     elemStiff_[1][3] = 0.0;
00243 
00244     elemStiff_[2][0] = 0.0;
00245     elemStiff_[2][1] = -1.0;
00246     elemStiff_[2][2] = 2.0;
00247     elemStiff_[2][3] = -1.0;
00248 
00249     elemStiff_[3][0] = -1.0;
00250     elemStiff_[3][1] = 0.0;
00251     elemStiff_[3][2] = -1.0;
00252     elemStiff_[3][3] = 2.0;
00253 }
00254 
00255 
00256 
00257 //-----------------------------------------------------------------------
00258 
00259 void Poisson_Elem::dumpToScreen() {
00260 
00261     int i;
00262     
00263     FEI_COUT << " globalElemID_ = " << (int)globalElemID_ << FEI_ENDL;
00264     FEI_COUT << " numElemRows_  = " << numElemRows_ << FEI_ENDL;
00265     FEI_COUT << " elemLength_  = " << elemLength_ << FEI_ENDL;
00266     FEI_COUT << " the " << numElemNodes_ << " nodes: ";
00267     for (i = 0; i < numElemNodes_; ++i) {
00268       FEI_COUT << "   " << (int)nodeList_[i];
00269     }
00270     FEI_COUT << FEI_ENDL;
00271     FEI_COUT << " the " << numElemRows_ << " load vector terms: ";
00272     for (i = 0; i < numElemRows_; ++i) {
00273       FEI_COUT << "   " << elemLoad_[i];
00274     }
00275     FEI_COUT << FEI_ENDL << FEI_ENDL;
00276     return;
00277 }
00278 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends