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

int XC::Collocation::domainChanged()
  {
    AnalysisModel *myModel = this->getAnalysisModelPtr();
    LinearSOE *theLinSOE = this->getLinearSOEPtr();
    const Vector &x = theLinSOE->getX();
    int size = x.Size();
    
    setRayleighDampingFactors();
    
    if(Ut.get().Size() != size)
      {
        // create the new
        Ut.resize(size);
        U.resize(size);
      }
    
    // now go through and populate U, by iterating through
    // the DOF_Groups and getting the last committed velocity and accel
    DOF_GrpIter &theDOFGroups = myModel->getDOFGroups();
    DOF_Group *dofGroupPtr= nullptr;
    
    while((dofGroupPtr = theDOFGroups()) != 0)
      {
        const XC::ID &id = dofGroupPtr->getID();
        const Vector &disp = dofGroupPtr->getCommittedDisp();	
        U.setDisp(id,disp);
        const Vector &vel = dofGroupPtr->getCommittedVel();
        U.setVel(id,vel);
        const Vector &accel = dofGroupPtr->getCommittedAccel();
        U.setAccel(id,accel);
        // NOTE WE CAN't DO TOGETHER BECAUSE DOF_GROUPS USING SINGLE VECTOR
      }
    
    return 0;
  }

int
DisplacementControl::update(const Vector &dU)
{

  if (theDofID == -1) {
    opserr << "DisplacementControl::newStep() - domainChanged has not been called\n";
    return -1;
  } 
    AnalysisModel *theModel = this->getAnalysisModel();
    LinearSOE *theLinSOE = this->getLinearSOE();    
    if (theModel == 0 || theLinSOE == 0) {
	opserr << "WARNING DisplacementControl::update() ";
	opserr << "No AnalysisModel or LinearSOE has been set\n";
	return -1;
    }

    (*deltaUbar) = dU; // have to do this as the SOE is gonna change
    double dUabar = (*deltaUbar)(theDofID);
    
    // determine dUhat    
    theLinSOE->setB(*phat);
    theLinSOE->solve();
    (*deltaUhat) = theLinSOE->getX();    

    double dUahat = (*deltaUhat)(theDofID);
    if (dUahat == 0.0) {
	opserr << "WARNING DisplacementControl::update() ";
	opserr << "dUahat is zero -- zero reference displacement at control node DOF\n";
	return -1;
    }
    
    // determine delta lambda(1) == dlambda    
    double dLambda = -dUabar/dUahat;
    
    // determine delta U(i)
    (*deltaU) = (*deltaUbar);    
    deltaU->addVector(1.0, *deltaUhat,dLambda);
    
    // update dU and dlambda
    (*deltaUstep) += *deltaU;
    deltaLambdaStep += dLambda;
    currentLambda += dLambda;

    // update the model
    theModel->incrDisp(*deltaU);    
    theModel->applyLoadDomain(currentLambda);    
    if (theModel->updateDomain() < 0) {
      opserr << "DisplacementControl::update - model failed to update for new dU\n";
      return -1;
    }
	
    
    // set the X soln in linearSOE to be deltaU for convergence Test
    theLinSOE->setX(*deltaU);

    numIncrLastStep++;

    return 0;
}

//! If the size of the LinearSOE has changed, the object deletes any old Vectors
//! created and then creates \f$6\f$ new Vector objects of size equal to {\em
//! theLinearSOE-\f$>\f$getNumEqn()}. There is a Vector object created to store
//! the current displacement, velocity and accelerations at times \f$t\f$ and
//! \f$t + \Delta t\f$. The response quantities at time \f$t + \Delta t\f$ are
//! then set by iterating over the DOF\_Group objects in the model and
//! obtaining their committed values. 
//! Returns \f$0\f$ if successful, otherwise a warning message and a negative
//! number is returned: \f$-1\f$ if no memory was available for constructing
//! the Vectors.
int XC::CentralDifference::domainChanged(void)
  {
    AnalysisModel *myModel = this->getAnalysisModelPtr();
    LinearSOE *theLinSOE = this->getLinearSOEPtr();
    const XC::Vector &x = theLinSOE->getX();
    int size = x.Size();
    
    setRayleighDampingFactors();

    if(Ut.get().Size() != size)
      {
        Utm1.resize(size);
        Ut.resize(size);
        U.resize(size);
      }
    
    // now go through and populate U, Udot and Udotdot by iterating through
    // the DOF_Groups and getting the last committed velocity and accel
    DOF_GrpIter &theDOFGroups = myModel->getDOFGroups();
    DOF_Group *dofGroupPtr= nullptr;
    
    while((dofGroupPtr = theDOFGroups()) != 0)
      {
        const ID &id = dofGroupPtr->getID();
        const int idSize = id.Size();
        
        const Vector &disp = dofGroupPtr->getCommittedDisp();	
        for(int i=0;i<idSize;i++)
          {
            const int loc = id(i);
            if(loc >= 0)
              { Utm1(loc) = disp(i); }
          }
        Ut.setDisp(id,disp);
        
        const Vector &vel = dofGroupPtr->getCommittedVel();
        U.setVel(id,vel);
        
        const XC::Vector &accel = dofGroupPtr->getCommittedAccel();	
        U.setAccel(id,accel);
        /** NOTE WE CAN't DO TOGETHER BECAUSE DOF_GROUPS USING SINGLE VECTOR ******
        for (int i=0; i < id.Size(); i++) {
          int loc = id(i);
          if (loc >= 0) {
 	    (*U)(loc) = disp(i);		
             (U.getDot())(loc) = vel(i);
          (U.getDotDot())(loc) = accel(i);
           }
         }
          ***/
      }
    
    std::cerr << getClassName() << "::" << __FUNCTION__
	      << "; WARNING: assuming Ut-1 = Ut\n";
    
    return 0;
  }

int CollocationHSFixedNumIter::commit(void)
{
    AnalysisModel *theModel = this->getAnalysisModel();
    if (theModel == 0)  {
        opserr << "WARNING CollocationHSFixedNumIter::commit() - no AnalysisModel set\n";
        return -1;
    }	  
    
    LinearSOE *theSOE = this->getLinearSOE();
    if (theSOE == 0)  {
        opserr << "WARNING CollocationHSFixedNumIter::commit() - no LinearSOE set\n";
        return -2;
    }
    
    if (theSOE->solve() < 0)  {
        opserr << "WARNING CollocationHSFixedNumIter::commit() - "
            << "the LinearSysOfEqn failed in solve()\n";
        return -3;
    }
    const Vector &deltaU = theSOE->getX();
    
    // determine the response at t+theta*deltaT
    U->addVector(1.0, deltaU, c1);
    
    Udot->addVector(1.0, deltaU, c2);
    
    Udotdot->addVector(1.0, deltaU, c3);
    
    // determine response quantities at t+deltaT
    Udotdot->addVector(1.0/theta, *Utdotdot, (theta-1.0)/theta);
    
    (*Udot) = *Utdot;
    double a1 = deltaT*(1.0 - gamma);
    double a2 = deltaT*gamma;
    Udot->addVector(1.0, *Utdotdot, a1);
    Udot->addVector(1.0, *Udotdot, a2);
    
    (*U) = *Ut;
    U->addVector(1.0, *Utdot, deltaT);
    double a3 = deltaT*deltaT*(0.5 - beta);
    double a4 = deltaT*deltaT*beta;
    U->addVector(1.0, *Utdotdot, a3);
    U->addVector(1.0, *Udotdot, a4);
    
    // update the response at the DOFs
    theModel->setResponse(*U,*Udot,*Udotdot);
    
    // set the time to be t+deltaT
    double time = theModel->getCurrentDomainTime();
    time += (1.0-theta)*deltaT;
    theModel->setCurrentDomainTime(time);
    
    return theModel->commitDomain();
}

int
ArcLength::newStep(void)
{
    // get pointers to AnalysisModel and LinearSOE
    AnalysisModel *theModel = this->getAnalysisModel();
    LinearSOE *theLinSOE = this->getLinearSOE();    
    if (theModel == 0 || theLinSOE == 0) {
	opserr << "WARNING ArcLength::newStep() ";
	opserr << "No AnalysisModel or LinearSOE has been set\n";
	return -1;
    }

    // get the current load factor
    currentLambda = theModel->getCurrentDomainTime();

    if (deltaLambdaStep < 0)
	signLastDeltaLambdaStep = -1;
    else
	signLastDeltaLambdaStep = +1;

    // determine dUhat
    this->formTangent();
    theLinSOE->setB(*phat);
    if (theLinSOE->solve() < 0) {
      opserr << "ArcLength::newStep(void) - failed in solver\n";
      return -1;
    }

    (*deltaUhat) = theLinSOE->getX();
    Vector &dUhat = *deltaUhat;
    
    // determine delta lambda(1) == dlambda
    double dLambda = sqrt(arcLength2/((dUhat^dUhat)+alpha2));
    dLambda *= signLastDeltaLambdaStep; // base sign of load change
                                        // on what was happening last step
    deltaLambdaStep = dLambda;
    currentLambda += dLambda;

    // determine delta U(1) == dU
    (*deltaU) = dUhat;
    (*deltaU) *= dLambda;
    (*deltaUstep) = (*deltaU);

    // update model with delta lambda and delta U
    theModel->incrDisp(*deltaU);    
    theModel->applyLoadDomain(currentLambda);    
    theModel->updateDomain();

    return 0;
}

int
Linear::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 *theAnalysisModel = this->getAnalysisModelPtr();
    LinearSOE  *theSOE = this->getLinearSOEptr();
    IncrementalIntegrator  *theIncIntegrator =
        this->getIncrementalIntegratorPtr();

    if ((theAnalysisModel == 0) || (theIncIntegrator ==0 ) || (theSOE == 0)) {
        opserr << "WARNING Linear::solveCurrentStep() -";
        opserr << "setLinks() has not been called.\n";
        return -5;
    }

    if (factorOnce != 2) {
        if (theIncIntegrator->formTangent(incrTangent) < 0) {
            opserr << "WARNING Linear::solveCurrentStep() -";
            opserr << "the Integrator failed in formTangent()\n";
            return -1;
        }
        if (factorOnce == 1)
            factorOnce = 2;
    }


    if (theIncIntegrator->formUnbalance() < 0) {
        opserr << "WARNING Linear::solveCurrentStep() -";
        opserr << "the Integrator failed in formUnbalance()\n";
        return -2;
    }

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

    const Vector &deltaU = theSOE->getX();

    if (theIncIntegrator->update(deltaU) < 0) {
        opserr << "WARNING Linear::solveCurrentStep() -";
        opserr << "the Integrator failed in update()\n";
        return -4;
    }

    return 0;
}

int 
InitialInterpolatedLineSearch::newStep(LinearSOE &theSOE)
{
  const Vector &dU = theSOE.getX();

  if (x == 0)
    x = new Vector(dU);

  if (x->Size() != dU.Size()) {
    delete x;
    x = new Vector(dU);
  }

  return 0;
}

int 
SecantLineSearch::newStep(LinearSOE &theSOE)
{
  const Vector &dU = theSOE.getX();

  if (x == 0)
    x = new Vector(dU);

  if (x->Size() != dU.Size()) {
    delete x;
    x = new Vector(dU);
  }

  return 0;
}

int
ArcLength1::update(const Vector &dU)
{
    AnalysisModel *theModel = this->getAnalysisModel();
    LinearSOE *theLinSOE = this->getLinearSOE();    
    if (theModel == 0 || theLinSOE == 0) {
	opserr << "WARNING ArcLength1::update() ";
	opserr << "No AnalysisModel or LinearSOE has been set\n";
	return -1;
    }

    (*deltaUbar) = dU; // have to do this as the SOE is gonna change

    // determine dUhat    
    theLinSOE->setB(*phat);
    theLinSOE->solve();
    (*deltaUhat) = theLinSOE->getX();    

    // determine delta lambda(i)
    double a = (*deltaUstep)^(*deltaUbar);
    double b = (*deltaUstep)^(*deltaUhat) + alpha2*deltaLambdaStep;
    if (b == 0) {
      opserr << "ArcLength1::update() - zero denominator,";
      opserr << " alpha was set to 0.0 and zero reference load\n";
      return -1;
    }
    double dLambda = -a/b;

    // determine delta U(i)
    (*deltaU) = (*deltaUbar);    
    deltaU->addVector(1.0, *deltaUhat,dLambda);
    
    // update dU and dlambda
    (*deltaUstep) += *deltaU;
    deltaLambdaStep += dLambda;
    currentLambda += dLambda;

    // update the model
    theModel->incrDisp(*deltaU);    
    theModel->applyLoadDomain(currentLambda);    
    theModel->updateDomain();
    
    // set the X soln in linearSOE to be deltaU for convergence Test
    theLinSOE->setX(*deltaU);

    return 0;
}

int
ArcLengthw::update(const Vector &dU)
{
    ofstream factor;
    factor.open("FS.dat",ios::app);

    factor<<"insideupdate"<<endln;
    //factor>>dU;
    factor<<"insideupdate1"<<endln;

    AnalysisModel *theModel = this->getAnalysisModel();
    LinearSOE *theLinSOE = this->getLinearSOE();    
    if (theModel == 0 || theLinSOE == 0) {
	opserr << "WARNING ArcLengthw::update() ";
	opserr << "No AnalysisModel or LinearSOE has been set\n";
	return -1;
    }

    (*deltaUbar) = dU; // have to do this as the SOE is gonna change

    // determine dUhat    
    theLinSOE->setB(*phat);
    theLinSOE->solve();
    (*deltaUhat) = theLinSOE->getX();    

    double dLambda = -((*phat)^(*deltaUbar))/((*phat)^(*deltaUhat));

    (*deltaU) = (*deltaUbar);    
    deltaU->addVector(1.0, *deltaUhat,dLambda);
    
    // update dU and dlambda
    (*deltaUstep) += *deltaU;
    deltaLambdaStep += dLambda;
    currentLambda += dLambda;

    // update the model
    theModel->incrDisp(*deltaU);    
    theModel->applyLoadDomain(currentLambda);    
    theModel->updateDomain();
    
    // set the X soln in linearSOE to be deltaU for convergence Test
    theLinSOE->setX(*deltaU);

    return 0;
}

int HHTHSFixedNumIter::commit(void)
{
    AnalysisModel *theModel = this->getAnalysisModel();
    if (theModel == 0)  {
        opserr << "WARNING HHTHSFixedNumIter::commit() - no AnalysisModel set\n";
        return -1;
    }
    
    LinearSOE *theSOE = this->getLinearSOE();
    if (theSOE == 0)  {
        opserr << "WARNING HHTHSFixedNumIter::commit() - no LinearSOE set\n";
        return -2;
    }
    
    if (this->formTangent(statusFlag) < 0)  {
        opserr << "WARNING HHTHSFixedNumIter::commit() - "
            << "the Integrator failed in formTangent()\n";	
        return -3;
    }
    
    if (theSOE->solve() < 0)  {
        opserr << "WARNING HHTHSFixedNumIter::commit() - "
            << "the LinearSysOfEqn failed in solve()\n";	
        return -4;
    }
    const Vector &deltaU = theSOE->getX();
    
    //  determine the response at t+deltaT
    U->addVector(1.0, deltaU, c1);
    
    Udot->addVector(1.0, deltaU, c2);
    
    Udotdot->addVector(1.0, deltaU, c3);
    
    // update the response at the DOFs
    theModel->setResponse(*U,*Udot,*Udotdot);
    
    // set the time to be t+deltaT
    double time = theModel->getCurrentDomainTime();
    time += (1.0-alphaF)*deltaT;
    theModel->setCurrentDomainTime(time);
    
    return theModel->commitDomain();
}

int 
AlgorithmIncrements::record(int cTag, double timeStamp)
{
  LinearSOE *theSOE = theAlgo->getLinearSOEptr();
  if (theSOE == 0)
    return 0;


  const Vector &B = theSOE->getB();
  const Vector &X = theSOE->getX();

  if (fileName != 0) {
    if (cTag == 0) {

      if (theFile)
	theFile.close();

      numRecord = 0;

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

  if (displayRecord == true)
    return this->plotData(X, B);

  return 0;
}

int HHTHSFixedNumIter::domainChanged()
{
    AnalysisModel *myModel = this->getAnalysisModel();
    LinearSOE *theLinSOE = this->getLinearSOE();
    const Vector &x = theLinSOE->getX();
    int size = x.Size();
    
    // if damping factors exist set them in the ele & node of the domain
    if (alphaM != 0.0 || betaK != 0.0 || betaKi != 0.0 || betaKc != 0.0)
        myModel->setRayleighDampingFactors(alphaM, betaK, betaKi, betaKc);
    
    // create the new Vector objects
    if (Ut == 0 || Ut->Size() != size)  {
        
        // delete the old
        if (Ut != 0)
            delete Ut;
        if (Utdot != 0)
            delete Utdot;
        if (Utdotdot != 0)
            delete Utdotdot;
        if (U != 0)
            delete U;
        if (Udot != 0)
            delete Udot;
        if (Udotdot != 0)
            delete Udotdot;
        if (Ualpha != 0)
            delete Ualpha;
        if (Ualphadot != 0)
            delete Ualphadot;
        if (Ualphadotdot != 0)
            delete Ualphadotdot;
        if (Utm1 != 0)
            delete Utm1;
        if (Utm2 != 0)
            delete Utm2;
        if (scaledDeltaU != 0)
            delete scaledDeltaU;
        
        // create the new
        Ut = new Vector(size);
        Utdot = new Vector(size);
        Utdotdot = new Vector(size);
        U = new Vector(size);
        Udot = new Vector(size);
        Udotdot = new Vector(size);
        Ualpha = new Vector(size);
        Ualphadot = new Vector(size);
        Ualphadotdot = new Vector(size);
        Utm1 = new Vector(size);
        Utm2 = new Vector(size);
        scaledDeltaU = new Vector(size);
        
        // check we obtained the new
        if (Ut == 0 || Ut->Size() != size ||
            Utdot == 0 || Utdot->Size() != size ||
            Utdotdot == 0 || Utdotdot->Size() != size ||
            U == 0 || U->Size() != size ||
            Udot == 0 || Udot->Size() != size ||
            Udotdot == 0 || Udotdot->Size() != size ||
            Ualpha == 0 || Ualpha->Size() != size ||
            Ualphadot == 0 || Ualphadot->Size() != size ||
            Ualphadotdot == 0 || Ualphadotdot->Size() != size ||
            Utm1 == 0 || Utm1->Size() != size ||
            Utm2 == 0 || Utm2->Size() != size ||
            scaledDeltaU == 0 || scaledDeltaU->Size() != size)  {
            
            opserr << "HHTHSFixedNumIter::domainChanged - ran out of memory\n";
            
            // delete the old
            if (Ut != 0)
                delete Ut;
            if (Utdot != 0)
                delete Utdot;
            if (Utdotdot != 0)
                delete Utdotdot;
            if (U != 0)
                delete U;
            if (Udot != 0)
                delete Udot;
            if (Udotdot != 0)
                delete Udotdot;
            if (Ualpha != 0)
                delete Ualpha;
            if (Ualphadot != 0)
                delete Ualphadot;
            if (Ualphadotdot != 0)
                delete Ualphadotdot;
            if (Utm1 != 0)
                delete Utm1;
            if (Utm2 != 0)
                delete Utm2;
            if (scaledDeltaU != 0)
                delete scaledDeltaU;
            
            Ut = 0; Utdot = 0; Utdotdot = 0;
            U = 0; Udot = 0; Udotdot = 0;
            Ualpha = 0; Ualphadot = 0; Ualphadotdot = 0;
            Utm1 = 0; Utm2 = 0; scaledDeltaU = 0;
//.........这里部分代码省略.........

int 
NewtonLineSearch::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 NewtonLineSearch::solveCurrentStep() - setLinks() has";
	opserr << " not been called - or no ConvergenceTest has been set\n";
	return -5;
    }	

    theLineSearch->newStep(*theSOE);

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

    if (theIntegrator->formUnbalance() < 0) {
      opserr << "WARNING NewtonLineSearch::solveCurrentStep() -";
      opserr << "the Integrator failed in formUnbalance()\n";	
      return -2;
    }	    

    int result = -1;
    do {

	//residual at this iteration before next solve 
	const Vector &Resid0 = theSOE->getB() ;
	
	//form the tangent
        if (theIntegrator->formTangent() < 0){
	    opserr << "WARNING NewtonLineSearch::solveCurrentStep() -";
	    opserr << "the Integrator failed in formTangent()\n";
	    return -1;
	}		    
	
	//solve 
	if (theSOE->solve() < 0) {
	    opserr << "WARNING NewtonLineSearch::solveCurrentStep() -";
	    opserr << "the LinearSysOfEqn failed in solve()\n";	
	    return -3;
	}	    


	//line search direction 
	const Vector &dx0 = theSOE->getX() ;

	//intial value of s
	double s0 = - (dx0 ^ Resid0) ; 

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

	if (theIntegrator->formUnbalance() < 0) {
	    opserr << "WARNING NewtonLineSearch::solveCurrentStep() -";
	    opserr << "the Integrator failed in formUnbalance()\n";	
	    return -2;
	}	

	// do a line search only if convergence criteria not met
	theOtherTest->start();
	result = theOtherTest->test();

	if (result < 1) {
	  //new residual 
	  const Vector &Resid = theSOE->getB() ;
	  
	  //new value of s 
	  double s = - ( dx0 ^ Resid ) ;
	  
	  if (theLineSearch != 0)
	    theLineSearch->search(s0, s, *theSOE, *theIntegrator);
	}

	this->record(0);
	  
	result = theTest->test();

    } while (result == -1);

    if (result == -2) {
      opserr << "NewtonLineSearch::solveCurrentStep() -";
      opserr << "the ConvergenceTest object failed in test()\n";
      return -3;
    }

    // note - if postive result we are returning what the convergence test returned
    // which should be the number of iterations
//.........这里部分代码省略.........

int 
NewtonRaphson::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 NewtonRaphson::solveCurrentStep() - setLinks() has";
	opserr << " not been called - or no ConvergenceTest has been set\n";
	return -5;
    }	

    if (theIntegrator->formUnbalance() < 0) {
      opserr << "WARNING NewtonRaphson::solveCurrentStep() -";
      opserr << "the Integrator failed in formUnbalance()\n";	
      return -2;
    }	    

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

    int result = -1;
    int count = 0;
    do {
      if (tangent == INITIAL_THEN_CURRENT_TANGENT) {
	if (count == 0) {
	  if (theIntegrator->formTangent(INITIAL_TANGENT) < 0){
	    opserr << "WARNING NewtonRaphson::solveCurrentStep() -";
	    opserr << "the Integrator failed in formTangent()\n";
	    return -1;
	  } 
	} else {
	  if (theIntegrator->formTangent(CURRENT_TANGENT) < 0){
	    opserr << "WARNING NewtonRaphson::solveCurrentStep() -";
	    opserr << "the Integrator failed in formTangent()\n";
	    return -1;
	  } 
	}	  
      } else {
	if (theIntegrator->formTangent(tangent) < 0){
	    opserr << "WARNING NewtonRaphson::solveCurrentStep() -";
	    opserr << "the Integrator failed in formTangent()\n";
	    return -1;
	}		    
      }

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

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

      if (theIntegrator->formUnbalance() < 0) {
	opserr << "WARNING NewtonRaphson::solveCurrentStep() -";
	opserr << "the Integrator failed in formUnbalance()\n";	
	return -2;
      }	

      result = theTest->test();
      this->record(count++);

    } while (result == -1);

    if (result == -2) {
      opserr << "NewtnRaphson::solveCurrentStep() -";
      opserr << "the ConvergenceTest object failed in test()\n";
      return -3;
    }

    // note - if postive result we are returning what the convergence test returned
    // which should be the number of iterations
    return result;
}