C++ LinearSOE::getNumEqn方法代码示例

int 
SecantAccelerator2::newStep(LinearSOE &theSOE)
{
  int newNumEqns = theSOE.getNumEqn();

  if (vOld != 0 && vOld->Size() != newNumEqns) {
    delete vOld;
    vOld = 0;
  }

  if (rOld != 0 && rOld->Size() != newNumEqns) {
    delete rOld;
    rOld = 0;
  }

  numEqns = newNumEqns;

  if (vOld == 0)
    vOld = new Vector(numEqns);

  if (rOld == 0)
    rOld = new Vector(numEqns);

  // Reset iteration counter
  iteration = 0;

  return 0;
}

int 
AlgorithmIncrements::playback(int cTag)
{

  if (fileName != 0) {

    LinearSOE *theSOE = theAlgo->getLinearSOEptr();
    if (theSOE == 0)
      return 0;

    Vector X(theSOE->getNumEqn());
    Vector B(theSOE->getNumEqn());

    if (theFile) {
      theFile.close();
    }

    ifstream aFile;
    aFile.open(fileName, ios::in);
    if (!aFile) {
      opserr << "WARNING - AlgorithmIncrements::playback()";
      opserr << " - could not open file " << fileName << endln;
    }
    for (int i = 0; i<numRecord; i++) {
      int ii;
      for (ii=0; X.Size(); ii++) theFile << X(ii);
      for (ii=0; X.Size(); ii++) theFile << B(ii);

      this->plotData(X,B);
      //      char c = getchar();
    }
  }

  return 0;
}

int KRAlphaExplicit::formTangent(int statFlag)
{
    statusFlag = statFlag;
    
    LinearSOE *theLinSOE = this->getLinearSOE();
    AnalysisModel *theModel = this->getAnalysisModel();
    if (theLinSOE == 0 || theModel == 0)  {
        opserr << "WARNING KRAlphaExplicit::formTangent() - ";
        opserr << "no LinearSOE or AnalysisModel has been set\n";
        return -1;
    }
    
    theLinSOE->zeroA();
    
    int size = theLinSOE->getNumEqn();
    ID id(size);
    for (int i=1; i<size; i++)  {
        id(i) = id(i-1) + 1;
    }
    if (theLinSOE->addA(*Mhat, id) < 0)  {
        opserr << "WARNING KRAlphaExplicit::formTangent() - ";
        opserr << "failed to add Mhat to A\n";
        return -2;
    }
    
    return 0;
}

int 
MillerAccelerator::newStep(LinearSOE &theSOE)
{
  int newNumEqns = theSOE.getNumEqn();

  if (newNumEqns != numEqns) {
    if (fData != 0) {
      delete [] fData;
      fData = 0;
    }
    if (work != 0) {
      delete [] work;
      work = 0;
    }
    numEqns = newNumEqns;
  }

  if (fData == 0)
    fData = new double [numEqns];

  // Make sure max dim <= num equations
  if (maxDimension > numEqns)
    maxDimension = numEqns;

  // Length of work array -- N*(2*MVEC+2)
  int lwork = numEqns*(2*maxDimension+2);
  
  if (work == 0)
    work = new double [lwork];

  // Reset iteration counter
  iteration = 1;
  dimension = (theTangent == CURRENT_TANGENT) ? maxDimension : 0;

  return 0;
}

int 
BFGS::solveCurrentStep(void)
{
 
    // set up some pointers and check they are valid
    // NOTE this could be taken away if we set Ptrs as protecetd in superclass

    AnalysisModel   *theAnaModel = this->getAnalysisModelPtr();

    IncrementalIntegrator *theIntegrator = this->getIncrementalIntegratorPtr();

    LinearSOE  *theSOE = this->getLinearSOEptr();

    if ((theAnaModel == 0) || (theIntegrator == 0) || (theSOE == 0)
	|| (theTest == 0)){
	opserr << "WARNING BFGS::solveCurrentStep() - setLinks() has";
	opserr << " not been called - or no ConvergenceTest has been set\n";
	return -5;
    }	

    // set itself as the ConvergenceTest objects EquiSolnAlgo
    theTest->setEquiSolnAlgo(*this);
    if (theTest->start() < 0) {
      opserr << "BFGS::solveCurrentStep() -";
      opserr << "the ConvergenceTest object failed in start()\n";
      return -3;
    }

    localTest->setEquiSolnAlgo(*this);

    if (rdotz == 0)
       rdotz = new double[numberLoops+3];

    if (sdotr == 0)
	sdotr = new double[numberLoops+3];


    int result = -1;
    int count = 0;
    do {

      // opserr << "      BFGS -- Forming New Tangent" << endln;

      //form the initial tangent
      if (theIntegrator->formTangent(tangent) < 0){
         opserr << "WARNING BFGS::solveCurrentStep() -";
         opserr << "the Integrator failed in formTangent()\n";
         return -1; 
      }

      //form the initial residual 
      if (theIntegrator->formUnbalance() < 0) {
        opserr << "WARNING BFGS::solveCurrentStep() -";
        opserr << "the Integrator failed in formUnbalance()\n";	
      }	    

      //solve
      if (theSOE->solve() < 0) {
	  opserr << "WARNING BFGS::solveCurrentStep() -";
	  opserr << "the LinearSysOfEqn failed in solve()\n";	
	  return -3;
	}	    

      //update
      if ( theIntegrator->update(theSOE->getX() ) < 0) {
	opserr << "WARNING BFGS::solveCurrentStep() -";
	opserr << "the Integrator failed in update()\n";	
	return -4;
      }	        


      //    int systemSize = ( theSOE->getB() ).Size();
      int systemSize = theSOE->getNumEqn( );

      //temporary vector
      if (temp == 0 )
	temp = new Vector(systemSize);

      //initial displacement increment
      if ( s[1] == 0 ) 
	s[1] = new Vector(systemSize);

      *s[1] = theSOE->getX( );

      if ( residOld == 0 ) 
	residOld = new Vector(systemSize);

      *residOld = theSOE->getB( ) ;
      *residOld *= (-1.0 );

      //form the residual again
      if (theIntegrator->formUnbalance() < 0) {
        opserr << "WARNING BFGS::solveCurrentStep() -";
        opserr << "the Integrator failed in formUnbalance()\n";	
      }	    

      if ( residNew == 0 ) 
	residNew = new Vector(systemSize);
 
      if ( du == 0 ) 
//.........这里部分代码省略.........

int 
KrylovAccelerator2::newStep(LinearSOE &theSOE)
{
  int newNumEqns = theSOE.getNumEqn();

  if (numEqns != newNumEqns) {
    if (v != 0) {
      for (int i = 0; i < maxDimension+1; i++)
	delete v[i];
      delete [] v;
      v = 0;
    }
    
    if (Av != 0) {
      for (int i = 0; i < maxDimension+1; i++)
	delete Av[i];
      delete [] Av;
      Av = 0;
    }
    
    if (AvData != 0) {
      delete [] AvData;
      AvData = 0;
    }
    
    if (rData != 0) {
      delete [] rData;
      rData = 0;
    }
    
    if (work != 0) {
      delete [] work;
      work = 0;
    }
  }

  numEqns = newNumEqns;
  if (maxDimension > numEqns)
    maxDimension = numEqns;

  if (v == 0) {
    v = new Vector*[maxDimension+1];
    for (int i = 0; i < maxDimension+1; i++)
      v[i] = new Vector(numEqns);
  }

  if (Av == 0) {
    Av = new Vector*[maxDimension+1];
    for (int i = 0; i < maxDimension+1; i++)
      Av[i] = new Vector(numEqns);
  }

  if (AvData == 0)
    AvData = new double [maxDimension*numEqns];

  if (rData == 0)
    // The LAPACK least squares subroutine overwrites the RHS vector
    // with the solution vector ... these vectors are not the same
    // size, so we need to use the max size
    rData = new double [(numEqns > maxDimension) ? numEqns : maxDimension];

  // Length of work vector should be >= 2*min(numEqns,maxDimension)
  // See dgels subroutine documentation
  lwork = 2 * ((numEqns < maxDimension) ? numEqns : maxDimension);
  if (lwork < 1)
    lwork = 1;

  if (work == 0)
    work = new double [lwork];

  // Reset dimension of subspace
  dimension = 0;

  return 0;
}

int KRAlphaExplicit::newStep(double _deltaT)
{
    updateCount = 0;
    
    if (beta == 0 || gamma == 0 )  {
        opserr << "WARNING KRAlphaExplicit::newStep() - error in variable\n";
        opserr << "gamma = " << gamma << " beta = " << beta << endln;
        return -1;
    }
    
    // get a pointer to the AnalysisModel
    AnalysisModel *theModel = this->getAnalysisModel();
    if (theModel == 0)  {
        opserr << "WARNING KRAlphaExplicit::newStep() - no AnalysisModel set\n";
        return -2;
    }
    
    if (initAlphaMatrices || _deltaT != deltaT)  {
        
        // update time step increment
        deltaT = _deltaT;
        if (deltaT <= 0.0)  {
            opserr << "WARNING KRAlphaExplicit::newStep() - error in variable\n";
            opserr << "dT = " << deltaT << endln;
            return -3;
        }
        
        // get the LinearSOE and the ConvergenceTest so we can switch back later
        LinearSOE *theLinSOE = this->getLinearSOE();
        ConvergenceTest *theTest = this->getConvergenceTest();
        
        // set up the FullLinearSOE (needed to compute the alpha matrices)
        int size = theLinSOE->getNumEqn();
        FullGenLinSolver *theFullLinSolver = new FullGenLinLapackSolver();
        LinearSOE *theFullLinSOE = new FullGenLinSOE(size, *theFullLinSolver);
        if (theFullLinSOE == 0)  {
            opserr << "WARNING KRAlphaExplicit::newStep() - failed to create FullLinearSOE\n";
            return -4;
        }
        theFullLinSOE->setLinks(*theModel);
        
        // now switch the SOE to the FullLinearSOE
        this->IncrementalIntegrator::setLinks(*theModel, *theFullLinSOE, theTest);
        
        // get a pointer to the A matrix of the FullLinearSOE
        const Matrix *tmp = theFullLinSOE->getA();
        if (tmp == 0)  {
            opserr << "WARNING KRAlphaExplicit::domainChanged() - ";
            opserr << "failed to get A matrix of FullGeneral LinearSOE\n";
            return -5;
        
        }
        
        // calculate the integration parameter matrices
        c1 = beta*deltaT*deltaT;
        c2 = gamma*deltaT;
        c3 = 1.0;
        this->TransientIntegrator::formTangent(INITIAL_TANGENT);
        Matrix A(*tmp);
        
        c1 *= (1.0 - alphaF);
        c2 *= (1.0 - alphaF);
        c3 = (1.0 -alphaM);
        this->TransientIntegrator::formTangent(INITIAL_TANGENT);
        Matrix B3(*tmp);
        
        // solve [M + gamma*deltaT*C + beta*deltaT^2*K]*[alpha3] = 
        // [alphaM*M + alphaF*gamma*deltaT*C + alphaF*beta*deltaT^2*K] for alpha3
        A.Solve(B3, *alpha3);
        
        c1 = 0.0;
        c2 = 0.0;
        c3 = 1.0;
        this->TransientIntegrator::formTangent(INITIAL_TANGENT);
        Matrix B1(*tmp);
        
        // solve [M + gamma*deltaT*C + beta*deltaT^2*K]*[alpha1] = [M] for alpha1
        A.Solve(B1, *alpha1);
        
        // calculate the effective mass matrix Mhat
        Mhat->addMatrix(0.0, B1, 1.0);
        Mhat->addMatrixProduct(1.0, B1, *alpha3, -1.0);
        
        // switch the SOE back to the user specified one
        this->IncrementalIntegrator::setLinks(*theModel, *theLinSOE, theTest);
        
        initAlphaMatrices = 0;
    }
    
    if (U == 0)  {
        opserr << "WARNING KRAlphaExplicit::newStep() - domainChange() failed or hasn't been called\n";
        return -6;
    }
    
    // set response at t to be that at t+deltaT of previous step
    (*Ut) = *U;
    (*Utdot) = *Udot;
    (*Utdotdot) = *Udotdot;
    
    // determine new response at time t+deltaT
//.........这里部分代码省略.........