00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034 #include <NewmarkHybridSimulation.h>
00035 #include <FE_Element.h>
00036 #include <LinearSOE.h>
00037 #include <AnalysisModel.h>
00038 #include <Vector.h>
00039 #include <DOF_Group.h>
00040 #include <DOF_GrpIter.h>
00041 #include <AnalysisModel.h>
00042 #include <Channel.h>
00043 #include <FEM_ObjectBroker.h>
00044 #include <ConvergenceTest.h>
00045
00046
00047 NewmarkHybridSimulation::NewmarkHybridSimulation()
00048 : TransientIntegrator(INTEGRATOR_TAGS_NewmarkHybridSimulation),
00049 gamma(0), beta(0),
00050 alphaM(0.0), betaK(0.0), betaKi(0.0), betaKc(0.0),
00051 theTest(0), rFact(1.0), c1(0.0), c2(0.0), c3(0.0),
00052 Ut(0), Utdot(0), Utdotdot(0), U(0), Udot(0), Udotdot(0)
00053 {
00054
00055 }
00056
00057 NewmarkHybridSimulation::NewmarkHybridSimulation(double _gamma, double _beta,
00058 ConvergenceTest &theT)
00059 : TransientIntegrator(INTEGRATOR_TAGS_NewmarkHybridSimulation),
00060 gamma(_gamma), beta(_beta),
00061 alphaM(0.0), betaK(0.0), betaKi(0.0), betaKc(0.0),
00062 theTest(&theT), rFact(1.0), c1(0.0), c2(0.0), c3(0.0),
00063 Ut(0), Utdot(0), Utdotdot(0), U(0), Udot(0), Udotdot(0)
00064 {
00065
00066 }
00067
00068 NewmarkHybridSimulation::NewmarkHybridSimulation(double _gamma, double _beta,
00069 ConvergenceTest &theT,
00070 double _alphaM, double _betaK, double _betaKi, double _betaKc)
00071 : TransientIntegrator(INTEGRATOR_TAGS_NewmarkHybridSimulation),
00072 gamma(_gamma), beta(_beta),
00073 alphaM(_alphaM), betaK(_betaK), betaKi(_betaKi), betaKc(_betaKc),
00074 theTest(&theT), rFact(1.0), c1(0.0), c2(0.0), c3(0.0),
00075 Ut(0), Utdot(0), Utdotdot(0), U(0), Udot(0), Udotdot(0)
00076 {
00077
00078 }
00079
00080 NewmarkHybridSimulation::~NewmarkHybridSimulation()
00081 {
00082
00083 if (Ut != 0)
00084 delete Ut;
00085 if (Utdot != 0)
00086 delete Utdot;
00087 if (Utdotdot != 0)
00088 delete Utdotdot;
00089 if (U != 0)
00090 delete U;
00091 if (Udot != 0)
00092 delete Udot;
00093 if (Udotdot != 0)
00094 delete Udotdot;
00095 }
00096
00097 int NewmarkHybridSimulation::newStep(double deltaT)
00098 {
00099 if (beta == 0 || gamma == 0) {
00100 opserr << "NewmarkHybridSimulation::newStep() - error in variable\n";
00101 opserr << "gamma = " << gamma << " beta = " << beta << endln;
00102 return -1;
00103 }
00104
00105 if (deltaT <= 0.0) {
00106 opserr << "NewmarkHybridSimulation::newStep() - error in variable\n";
00107 opserr << "dT = " << deltaT << endln;
00108 return -2;
00109 }
00110
00111
00112 AnalysisModel *theModel = this->getAnalysisModelPtr();
00113
00114
00115 c1 = 1.0;
00116 c2 = gamma/(beta*deltaT);
00117 c3 = 1.0/(beta*deltaT*deltaT);
00118
00119 if (U == 0) {
00120 opserr << "NewmarkHybridSimulation::newStep() - domainChange() failed or hasn't been called\n";
00121 return -3;
00122 }
00123
00124
00125 (*Ut) = *U;
00126 (*Utdot) = *Udot;
00127 (*Utdotdot) = *Udotdot;
00128
00129
00130 double time = theModel->getCurrentDomainTime();
00131 time += deltaT;
00132 theModel->applyLoadDomain(time);
00133
00134
00135 double a1 = (1.0 - gamma/beta);
00136 double a2 = (deltaT)*(1.0 - 0.5*gamma/beta);
00137 Udot->addVector(a1, *Utdotdot, a2);
00138
00139 double a3 = -1.0/(beta*deltaT);
00140 double a4 = 1.0 - 0.5/beta;
00141 Udotdot->addVector(a4, *Utdot, a3);
00142
00143
00144 theModel->setVel(*Udot);
00145 theModel->setAccel(*Udotdot);
00146
00147 return 0;
00148 }
00149
00150 int NewmarkHybridSimulation::revertToLastStep()
00151 {
00152
00153 if (U != 0) {
00154 (*U) = *Ut;
00155 (*Udot) = *Utdot;
00156 (*Udotdot) = *Utdotdot;
00157 }
00158
00159 return 0;
00160 }
00161
00162 int NewmarkHybridSimulation::formEleTangent(FE_Element *theEle)
00163 {
00164 theEle->zeroTangent();
00165
00166 if (statusFlag == CURRENT_TANGENT) {
00167 theEle->addKtToTang(c1);
00168 theEle->addCtoTang(c2);
00169 theEle->addMtoTang(c3);
00170 } else if (statusFlag == INITIAL_TANGENT) {
00171 theEle->addKiToTang(c1);
00172 theEle->addCtoTang(c2);
00173 theEle->addMtoTang(c3);
00174 }
00175
00176 return 0;
00177 }
00178
00179 int NewmarkHybridSimulation::formNodTangent(DOF_Group *theDof)
00180 {
00181 theDof->zeroTangent();
00182
00183 theDof->addCtoTang(c2);
00184 theDof->addMtoTang(c3);
00185
00186 return 0;
00187 }
00188
00189 int NewmarkHybridSimulation::domainChanged()
00190 {
00191 AnalysisModel *myModel = this->getAnalysisModelPtr();
00192 LinearSOE *theLinSOE = this->getLinearSOEPtr();
00193 const Vector &x = theLinSOE->getX();
00194 int size = x.Size();
00195
00196
00197 if (alphaM != 0.0 || betaK != 0.0 || betaKi != 0.0 || betaKc != 0.0)
00198 myModel->setRayleighDampingFactors(alphaM, betaK, betaKi, betaKc);
00199
00200
00201 if (Ut == 0 || Ut->Size() != size) {
00202
00203
00204 if (Ut != 0)
00205 delete Ut;
00206 if (Utdot != 0)
00207 delete Utdot;
00208 if (Utdotdot != 0)
00209 delete Utdotdot;
00210 if (U != 0)
00211 delete U;
00212 if (Udot != 0)
00213 delete Udot;
00214 if (Udotdot != 0)
00215 delete Udotdot;
00216
00217
00218 Ut = new Vector(size);
00219 Utdot = new Vector(size);
00220 Utdotdot = new Vector(size);
00221 U = new Vector(size);
00222 Udot = new Vector(size);
00223 Udotdot = new Vector(size);
00224
00225
00226 if (Ut == 0 || Ut->Size() != size ||
00227 Utdot == 0 || Utdot->Size() != size ||
00228 Utdotdot == 0 || Utdotdot->Size() != size ||
00229 U == 0 || U->Size() != size ||
00230 Udot == 0 || Udot->Size() != size ||
00231 Udotdot == 0 || Udotdot->Size() != size) {
00232
00233 opserr << "NewmarkHybridSimulation::domainChanged - ran out of memory\n";
00234
00235
00236 if (Ut != 0)
00237 delete Ut;
00238 if (Utdot != 0)
00239 delete Utdot;
00240 if (Utdotdot != 0)
00241 delete Utdotdot;
00242 if (U != 0)
00243 delete U;
00244 if (Udot != 0)
00245 delete Udot;
00246 if (Udotdot != 0)
00247 delete Udotdot;
00248
00249 Ut = 0; Utdot = 0; Utdotdot = 0;
00250 U = 0; Udot = 0; Udotdot = 0;
00251
00252 return -1;
00253 }
00254 }
00255
00256
00257
00258 DOF_GrpIter &theDOFs = myModel->getDOFs();
00259 DOF_Group *dofPtr;
00260
00261 while ((dofPtr = theDOFs()) != 0) {
00262 const ID &id = dofPtr->getID();
00263 int idSize = id.Size();
00264
00265 int i;
00266 const Vector &disp = dofPtr->getCommittedDisp();
00267 for (i=0; i < idSize; i++) {
00268 int loc = id(i);
00269 if (loc >= 0) {
00270 (*U)(loc) = disp(i);
00271 }
00272 }
00273
00274 const Vector &vel = dofPtr->getCommittedVel();
00275 for (i=0; i < idSize; i++) {
00276 int loc = id(i);
00277 if (loc >= 0) {
00278 (*Udot)(loc) = vel(i);
00279 }
00280 }
00281
00282 const Vector &accel = dofPtr->getCommittedAccel();
00283 for (i=0; i < idSize; i++) {
00284 int loc = id(i);
00285 if (loc >= 0) {
00286 (*Udotdot)(loc) = accel(i);
00287 }
00288 }
00289 }
00290
00291 return 0;
00292 }
00293
00294 int NewmarkHybridSimulation::update(const Vector &deltaU)
00295 {
00296 AnalysisModel *theModel = this->getAnalysisModelPtr();
00297 if (theModel == 0) {
00298 opserr << "WARNING NewmarkHybridSimulation::update() - no AnalysisModel set\n";
00299 return -1;
00300 }
00301
00302
00303 if (Ut == 0) {
00304 opserr << "WARNING NewmarkHybridSimulation::update() - domainChange() failed or not called\n";
00305 return -2;
00306 }
00307
00308
00309 if (deltaU.Size() != U->Size()) {
00310 opserr << "WARNING NewmarkHybridSimulation::update() - Vectors of incompatible size ";
00311 opserr << " expecting " << U->Size() << " obtained " << deltaU.Size() << endln;
00312 return -3;
00313 }
00314
00315
00316 rFact = 1.0/(theTest->getMaxNumTests() - theTest->getNumTests() + 1.0);
00317
00318
00319 (*U) += rFact*deltaU;
00320
00321 Udot->addVector(1.0, deltaU, rFact*c2);
00322
00323 Udotdot->addVector(1.0, deltaU, rFact*c3);
00324
00325
00326 theModel->setResponse(*U,*Udot,*Udotdot);
00327 if (theModel->updateDomain() < 0) {
00328 opserr << "NewmarkHybridSimulation::update() - failed to update the domain\n";
00329 return -4;
00330 }
00331
00332 return 0;
00333 }
00334
00335 int NewmarkHybridSimulation::sendSelf(int cTag, Channel &theChannel)
00336 {
00337 Vector data(8);
00338 data(0) = gamma;
00339 data(1) = beta;
00340 data(2) = theTest->getClassTag();
00341 data(3) = theTest->getDbTag();
00342 data(4) = alphaM;
00343 data(5) = betaK;
00344 data(6) = betaKi;
00345 data(7) = betaKc;
00346
00347 if (theChannel.sendVector(this->getDbTag(), cTag, data) < 0) {
00348 opserr << "WARNING NewmarkHybridSimulation::sendSelf() - could not send data\n";
00349 return -1;
00350 }
00351
00352 if (theTest->sendSelf(cTag, theChannel) < 0) {
00353 opserr << "WARNING NewmarkHybridSimulation::sendSelf() - failed to send CTest object\n";
00354 return -1;
00355 }
00356
00357 return 0;
00358 }
00359
00360 int NewmarkHybridSimulation::recvSelf(int cTag, Channel &theChannel, FEM_ObjectBroker &theBroker)
00361 {
00362 Vector data(8);
00363 if (theChannel.recvVector(this->getDbTag(), cTag, data) < 0) {
00364 opserr << "WARNING NewmarkHybridSimulation::recvSelf() - could not receive data\n";
00365 gamma = 0.5; beta = 0.25;
00366 return -1;
00367 }
00368
00369 gamma = data(0);
00370 beta = data(1);
00371 int ctType = int(data(2));
00372 int ctDb = int(data(3));
00373 alphaM = data(4);
00374 betaK = data(5);
00375 betaKi = data(6);
00376 betaKc = data(7);
00377
00378 theTest = theBroker.getNewConvergenceTest(ctType);
00379 theTest->setDbTag(ctDb);
00380 if (theTest->recvSelf(cTag, theChannel, theBroker) < 0) {
00381 opserr << "WARNING NewmarkHybridSimulation::recvSelf() - failed to recv CTest object\n";
00382 return -1;
00383 }
00384
00385 return 0;
00386 }
00387
00388 void NewmarkHybridSimulation::Print(OPS_Stream &s, int flag)
00389 {
00390 AnalysisModel *theModel = this->getAnalysisModelPtr();
00391 if (theModel != 0) {
00392 double currentTime = theModel->getCurrentDomainTime();
00393 s << "\t NewmarkHybridSimulation - currentTime: " << currentTime;
00394 s << " gamma: " << gamma << " beta: " << beta << endln;
00395 s << " c1: " << c1 << " c2: " << c2 << " c3: " << c3 << endln;
00396 s << " Rayleigh Damping - alphaM: " << alphaM;
00397 s << " betaK: " << betaK << " betaKi: " << betaKi << endln;
00398 } else
00399 s << "\t NewmarkHybridSimulation - no associated AnalysisModel\n";
00400 }