OpenSees » SRC » system_of_eqn » linearSOE » sparseSYM

SymSparseLinSolver.cpp

Go to the documentation of this file.
00001 // File: ~/system_of_eqn/linearSOE/symLinSolver/SymSparseLinSolver.C
00002 //
00003 // Written: Jun Peng  (junpeng@stanford.edu)
00004 //          Prof. Kincho H. Law
00005 //          Stanford University
00006 // Created: 12/98
00007 // Revision: A
00008 //
00009 // Description: This file contains the class definition for 
00010 // SymSparseinSolver. It solves the SymSparseLinSOE object by calling
00011 // some "C" functions. The solver used here is generalized sparse
00012 // solver. The user can choose three different ordering schema.
00013 //
00014 // What: "@(#) SymSparseLinSolver.C, revA"
00015 
00016 
00017 #include "SymSparseLinSOE.h"
00018 #include "SymSparseLinSolver.h"
00019 #include <f2c.h>
00020 #include <math.h>
00021 #include <Channel.h>
00022 #include <FEM_ObjectBroker.h>
00023 
00024 extern "C" {
00025   #include "FeStructs.h"
00026 }
00027 
00028 
00029 
00030 SymSparseLinSolver::SymSparseLinSolver()
00031 :LinearSOESolver(SOLVER_TAGS_SymSparseLinSolver),
00032  theSOE(0)
00033 {
00034     // nothing to do.
00035 }
00036 
00037 
00038 SymSparseLinSolver::~SymSparseLinSolver()
00039 { 
00040     // nothing to do.
00041 }
00042 
00043 
00044 extern "C" int pfsfct(int neqns, double *diag, double **penv, int nblks, int *xblk,
00045                       OFFDBLK **begblk, OFFDBLK *first, int *rowblks);
00046 
00047 extern "C" void pfsslv(int neqns, double *diag, double **penv, int nblks,
00048                        int *xblk, double *rhs, OFFDBLK **begblk);
00049 
00050 
00051 int
00052 SymSparseLinSolver::solve(void)
00053 { 
00054     if (theSOE == 0) {
00055         opserr << "WARNING SymSparseLinSolver::solve(void)- ";
00056         opserr << " No LinearSOE object has been set\n";
00057         return -1;
00058     }
00059 
00060     int      nblks = theSOE->nblks;
00061     int      *xblk = theSOE->xblk;
00062     int      *invp = theSOE->invp;
00063     double   *diag = theSOE->diag;
00064     double   **penv = theSOE->penv;
00065     int      *rowblks = theSOE->rowblks;
00066     OFFDBLK  **begblk = theSOE->begblk;
00067     OFFDBLK  *first = theSOE->first;
00068 
00069     int neq = theSOE->size;
00070 
00071     // check for quick return
00072     if (neq == 0)
00073         return 0;
00074 
00075     // first copy B into X
00076 
00077     for (int i=0; i<neq; i++) {
00078         theSOE->X[i] = theSOE->B[i];
00079     }
00080     double *Xptr = theSOE->X;
00081 
00082     if (theSOE->factored == false) {
00083 
00084         //factor the matrix
00085         //call the "C" function to do the numerical factorization.
00086         int factor;
00087         factor = pfsfct(neq, diag, penv, nblks, xblk, begblk, first, rowblks);
00088         if (factor > 0) {
00089             opserr << "In SymSparseLinSolver: error in factorization.\n";
00090             return -1;
00091         }
00092         theSOE->factored = true;
00093     }
00094 
00095     // do forward and backward substitution.
00096     // call the "C" function.
00097 
00098     pfsslv(neq, diag, penv, nblks, xblk, Xptr, begblk);
00099 
00100     // Since the X we get by solving AX=B is P*X, we need to reordering
00101     // the Xptr to ge the wanted X.
00102 
00103     double *tempX = new double[neq];
00104     if (tempX == 0) {
00105         opserr << "WARNING SymSparseLinSover::SymSparseLinSolver :";
00106         opserr << " ran out of memory for vectors (tempX) ";
00107         return -1;
00108     } 
00109 
00110     for (int m=0; m<neq; m++) {
00111         tempX[m] = Xptr[invp[m]];
00112     }
00113         
00114     for (int k=0; k<neq; k++) {
00115         Xptr[k] = tempX[k];
00116     }
00117         
00118     delete [] tempX;
00119     return 0;
00120 }
00121 
00122 
00123 int
00124 SymSparseLinSolver::setSize()
00125 {
00126     // nothing to do
00127     return 0;
00128 }
00129 
00130 
00131 int
00132 SymSparseLinSolver::setLinearSOE(SymSparseLinSOE &theLinearSOE)
00133 {
00134     theSOE = &theLinearSOE;
00135     return 0;
00136 }
00137 
00138 
00139 int
00140 SymSparseLinSolver::sendSelf(int cTAg, Channel &theChannel)
00141 {
00142     // doing nothing
00143     return 0;
00144 }
00145 
00146 
00147 int
00148 SymSparseLinSolver::recvSelf(int cTag,
00149                              Channel &theChannel, FEM_ObjectBroker &theBroker)
00150 {
00151     // nothing to do
00152     return 0;
00153 }
00154 
00155 
00156 
00157

Generated on Mon Oct 23 15:05:29 2006 for OpenSees by doxygen 1.5.0