FEI Version of the Day
test_LinearSystem.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 #include <fei_macros.hpp>
00044 
00045 #include <test_utils/test_LinearSystem.hpp>
00046 #include <test_utils/test_VectorSpace.hpp>
00047 #include <test_utils/test_MatrixGraph.hpp>
00048 #include <snl_fei_Factory.hpp>
00049 #include <fei_Vector.hpp>
00050 #include <fei_Matrix.hpp>
00051 #include <snl_fei_LinearSystem_General.hpp>
00052 
00053 #include <test_utils/LibraryFactory.hpp>
00054 
00055 #ifdef HAVE_FEI_AZTECOO
00056 #include <fei_Aztec_LinSysCore.hpp>
00057 #endif
00058 
00059 #undef fei_file
00060 #define fei_file "test_LinearSystem.cpp"
00061 #include <fei_ErrMacros.hpp>
00062 
00063 test_LinearSystem::test_LinearSystem(MPI_Comm comm)
00064  : tester(comm)
00065 {
00066 }
00067 
00068 test_LinearSystem::~test_LinearSystem()
00069 {
00070 }
00071 
00072 int test_LinearSystem::runtests()
00073 {
00074   CHK_ERR( test1() );
00075   CHK_ERR( test2() );
00076   CHK_ERR( test3() );
00077   CHK_ERR( test4() );
00078   CHK_ERR( test5() );
00079 
00080   return(0);
00081 }
00082 
00083 int test_LinearSystem::test1()
00084 {
00085   return(0);
00086 }
00087 
00088 int test_LinearSystem::test2()
00089 {
00090 #ifdef HAVE_FEI_AZTECOO
00091   fei::SharedPtr<testData> testdata(new testData(localProc_, numProcs_));
00092   std::vector<int>& fieldIDs = testdata->fieldIDs;
00093   std::vector<int>& idTypes = testdata->idTypes;
00094   std::vector<int>& ids = testdata->ids;
00095 
00096   fei::SharedPtr<LinearSystemCore> az_lsc(new fei_trilinos::Aztec_LinSysCore(comm_));
00097 
00098   char* param = new char[64];
00099   sprintf(param,"debugOutput .");
00100 
00101   CHK_ERR( az_lsc->parameters(1, &param) );
00102   delete [] param;
00103 
00104   fei::SharedPtr<fei::Factory> factory(new snl_fei::Factory(comm_, az_lsc));
00105 
00106   fei::SharedPtr<fei::VectorSpace> vectorSpacePtr =
00107     test_VectorSpace::create_VectorSpace(comm_,
00108            testdata.get(), localProc_, numProcs_,
00109            false, false, "U_LS2", factory);
00110 
00111   fei::SharedPtr<fei::MatrixGraph> matrixGraphPtr =
00112     test_MatrixGraph::create_MatrixGraph(testdata.get(), localProc_, numProcs_,
00113            false, false, "U_LS2", vectorSpacePtr,
00114            factory, path_);
00115 
00116   std::vector<int> crIDTypes(2);
00117   std::vector<int> crFieldIDs(2);
00118   crIDTypes[0] = idTypes[0]; crIDTypes[1] = idTypes[0];
00119   crFieldIDs[0] = fieldIDs[0]; crFieldIDs[1] = fieldIDs[0];
00120 
00121   CHK_ERR( matrixGraphPtr->initLagrangeConstraint(0, idTypes[1],
00122               2, //numIDs
00123               &crIDTypes[0],
00124               &(ids[1]),
00125               &crFieldIDs[0]) );
00126 
00127   CHK_ERR( matrixGraphPtr->initComplete() );
00128 
00129   fei::SharedPtr<fei::Vector> vec_lsc = factory->createVector(vectorSpacePtr);
00130 
00131   fei::SharedPtr<fei::Vector> vec_lsc2 = factory->createVector(vectorSpacePtr, true);
00132 
00133   fei::SharedPtr<fei::Matrix> mat_lsc = factory->createMatrix(matrixGraphPtr);
00134 
00135   fei::SharedPtr<fei::LinearSystem> linsys = factory->createLinearSystem(matrixGraphPtr);
00136   linsys->setMatrix(mat_lsc);
00137   linsys->setSolutionVector(vec_lsc2);
00138   linsys->setRHS(vec_lsc);
00139 
00140   int blockID=0;
00141   int numIndices = matrixGraphPtr->getConnectivityNumIndices(blockID);
00142 
00143   std::vector<int> indicesArray(numIndices);
00144   int* indicesPtr = &indicesArray[0];
00145 
00146   int checkNumIndices = 0;
00147   CHK_ERR( matrixGraphPtr->getConnectivityIndices(blockID, 0,
00148                numIndices, indicesPtr,
00149                checkNumIndices) );
00150 
00151   std::vector<double> data(ids.size(), 1.0);
00152   double* dptr = &data[0];
00153   std::vector<double*> coefPtrs(ids.size());
00154   std::vector<double> crdata(2);
00155   crdata[0] = 1.0;
00156   crdata[1] = -1.0;
00157 
00158   for(unsigned ii=0; ii<ids.size(); ++ii) coefPtrs[ii] = dptr;
00159 
00160   CHK_ERR( mat_lsc->sumIn(numIndices, indicesPtr, numIndices, indicesPtr,
00161         &coefPtrs[0]) );
00162 
00163   CHK_ERR( vec_lsc->sumInFieldData(fieldIDs[0], idTypes[0],
00164             ids.size(), &ids[0],
00165             &data[0]) );
00166 
00167   CHK_ERR( linsys->loadLagrangeConstraint(0, &crdata[0], 0.0) );
00168 
00169   CHK_ERR( mat_lsc->gatherFromOverlap() );
00170 
00171   CHK_ERR( az_lsc->matrixLoadComplete() );
00172 
00173    CHK_ERR( linsys->loadComplete() );
00174 
00175   std::vector<int> crindices;
00176   linsys->getConstrainedEqns(crindices);
00177   if (crindices.size() != 2) {
00178     ERReturn(-7);
00179   }
00180 
00181   CHK_ERR( az_lsc->writeSystem("U_LS2") );
00182 
00183   MPI_Barrier(comm_);
00184 #endif
00185 
00186   return(0);
00187 }
00188 
00189 int test_LinearSystem::test3()
00190 {
00191 #ifdef HAVE_FEI_AZTECOO
00192   fei::SharedPtr<testData> testdata(new testData(localProc_, numProcs_));
00193   std::vector<int>& fieldIDs = testdata->fieldIDs;
00194   std::vector<int>& idTypes = testdata->idTypes;
00195   std::vector<int>& ids = testdata->ids;
00196 
00197   fei::SharedPtr<LinearSystemCore> az_lsc(new fei_trilinos::Aztec_LinSysCore(comm_));
00198 
00199   char* param = new char[64];
00200   sprintf(param,"debugOutput .");
00201 
00202   CHK_ERR( az_lsc->parameters(1, &param) );
00203 
00204   fei::SharedPtr<fei::Factory> factory(new snl_fei::Factory(comm_, az_lsc));
00205 
00206   fei::SharedPtr<fei::VectorSpace> vectorSpacePtr =
00207     test_VectorSpace::create_VectorSpace(comm_,
00208            testdata.get(), localProc_, numProcs_,
00209            false, false, "U_LS3", factory);
00210 
00211   fei::SharedPtr<fei::MatrixGraph> matrixGraphPtr =
00212     test_MatrixGraph::create_MatrixGraph(testdata.get(), localProc_, numProcs_,
00213            false, false, "U_LS3", vectorSpacePtr,
00214            factory, path_);
00215 
00216   std::vector<int> crIDTypes(2);
00217   std::vector<int> crFieldIDs(2);
00218   crIDTypes[0] = idTypes[0]; crIDTypes[1] = idTypes[0];
00219   crFieldIDs[0] = fieldIDs[0]; crFieldIDs[1] = fieldIDs[0];
00220 
00221   CHK_ERR( matrixGraphPtr->initPenaltyConstraint(0, idTypes[1],
00222               2, //numIDs
00223               &crIDTypes[0],
00224               &(ids[1]),
00225               &crFieldIDs[0]) );
00226 
00227   CHK_ERR( matrixGraphPtr->initComplete() );
00228 
00229   fei::SharedPtr<fei::Vector> vec_lsc = factory->createVector(vectorSpacePtr);
00230 
00231   fei::SharedPtr<fei::Vector> vec_lsc2 = factory->createVector(vectorSpacePtr, true);
00232 
00233   fei::SharedPtr<fei::Matrix> mat_lsc = factory->createMatrix(matrixGraphPtr);
00234 
00235   fei::SharedPtr<fei::LinearSystem> linsys = factory->createLinearSystem(matrixGraphPtr);
00236   CHK_ERR( linsys->parameters(1, &param) );
00237   delete [] param;
00238 
00239   linsys->setMatrix(mat_lsc);
00240   linsys->setSolutionVector(vec_lsc2);
00241   linsys->setRHS(vec_lsc);
00242 
00243   int blockID=0;
00244   int numIndices = matrixGraphPtr->getConnectivityNumIndices(blockID);
00245 
00246   std::vector<int> indicesArray(numIndices);
00247   int* indicesPtr = &indicesArray[0];
00248 
00249   int checkNumIndices = 0;
00250   CHK_ERR( matrixGraphPtr->getConnectivityIndices(blockID, 0,
00251                numIndices, indicesPtr,
00252                checkNumIndices) );
00253 
00254   std::vector<double> data(ids.size(), 1.0);
00255   double* dptr = &data[0];
00256   std::vector<double*> coefPtrs(ids.size());
00257   std::vector<double> crdata(2);
00258   crdata[0] = 1.0;
00259   crdata[1] = -1.0;
00260 
00261   for(unsigned ii=0; ii<ids.size(); ++ii) coefPtrs[ii] = dptr;
00262 
00263   CHK_ERR( mat_lsc->sumIn(numIndices, indicesPtr, numIndices, indicesPtr,
00264         &coefPtrs[0]) );
00265 
00266   CHK_ERR( vec_lsc->sumInFieldData(fieldIDs[0], idTypes[0],
00267             ids.size(), &ids[0],
00268             &data[0]) );
00269 
00270   CHK_ERR( linsys->loadPenaltyConstraint(0, &crdata[0], 100.0, 0.0) );
00271 
00272   CHK_ERR( mat_lsc->gatherFromOverlap() );
00273 
00274   CHK_ERR( az_lsc->matrixLoadComplete() );
00275 
00276    CHK_ERR( linsys->loadComplete() );
00277 
00278   CHK_ERR( az_lsc->writeSystem("U_LS3") );
00279 
00280   MPI_Barrier(comm_);
00281 #endif
00282 
00283   return(0);
00284 }
00285 
00286 int test_LinearSystem::test4()
00287 {
00288 #ifdef HAVE_FEI_AZTECOO
00289   fei::SharedPtr<testData> testdata(new testData(localProc_, numProcs_));
00290   std::vector<int>& fieldIDs = testdata->fieldIDs;
00291   std::vector<int>& idTypes = testdata->idTypes;
00292   std::vector<int>& ids = testdata->ids;
00293 
00294   fei::SharedPtr<LinearSystemCore> az_lsc(new fei_trilinos::Aztec_LinSysCore(comm_));
00295 
00296   char* param = new char[64];
00297   sprintf(param,"debugOutput .");
00298 
00299   CHK_ERR( az_lsc->parameters(1, &param) );
00300   delete [] param;
00301 
00302   fei::SharedPtr<fei::Factory> factory(new snl_fei::Factory(comm_, az_lsc));
00303 
00304   fei::SharedPtr<fei::VectorSpace> vectorSpacePtr =
00305     test_VectorSpace::create_VectorSpace(comm_,
00306            testdata.get(), localProc_, numProcs_,
00307            false, false, "U_LS4", factory);
00308 
00309   fei::SharedPtr<fei::MatrixGraph> matrixGraphPtr =
00310     test_MatrixGraph::create_MatrixGraph(testdata.get(), localProc_, numProcs_,
00311            false, false, "U_LS4", vectorSpacePtr,
00312            factory, path_);
00313 
00314   std::vector<int> crIDTypes(2);
00315   std::vector<int> crFieldIDs(2);
00316   crIDTypes[0] = idTypes[0]; crIDTypes[1] = idTypes[0];
00317   crFieldIDs[0] = fieldIDs[0]; crFieldIDs[1] = fieldIDs[0];
00318 
00319   CHK_ERR( matrixGraphPtr->initLagrangeConstraint(0, idTypes[1],
00320               2, //numIDs
00321               &crIDTypes[0],
00322               &(ids[1]),
00323               &crFieldIDs[0]) );
00324 
00325   CHK_ERR( matrixGraphPtr->initComplete() );
00326 
00327   fei::SharedPtr<fei::Vector> vec_lsc = factory->createVector(vectorSpacePtr);
00328 
00329   fei::SharedPtr<fei::Vector> vec_lsc2 = factory->createVector(vectorSpacePtr, true);
00330 
00331   fei::SharedPtr<fei::Matrix> mat_lsc = factory->createMatrix(matrixGraphPtr);
00332 
00333   fei::SharedPtr<fei::LinearSystem> linsys = factory->createLinearSystem(matrixGraphPtr);
00334   linsys->setMatrix(mat_lsc);
00335   linsys->setSolutionVector(vec_lsc2);
00336   linsys->setRHS(vec_lsc);
00337 
00338   int blockID=0;
00339   int numIndices = matrixGraphPtr->getConnectivityNumIndices(blockID);
00340 
00341   std::vector<int> indicesArray(numIndices);
00342   int* indicesPtr = &indicesArray[0];
00343 
00344   int checkNumIndices = 0;
00345   CHK_ERR( matrixGraphPtr->getConnectivityIndices(blockID, 0,
00346                numIndices, indicesPtr,
00347                checkNumIndices) );
00348 
00349   std::vector<double> data(ids.size(), 1.0);
00350   double* dptr = &data[0];
00351   std::vector<double*> coefPtrs(ids.size());
00352   std::vector<double> crdata(2);
00353   crdata[0] = 1.0;
00354   crdata[1] = -1.0;
00355 
00356   for(unsigned ii=0; ii<ids.size(); ++ii) coefPtrs[ii] = dptr;
00357 
00358   CHK_ERR( mat_lsc->sumIn(numIndices, indicesPtr, numIndices, indicesPtr,
00359         &coefPtrs[0]) );
00360 
00361   CHK_ERR( vec_lsc->sumInFieldData(fieldIDs[0], idTypes[0],
00362            ids.size(), &ids[0], &data[0]) );
00363 
00364   CHK_ERR( linsys->loadLagrangeConstraint(0, &crdata[0], 0.0) );
00365 
00366   CHK_ERR( mat_lsc->gatherFromOverlap() );
00367 
00368   CHK_ERR( az_lsc->matrixLoadComplete() );
00369 
00370   CHK_ERR( linsys->loadComplete() );
00371 
00372   CHK_ERR( az_lsc->writeSystem("U_LS4") );
00373 
00374   MPI_Barrier(comm_);
00375 #endif
00376 
00377   return(0);
00378 }
00379 
00380 int test_LinearSystem::test5()
00381 {
00382   return(0);
00383 }
00384 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends