Pinching4Material.cppGo to the documentation of this file.00001 /* ****************************************************************** ** 00002 ** OpenSees - Open System for Earthquake Engineering Simulation ** 00003 ** Pacific Earthquake Engineering Research Center ** 00004 ** ** 00005 ** ** 00006 ** (C) Copyright 1999, The Regents of the University of California ** 00007 ** All Rights Reserved. ** 00008 ** ** 00009 ** Commercial use of this program without express permission of the ** 00010 ** University of California, Berkeley, is strictly prohibited. See ** 00011 ** file 'COPYRIGHT' in main directory for information on usage and ** 00012 ** redistribution, and for a DISCLAIMER OF ALL WARRANTIES. ** 00013 ** ** 00014 ** Developed by: ** 00015 ** Frank McKenna (fmckenna@ce.berkeley.edu) ** 00016 ** Gregory L. Fenves (fenves@ce.berkeley.edu) ** 00017 ** Filip C. Filippou (filippou@ce.berkeley.edu) ** 00018 ** ** 00019 ** ****************************************************************** */ 00020 00021 // $Revision: 1.2 $ 00022 // $Date: 2004/10/06 19:21:12 $ 00023 // $Source: /usr/local/cvs/OpenSees/SRC/material/uniaxial/Pinching4Material.cpp,v $ 00024 00025 00026 // Written: NM (nmitra@u.washington.edu) 00027 // Created: January 2002 00028 // Updates: September 2004 00029 // 00030 // Description: This file contains the class implementation for 00031 // Pinching material which is defined by 4 points on the positive and 00032 // negative envelopes and a bunch of damage parameters. The material accounts for 00033 // 3 types of damage rules : Strength degradation, Stiffness degradation, 00034 // unloading stiffness degradation. 00035 // Updates: damage calculations and several bug fixes 00036 00037 00038 #include <Pinching4Material.h> 00039 #include <OPS_Globals.h> 00040 #include <stdlib.h> 00041 #include <math.h> 00042 #include <float.h> 00043 #include <OPS_Stream.h> 00044 00045 #include <string.h> 00046 00047 00048 Pinching4Material::Pinching4Material(int tag, 00049 double f1p, double d1p, double f2p, double d2p, 00050 double f3p, double d3p, double f4p, double d4p, 00051 double f1n, double d1n, double f2n, double d2n, 00052 double f3n, double d3n, double f4n, double d4n, 00053 double mdp, double mfp, double msp, 00054 double mdn, double mfn, double msn, 00055 double gk1, double gk2, double gk3, double gk4, double gklim, 00056 double gd1, double gd2, double gd3, double gd4, double gdlim, 00057 double gf1, double gf2, double gf3, double gf4, double gflim, double ge, int dc): 00058 UniaxialMaterial(tag, MAT_TAG_Pinching4), 00059 stress1p(f1p), strain1p(d1p), stress2p(f2p), strain2p(d2p), 00060 stress3p(f3p), strain3p(d3p), stress4p(f4p), strain4p(d4p), 00061 stress1n(f1n), strain1n(d1n), stress2n(f2n), strain2n(d2n), 00062 stress3n(f3n), strain3n(d3n), stress4n(f4n), strain4n(d4n), 00063 envlpPosStress(6), envlpPosStrain(6), envlpNegStress(6), envlpNegStrain(6), tagMat(tag), 00064 gammaK1(gk1), gammaK2(gk2), gammaK3(gk3), gammaK4(gk4), gammaKLimit(gklim), 00065 gammaD1(gd1), gammaD2(gd2), gammaD3(gd3), gammaD4(gd4), gammaDLimit(gdlim), 00066 gammaF1(gf1), gammaF2(gf2), gammaF3(gf3), gammaF4(gf4), gammaFLimit(gflim), 00067 gammaE(ge), TnCycle(0.0), CnCycle(0.0), DmgCyc(dc), 00068 rDispP(mdp), rForceP(mfp), uForceP(msp), rDispN(mdn), rForceN(mfn), uForceN(msn), 00069 state3Stress(4), state3Strain(4), state4Stress(4), state4Strain(4), 00070 envlpPosDamgdStress(6), envlpNegDamgdStress(6) 00071 { 00072 bool error = false; 00073 00074 // Positive backbone parameters 00075 if (strain1p <= 0.0) 00076 error = true; 00077 00078 if (strain2p <= 0.0) 00079 error = true; 00080 00081 if (strain3p <= 0.0) 00082 error = true; 00083 00084 if (strain4p <= 0.0) 00085 error = true; 00086 00087 // Negative backbone parameters 00088 if (strain1n >= 0.0) 00089 error = true; 00090 00091 if (strain2n >= 0.0) 00092 error = true; 00093 00094 if (strain3n >= 0.0) 00095 error = true; 00096 00097 if (strain4n >= 0.0) 00098 error = true; 00099 00100 if (error){ 00101 opserr << "ERROR: -- input backbone is not unique (one-to-one) , Pinching4Material::Pinching4Material" << "\a"; 00102 } 00103 00104 envlpPosStress.Zero(); envlpPosStrain.Zero(); envlpNegStress.Zero(); envlpNegStrain.Zero(); 00105 energyCapacity = 0.0; kunload = 0.0; elasticStrainEnergy = 0.0; 00106 00107 // set envelope slopes 00108 this->SetEnvelope(); 00109 00110 envlpPosDamgdStress = envlpPosStress; envlpNegDamgdStress = envlpNegStress; 00111 state3Stress.Zero(); state3Strain.Zero(); state4Stress.Zero(); state4Strain.Zero(); 00112 // Initialize history variables 00113 this->revertToStart(); 00114 this->revertToLastCommit(); 00115 } 00116 00117 Pinching4Material::Pinching4Material(int tag, 00118 double f1p, double d1p, double f2p, double d2p, 00119 double f3p, double d3p, double f4p, double d4p, 00120 double mdp, double mfp, double msp, 00121 double gk1, double gk2, double gk3, double gk4, double gklim, 00122 double gd1, double gd2, double gd3, double gd4, double gdlim, 00123 double gf1, double gf2, double gf3, double gf4, double gflim, double ge, int dc): 00124 UniaxialMaterial(tag, MAT_TAG_Pinching4), 00125 stress1p(f1p), strain1p(d1p), stress2p(f2p), strain2p(d2p), 00126 stress3p(f3p), strain3p(d3p), stress4p(f4p), strain4p(d4p), 00127 envlpPosStress(6), envlpPosStrain(6), envlpNegStress(6), envlpNegStrain(6), tagMat(tag), 00128 gammaK1(gk1), gammaK2(gk2), gammaK3(gk3), gammaK4(gk4), gammaKLimit(gklim), 00129 gammaD1(gd1), gammaD2(gd2), gammaD3(gd3), gammaD4(gd4), gammaDLimit(gdlim), 00130 gammaF1(gf1), gammaF2(gf2), gammaF3(gf3), gammaF4(gf4), gammaFLimit(gflim), 00131 gammaE(ge), TnCycle(0.0), CnCycle(0.0), DmgCyc(dc), 00132 rDispP(mdp), rForceP(mfp), uForceP(msp), 00133 state3Stress(4), state3Strain(4), state4Stress(4), state4Strain(4), 00134 envlpPosDamgdStress(6), envlpNegDamgdStress(6) 00135 00136 { 00137 bool error = false; 00138 00139 // Positive backbone parameters 00140 if (strain1p <= 0.0) 00141 error = true; 00142 00143 if (strain2p <= 0.0) 00144 error = true; 00145 00146 if (strain3p <= 0.0) 00147 error = true; 00148 00149 if (strain4p <= 0.0) 00150 error = true; 00151 00152 if (error){ 00153 opserr << "ERROR: -- input backbone is not unique (one-to-one) , Pinching4Material::Pinching4Material" << "\a"; 00154 } 00155 00156 strain1n = -strain1p; stress1n = -stress1p; strain2n = -strain2p; stress2n = -stress2p; 00157 strain3n = -strain3p; stress3n = -stress3p; strain4n = -strain4p; stress4n = -stress4p; 00158 rDispN = rDispP; rForceN = rForceP; uForceN = uForceP; 00159 00160 envlpPosStress.Zero(); envlpPosStrain.Zero(); envlpNegStress.Zero(); envlpNegStrain.Zero(); 00161 energyCapacity = 0.0; kunload = 0.0; elasticStrainEnergy = 0.0; 00162 state3Stress.Zero(); state3Strain.Zero(); state4Stress.Zero(); state4Strain.Zero(); 00163 // set envelope slopes 00164 this->SetEnvelope(); 00165 00166 envlpPosDamgdStress = envlpPosStress; envlpNegDamgdStress = envlpNegStress; 00167 00168 // Initialize history variables 00169 this->revertToStart(); 00170 this->revertToLastCommit(); 00171 } 00172 00173 00174 Pinching4Material::Pinching4Material(): 00175 UniaxialMaterial(0, MAT_TAG_Pinching4), 00176 stress1p(0.0), strain1p(0.0), stress2p(0.0), strain2p(0.0), 00177 stress3p(0.0), strain3p(0.0), stress4p(0.0), strain4p(0.0), 00178 stress1n(0.0), strain1n(0.0), stress2n(0.0), strain2n(0.0), 00179 stress3n(0.0), strain3n(0.0), stress4n(0.0), strain4n(0.0), 00180 gammaK1(0.0), gammaK2(0.0), gammaK3(0.0), gammaKLimit(0.0), 00181 gammaD1(0.0), gammaD2(0.0), gammaD3(0.0), gammaDLimit(0.0), 00182 gammaF1(0.0), gammaF2(0.0), gammaF3(0.0), gammaFLimit(0.0), gammaE(0.0), 00183 rDispP(0.0), rForceP(0.0), uForceP(0.0), rDispN(0.0), rForceN(0.0), uForceN(0.0) 00184 { 00185 00186 } 00187 00188 Pinching4Material::~Pinching4Material() 00189 { 00190 00191 } 00192 00193 int Pinching4Material::setTrialStrain(double strain, double CstrainRate) 00194 { 00195 00196 Tstate = Cstate; 00197 Tenergy = Cenergy; 00198 Tstrain = strain; 00199 lowTstateStrain = lowCstateStrain; 00200 hghTstateStrain = hghCstateStrain; 00201 lowTstateStress = lowCstateStress; 00202 hghTstateStress = hghCstateStress; 00203 TminStrainDmnd = CminStrainDmnd; 00204 TmaxStrainDmnd = CmaxStrainDmnd; 00205 TgammaF = CgammaF; 00206 TgammaK = CgammaK; 00207 TgammaD = CgammaD; 00208 00209 dstrain = Tstrain - Cstrain; 00210 if (dstrain<1e-12 && dstrain>-1e-12){ 00211 dstrain = 0.0; 00212 } 00213 00214 00215 // determine new state if there is a change in state 00216 getstate(Tstrain,dstrain); 00217 00218 switch (Tstate) 00219 { 00220 00221 case 0: 00222 Ttangent = envlpPosStress(0)/envlpPosStrain(0); 00223 Tstress = Ttangent*Tstrain; 00224 break; 00225 case 1: 00226 Tstress = posEnvlpStress(strain); 00227 Ttangent = posEnvlpTangent(strain); 00228 break; 00229 case 2: 00230 Ttangent = negEnvlpTangent(strain); 00231 Tstress = negEnvlpStress(strain); 00232 break; 00233 case 3: 00234 kunload = (hghTstateStrain<0.0) ? kElasticNegDamgd:kElasticPosDamgd; 00235 state3Strain(0) = lowTstateStrain; 00236 state3Strain(3) = hghTstateStrain; 00237 state3Stress(0) = lowTstateStress; 00238 state3Stress(3) = hghTstateStress; 00239 00240 getState3(state3Strain,state3Stress,kunload); 00241 Ttangent = Envlp3Tangent(state3Strain,state3Stress,strain); 00242 Tstress = Envlp3Stress(state3Strain,state3Stress,strain); 00243 00244 //Print(opserr,1); 00245 break; 00246 case 4: 00247 kunload = (lowTstateStrain<0.0) ? kElasticNegDamgd:kElasticPosDamgd; 00248 state4Strain(0) = lowTstateStrain; 00249 state4Strain(3) = hghTstateStrain; 00250 state4Stress(0) = lowTstateStress; 00251 state4Stress(3) = hghTstateStress; 00252 00253 getState4(state4Strain,state4Stress,kunload); 00254 Ttangent = Envlp4Tangent(state4Strain,state4Stress,strain); 00255 Tstress = Envlp4Stress(state4Strain,state4Stress,strain); 00256 break; 00257 } 00258 00259 double denergy = 0.5*(Tstress+Cstress)*dstrain; 00260 elasticStrainEnergy = (Tstrain>0.0) ? 0.5*Tstress/kElasticPosDamgd*Tstress:0.5*Tstress/kElasticNegDamgd*Tstress; 00261 00262 Tenergy = Cenergy + denergy; 00263 00264 updateDmg(Tstrain,dstrain); 00265 return 0; 00266 } 00267 00268 double Pinching4Material::getStrain(void) 00269 { 00270 return Tstrain; 00271 } 00272 00273 double Pinching4Material::getStress(void) 00274 { 00275 return Tstress; 00276 } 00277 00278 double Pinching4Material::getTangent(void) 00279 { 00280 return Ttangent; 00281 } 00282 00283 double Pinching4Material::getInitialTangent(void) 00284 { 00285 return envlpPosStress(0)/envlpPosStrain(0); 00286 } 00287 00288 int Pinching4Material::commitState(void) 00289 { 00290 Cstate = Tstate; 00291 00292 if (dstrain>1e-12||dstrain<-(1e-12)) { 00293 CstrainRate = dstrain;} 00294 else { 00295 CstrainRate = TstrainRate;} 00296 00297 lowCstateStrain = lowTstateStrain; 00298 lowCstateStress = lowTstateStress; 00299 hghCstateStrain = hghTstateStrain; 00300 hghCstateStress = hghTstateStress; 00301 CminStrainDmnd = TminStrainDmnd; 00302 CmaxStrainDmnd = TmaxStrainDmnd; 00303 Cenergy = Tenergy; 00304 00305 Cstress = Tstress; 00306 Cstrain = Tstrain; 00307 00308 CgammaK = TgammaK; 00309 CgammaD = TgammaD; 00310 CgammaF = TgammaF; 00311 00312 // define adjusted strength and stiffness parameters 00313 kElasticPosDamgd = kElasticPos*(1 - gammaKUsed); 00314 kElasticNegDamgd = kElasticNeg*(1 - gammaKUsed); 00315 00316 uMaxDamgd = TmaxStrainDmnd*(1 + CgammaD); 00317 uMinDamgd = TminStrainDmnd*(1 + CgammaD); 00318 00319 envlpPosDamgdStress = envlpPosStress*(1-gammaFUsed); 00320 envlpNegDamgdStress = envlpNegStress*(1-gammaFUsed); 00321 00322 CnCycle = TnCycle; // number of cycles of loading 00323 00324 return 0; 00325 } 00326 00327 int Pinching4Material::revertToLastCommit(void) 00328 { 00329 00330 Tstate = Cstate; 00331 00332 TstrainRate = CstrainRate; 00333 00334 lowTstateStrain = lowCstateStrain; 00335 lowTstateStress = lowCstateStress; 00336 hghTstateStrain = hghCstateStrain; 00337 hghTstateStress = hghCstateStress; 00338 TminStrainDmnd = CminStrainDmnd; 00339 TmaxStrainDmnd = CmaxStrainDmnd; 00340 Tenergy = Cenergy; 00341 00342 Tstrain = Cstrain; Tstress = Cstress; 00343 00344 TgammaD = CgammaD; 00345 TgammaK = CgammaK; 00346 TgammaF = CgammaF; 00347 00348 TnCycle = CnCycle; 00349 00350 return 0; 00351 } 00352 00353 int Pinching4Material::revertToStart(void) 00354 { 00355 Cstate = 0; 00356 Cstrain = 0.0; 00357 Cstress = 0.0; 00358 CstrainRate = 0.0; 00359 lowCstateStrain = envlpNegStrain(0); 00360 lowCstateStress = envlpNegStress(0); 00361 hghCstateStrain = envlpPosStrain(0); 00362 hghCstateStress = envlpPosStress(0); 00363 CminStrainDmnd = envlpNegStrain(1); 00364 CmaxStrainDmnd = envlpPosStrain(1); 00365 Cenergy = 0.0; 00366 CgammaK = 0.0; 00367 CgammaD = 0.0; 00368 CgammaF = 0.0; 00369 CnCycle = 0.0; 00370 00371 Ttangent = envlpPosStress(0)/envlpPosStrain(0); 00372 dstrain = 0.0; 00373 gammaKUsed = 0.0; 00374 gammaFUsed = 0.0; 00375 00376 kElasticPosDamgd = kElasticPos; 00377 kElasticNegDamgd = kElasticNeg; 00378 uMaxDamgd = CmaxStrainDmnd; 00379 uMinDamgd = CminStrainDmnd; 00380 00381 return 0; 00382 } 00383 00384 UniaxialMaterial* Pinching4Material::getCopy(void) 00385 { 00386 Pinching4Material *theCopy = new Pinching4Material (this->getTag(), 00387 stress1p, strain1p, stress2p, strain2p, stress3p, strain3p, stress4p, strain4p, 00388 stress1n, strain1n, stress2n, strain2n, stress3n, strain3n, stress4n, strain4n, 00389 rDispP, rForceP, uForceP, rDispN, rForceN, uForceN,gammaK1,gammaK2,gammaK3,gammaK4, 00390 gammaKLimit,gammaD1,gammaD2,gammaD3,gammaD4,gammaDLimit,gammaF1,gammaF2,gammaF3,gammaF4, 00391 gammaFLimit,gammaE,DmgCyc); 00392 00393 theCopy->rDispN = rDispN; 00394 theCopy->rDispP = rDispP; 00395 theCopy->rForceN = rForceN; 00396 theCopy->rForceP = rForceP; 00397 theCopy->uForceN = uForceN; 00398 theCopy->uForceP = uForceP; 00399 00400 // Trial state variables 00401 theCopy->Tstress = Tstress; 00402 theCopy->Tstrain = Tstrain; 00403 theCopy->Ttangent = Ttangent; 00404 00405 // Coverged material history parameters 00406 theCopy->Cstate = Cstate; 00407 theCopy->Cstrain = Cstrain; 00408 theCopy->Cstress = Cstress; 00409 theCopy->CstrainRate = CstrainRate; 00410 00411 theCopy->lowCstateStrain = lowCstateStrain; 00412 theCopy->lowCstateStress = lowCstateStress; 00413 theCopy->hghCstateStrain = hghCstateStrain; 00414 theCopy->hghCstateStress = hghCstateStress; 00415 theCopy->CminStrainDmnd = CminStrainDmnd; 00416 theCopy->CmaxStrainDmnd = CmaxStrainDmnd; 00417 theCopy->Cenergy = Cenergy; 00418 theCopy->CgammaK = CgammaK; 00419 theCopy->CgammaD = CgammaD; 00420 theCopy->CgammaF = CgammaF; 00421 theCopy->CnCycle = CnCycle; 00422 theCopy->gammaKUsed = gammaKUsed; 00423 theCopy->gammaFUsed = gammaFUsed; 00424 theCopy->DmgCyc = DmgCyc; 00425 00426 // trial material history parameters 00427 theCopy->Tstate = Tstate; 00428 theCopy->dstrain = dstrain; 00429 theCopy->lowTstateStrain = lowTstateStrain; 00430 theCopy->lowTstateStress = lowTstateStress; 00431 theCopy->hghTstateStrain = hghTstateStrain; 00432 theCopy->hghTstateStress = hghTstateStress; 00433 theCopy->TminStrainDmnd = TminStrainDmnd; 00434 theCopy->TmaxStrainDmnd = TmaxStrainDmnd; 00435 theCopy->Tenergy = Tenergy; 00436 theCopy->TgammaK = TgammaK; 00437 theCopy->TgammaD = TgammaD; 00438 theCopy->TgammaF = TgammaF; 00439 theCopy->TnCycle = TnCycle; 00440 00441 // Strength and stiffness parameters 00442 theCopy->kElasticPos = kElasticPos; 00443 theCopy->kElasticNeg = kElasticNeg; 00444 theCopy->kElasticPosDamgd = kElasticPosDamgd; 00445 theCopy->kElasticNegDamgd = kElasticNegDamgd; 00446 theCopy->uMaxDamgd = uMaxDamgd; 00447 theCopy->uMinDamgd = uMinDamgd; 00448 00449 for (int i = 0; i<6; i++) 00450 { 00451 theCopy->envlpPosStrain(i) = envlpPosStrain(i); 00452 theCopy->envlpPosStress(i) = envlpPosStress(i); 00453 theCopy->envlpNegStrain(i) = envlpNegStrain(i); 00454 theCopy->envlpNegStress(i) = envlpNegStress(i); 00455 theCopy->envlpNegDamgdStress(i) = envlpNegDamgdStress(i); 00456 theCopy->envlpPosDamgdStress(i) = envlpPosDamgdStress(i); 00457 } 00458 00459 for (int j = 0; j<4; j++) 00460 { 00461 theCopy->state3Strain(j) = state3Strain(j); 00462 theCopy->state3Stress(j) = state3Stress(j); 00463 theCopy->state4Strain(j) = state4Strain(j); 00464 theCopy->state4Stress(j) = state4Stress(j); 00465 } 00466 00467 theCopy->energyCapacity = energyCapacity; 00468 theCopy->kunload = kunload; 00469 theCopy->elasticStrainEnergy = elasticStrainEnergy; 00470 00471 return theCopy; 00472 } 00473 00474 int Pinching4Material::sendSelf(int commitTag, Channel &theChannel) 00475 { 00476 return -1; 00477 } 00478 00479 int Pinching4Material::recvSelf(int commitTag, Channel &theChannel, 00480 FEM_ObjectBroker & theBroker) 00481 { 00482 return -1; 00483 } 00484 00485 void Pinching4Material::Print(OPS_Stream &s, int flag) 00486 { 00487 s << "Pinching4Material, tag: " << this-> getTag() << endln; 00488 s << "strain: " << Tstrain << endln; 00489 s << "stress: " << Tstress << endln; 00490 s << "state: " << Tstate << endln; 00491 } 00492 00493 void Pinching4Material::SetEnvelope(void) 00494 { 00495 double kPos = stress1p/strain1p; 00496 double kNeg = stress1n/strain1n; 00497 double k = (kPos>kNeg) ? kPos:kNeg; 00498 double u = (strain1p>-strain1n) ? 1e-4*strain1p:-1e-4*strain1n; 00499 00500 envlpPosStrain(0) = u; 00501 envlpPosStress(0) = u*k; 00502 envlpNegStrain(0) = -u; 00503 envlpNegStress(0) = -u*k; 00504 00505 envlpPosStrain(1) = strain1p; 00506 envlpPosStrain(2) = strain2p; 00507 envlpPosStrain(3) = strain3p; 00508 envlpPosStrain(4) = strain4p; 00509 00510 envlpNegStrain(1) = strain1n; 00511 envlpNegStrain(2) = strain2n; 00512 envlpNegStrain(3) = strain3n; 00513 envlpNegStrain(4) = strain4n; 00514 00515 envlpPosStress(1) = stress1p; 00516 envlpPosStress(2) = stress2p; 00517 envlpPosStress(3) = stress3p; 00518 envlpPosStress(4) = stress4p; 00519 00520 envlpNegStress(1) = stress1n; 00521 envlpNegStress(2) = stress2n; 00522 envlpNegStress(3) = stress3n; 00523 envlpNegStress(4) = stress4n; 00524 00525 double k1 = (stress4p - stress3p)/(strain4p - strain3p); 00526 double k2 = (stress4n - stress3n)/(strain4n - strain3n); 00527 00528 00529 envlpPosStrain(5) = 1e+6*strain4p; 00530 envlpPosStress(5) = (k1>0.0)? stress4p+k1*(envlpPosStrain(5) - strain4p):stress4p*1.1; 00531 envlpNegStrain(5) = 1e+6*strain4n; 00532 envlpNegStress(5) = (k2>0.0)? stress4n+k2*(envlpNegStrain(5) - strain4n):stress4n*1.1; 00533 00534 // define crtical material properties 00535 kElasticPos = envlpPosStress(1)/envlpPosStrain(1); 00536 kElasticNeg = envlpNegStress(1)/envlpNegStrain(1); 00537 00538 double energypos = 0.5*envlpPosStrain(0)*envlpPosStress(0); 00539 00540 for (int jt = 0; jt<4; jt++){ 00541 energypos += 0.5*(envlpPosStress(jt) + envlpPosStress(jt+1))*(envlpPosStrain(jt+1)-envlpPosStrain(jt)); 00542 } 00543 00544 double energyneg = 0.5*envlpNegStrain(0)*envlpNegStress(0); 00545 00546 for (int jy = 0; jy<4; jy++){ 00547 energyneg += 0.5*(envlpNegStress(jy) + envlpNegStress(jy+1))*(envlpNegStrain(jy+1)-envlpNegStrain(jy)); 00548 } 00549 00550 double max_energy = (energypos>energyneg) ? energypos:energyneg; 00551 00552 energyCapacity = gammaE*max_energy; 00553 00554 00555 } 00556 00557 void Pinching4Material::getstate(double u,double du) 00558 { 00559 int cid = 0; 00560 int cis = 0; 00561 int newState = 0; 00562 if (du*CstrainRate<=0.0){ 00563 cid = 1; 00564 } 00565 if (u<lowTstateStrain || u>hghTstateStrain || cid) { 00566 if (Tstate == 0) { 00567 if (u>hghTstateStrain) { 00568 cis = 1; 00569 newState = 1; 00570 lowTstateStrain = envlpPosStrain(0); 00571 lowTstateStress = envlpPosStress(0); 00572 hghTstateStrain = envlpPosStrain(5); 00573 hghTstateStress = envlpPosStress(5); 00574 } 00575 else if (u<lowTstateStrain){ 00576 cis = 1; 00577 newState = 2; 00578 lowTstateStrain = envlpNegStrain(5); 00579 lowTstateStress = envlpNegStress(5); 00580 hghTstateStrain = envlpNegStrain(0); 00581 hghTstateStress = envlpNegStress(0); 00582 } 00583 } 00584 else if (Tstate==1 && du<0.0) { 00585 cis = 1; 00586 if (Cstrain>TmaxStrainDmnd) { 00587 TmaxStrainDmnd = u - du; 00588 } 00589 if (TmaxStrainDmnd<uMaxDamgd) { 00590 TmaxStrainDmnd = uMaxDamgd; 00591 } 00592 if (u<uMinDamgd) { 00593 newState = 2; 00594 gammaFUsed = CgammaF; 00595 for (int i=0; i<=5; i++) { 00596 envlpNegDamgdStress(i) = envlpNegStress(i)*(1.0-gammaFUsed); 00597 } 00598 lowTstateStrain = envlpNegStrain(5); 00599 lowTstateStress = envlpNegStress(5); 00600 hghTstateStrain = envlpNegStrain(0); 00601 hghTstateStress = envlpNegStress(0); 00602 } 00603 else { 00604 newState = 3; 00605 lowTstateStrain = uMinDamgd; 00606 gammaFUsed = CgammaF; 00607 for (int i=0; i<=5; i++) { 00608 envlpNegDamgdStress(i) = envlpNegStress(i)*(1.0-gammaFUsed); 00609 } 00610 lowTstateStress = negEnvlpStress(uMinDamgd); 00611 hghTstateStrain = Cstrain; 00612 hghTstateStress = Cstress; 00613 } 00614 gammaKUsed = CgammaK; 00615 kElasticPosDamgd = kElasticPos*(1.0-gammaKUsed); 00616 } 00617 else if (Tstate ==2 && du>0.0){ 00618 cis = 1; 00619 if (Cstrain<TminStrainDmnd) { 00620 TminStrainDmnd = Cstrain; 00621 } 00622 if (TminStrainDmnd>uMinDamgd) { 00623 TminStrainDmnd = uMinDamgd; 00624 } 00625 if (u>uMaxDamgd) { 00626 newState = 1; 00627 gammaFUsed = CgammaF; 00628 for (int i=0; i<=5; i++) { 00629 envlpPosDamgdStress(i) = envlpPosStress(i)*(1.0-gammaFUsed); 00630 } 00631 lowTstateStrain = envlpPosStrain(0); 00632 lowTstateStress = envlpPosStress(0); 00633 hghTstateStrain = envlpPosStrain(5); 00634 hghTstateStress = envlpPosStress(5); 00635 } 00636 else { 00637 newState = 4; 00638 lowTstateStrain = Cstrain; 00639 lowTstateStress = Cstress; 00640 hghTstateStrain = uMaxDamgd; 00641 gammaFUsed = CgammaF; 00642 for (int i=0; i<=5; i++) { 00643 envlpPosDamgdStress(i) = envlpPosStress(i)*(1.0-gammaFUsed); 00644 } 00645 hghTstateStress = posEnvlpStress(uMaxDamgd); 00646 } 00647 gammaKUsed = CgammaK; 00648 kElasticNegDamgd = kElasticNeg*(1.0-gammaKUsed); 00649 } 00650 else if (Tstate ==3) { 00651 if (u<lowTstateStrain){ 00652 cis = 1; 00653 newState = 2; 00654 lowTstateStrain = envlpNegStrain(5); 00655 hghTstateStrain = envlpNegStrain(0); 00656 lowTstateStress = envlpNegDamgdStress(5); 00657 hghTstateStress = envlpNegDamgdStress(0); 00658 } 00659 else if (u>uMaxDamgd && du>0.0) { 00660 cis = 1; 00661 newState = 1; 00662 lowTstateStrain = envlpPosStrain(0); 00663 lowTstateStress = envlpPosStress(0); 00664 hghTstateStrain = envlpPosStrain(5); 00665 hghTstateStress = envlpPosStress(5); 00666 } 00667 else if (du>0.0) { 00668 cis = 1; 00669 newState = 4; 00670 lowTstateStrain = Cstrain; 00671 lowTstateStress = Cstress; 00672 hghTstateStrain = uMaxDamgd; 00673 gammaFUsed = CgammaF; 00674 for (int i=0; i<=5; i++) { 00675 envlpPosDamgdStress(i) = envlpPosStress(i)*(1.0-gammaFUsed); 00676 } 00677 hghTstateStress = posEnvlpStress(uMaxDamgd); 00678 gammaKUsed = CgammaK; 00679 kElasticNegDamgd = kElasticNeg*(1.0-gammaKUsed); 00680 } 00681 } 00682 else if (Tstate == 4){ 00683 if (u>hghTstateStrain){ 00684 cis = 1; 00685 newState = 1; 00686 lowTstateStrain = envlpPosStrain(0); 00687 lowTstateStress = envlpPosDamgdStress(0); 00688 hghTstateStrain = envlpPosStrain(5); 00689 hghTstateStress = envlpPosDamgdStress(5); 00690 } 00691 else if (u<uMinDamgd && du <0.0) { 00692 cis = 1; 00693 newState = 2; 00694 lowTstateStrain = envlpNegStrain(5); 00695 lowTstateStress = envlpNegDamgdStress(5); 00696 hghTstateStrain = envlpNegStrain(0); 00697 hghTstateStress = envlpNegDamgdStress(0); 00698 } 00699 else if (du<0.0) { 00700 cis = 1; 00701 newState = 3; 00702 lowTstateStrain = uMinDamgd; 00703 gammaFUsed = CgammaF; 00704 for (int i=0; i<=5; i++) { 00705 envlpNegDamgdStress(i) = envlpNegStress(i)*(1.0-gammaFUsed); 00706 } 00707 lowTstateStress = negEnvlpStress(uMinDamgd); 00708 hghTstateStrain = Cstrain; 00709 hghTstateStress = Cstress; 00710 gammaKUsed = CgammaK; 00711 kElasticPosDamgd = kElasticPos*(1.0-gammaKUsed); 00712 } 00713 } 00714 } 00715 if (cis) { 00716 Tstate = newState; 00717 00718 } 00719 } 00720 00721 double Pinching4Material::posEnvlpStress(double u) 00722 { 00723 double k = 0.0; 00724 int i = 0; 00725 double f = 0.0; 00726 while (k==0.0 && i<=4){ 00727 00728 if (u<=envlpPosStrain(i+1)){ 00729 k = (envlpPosDamgdStress(i+1)-envlpPosDamgdStress(i))/(envlpPosStrain(i+1)-envlpPosStrain(i)); 00730 f = envlpPosDamgdStress(i) + (u-envlpPosStrain(i))*k; 00731 } 00732 i++; 00733 } 00734 00735 00736 if (k==0.0){ 00737 k = (envlpPosDamgdStress(5) - envlpPosDamgdStress(4))/(envlpPosStrain(5) - envlpPosStrain(4)); 00738 f = envlpPosDamgdStress(5) + k*(u-envlpPosStrain(5)); 00739 } 00740 00741 return f; 00742 00743 } 00744 00745 double Pinching4Material::posEnvlpTangent(double u) 00746 { 00747 double k = 0.0; 00748 int i = 0; 00749 while (k==0.0 && i<=4){ 00750 00751 if (u<=envlpPosStrain(i+1)){ 00752 k = (envlpPosDamgdStress(i+1)-envlpPosDamgdStress(i))/(envlpPosStrain(i+1)-envlpPosStrain(i)); 00753 } 00754 i++; 00755 } 00756 00757 if (k==0.0){ 00758 k = (envlpPosDamgdStress(5) - envlpPosDamgdStress(4))/(envlpPosStrain(5) - envlpPosStrain(4)); 00759 } 00760 00761 return k; 00762 00763 } 00764 00765 00766 double Pinching4Material::negEnvlpStress(double u) 00767 { 00768 double k = 0.0; 00769 int i = 0; 00770 double f = 0.0; 00771 while (k==0.0 && i<=4){ 00772 if (u>=envlpNegStrain(i+1)){ 00773 k = (envlpNegDamgdStress(i)-envlpNegDamgdStress(i+1))/(envlpNegStrain(i)-envlpNegStrain(i+1)); 00774 f = envlpNegDamgdStress(i+1)+(u-envlpNegStrain(i+1))*k; 00775 } 00776 i++; 00777 } 00778 00779 if (k==0.0){ 00780 k = (envlpNegDamgdStress(4) - envlpNegDamgdStress(5))/(envlpNegStrain(4)-envlpNegStrain(5)); 00781 f = envlpNegDamgdStress(5) + k*(u-envlpNegStrain(5)); 00782 } 00783 return f; 00784 00785 } 00786 00787 double Pinching4Material::negEnvlpTangent(double u) 00788 { 00789 double k = 0.0; 00790 int i = 0; 00791 while (k==0.0 && i<=4){ 00792 00793 if (u>=envlpNegStrain(i+1)){ 00794 k = (envlpNegDamgdStress(i)-envlpNegDamgdStress(i+1))/(envlpNegStrain(i)-envlpNegStrain(i+1)); 00795 } 00796 i++; 00797 } 00798 00799 if (k==0.0){ 00800 k = (envlpNegDamgdStress(4) - envlpNegDamgdStress(5))/(envlpNegStrain(4)-envlpNegStrain(5)); 00801 } 00802 return k; 00803 00804 } 00805 00806 00807 void Pinching4Material::getState3(Vector& state3Strain, Vector& state3Stress, double kunload) 00808 { 00809 00810 double kmax = (kunload>kElasticNegDamgd) ? kunload:kElasticNegDamgd; 00811 00812 if (state3Strain(0)*state3Strain(3) <0.0){ 00813 // trilinear unload reload path expected, first define point for reloading 00814 state3Strain(1) = lowTstateStrain*rDispN; 00815 if (rForceN-uForceN > 1e-8) { 00816 state3Stress(1) = lowTstateStress*rForceN; 00817 } 00818 else { 00819 if (TminStrainDmnd < envlpNegStrain(3)) { 00820 double st1 = lowTstateStress*uForceN*(1.0+1e-6); 00821 double st2 = envlpNegDamgdStress(4)*(1.0+1e-6); 00822 state3Stress(1) = (st1<st2) ? st1:st2; 00823 } 00824 else { 00825 double st1 = envlpNegDamgdStress(3)*uForceN*(1.0+1e-6); 00826 double st2 = envlpNegDamgdStress(4)*(1.0+1e-6); 00827 state3Stress(1) = (st1<st2) ? st1:st2; 00828 } 00829 } 00830 // if reload stiffness exceeds unload stiffness, reduce reload stiffness to make it equal to unload stiffness 00831 if ((state3Stress(1)-state3Stress(0))/(state3Strain(1)-state3Strain(0)) > kElasticNegDamgd) { 00832 state3Strain(1) = lowTstateStrain + (state3Stress(1)-state3Stress(0))/kElasticNegDamgd; 00833 } 00834 // check that reloading point is not behind point 4 00835 if (state3Strain(1)>state3Strain(3)) { 00836 // path taken to be a straight line between points 1 and 4 00837 double du = state3Strain(3) - state3Strain(0); 00838 double df = state3Stress(3) - state3Stress(0); 00839 state3Strain(1) = state3Strain(0) + 0.33*du; 00840 state3Strain(2) = state3Strain(0) + 0.67*du; 00841 state3Stress(1) = state3Stress(0) + 0.33*df; 00842 state3Stress(2) = state3Stress(0) + 0.67*df; 00843 } 00844 else { 00845 if (TminStrainDmnd < envlpNegStrain(3)) { 00846 state3Stress(2) = uForceN*envlpNegDamgdStress(4); 00847 } 00848 else { 00849 state3Stress(2) = uForceN*envlpNegDamgdStress(3); 00850 } 00851 state3Strain(2) = hghTstateStrain - (hghTstateStress-state3Stress(2))/kunload; 00852 00853 if (state3Strain(2) > state3Strain(3)) { 00854 // point 3 should be along a line between 2 and 4 00855 double du = state3Strain(3) - state3Strain(1); 00856 double df = state3Stress(3) - state3Stress(1); 00857 state3Strain(2) = state3Strain(1) + 0.5*du; 00858 state3Stress(2) = state3Stress(1) + 0.5*df; 00859 } 00860 else if ((state3Stress(2) - state3Stress(1))/(state3Strain(2) - state3Strain(1)) > kmax) { 00861 // linear unload-reload path expected 00862 double du = state3Strain(3) - state3Strain(0); 00863 double df = state3Stress(3) - state3Stress(0); 00864 state3Strain(1) = state3Strain(0) + 0.33*du; 00865 state3Strain(2) = state3Strain(0) + 0.67*du; 00866 state3Stress(1) = state3Stress(0) + 0.33*df; 00867 state3Stress(2) = state3Stress(0) + 0.67*df; 00868 } 00869 else if ((state3Strain(2) < state3Strain(1))||((state3Stress(2)-state3Stress(1))/(state3Strain(2)-state3Strain(1))<0)) { 00870 if (state3Strain(2)<0.0) { 00871 // pt 3 should be along a line between 2 and 4 00872 double du = state3Strain(3)-state3Strain(1); 00873 double df = state3Stress(3)-state3Stress(1); 00874 state3Strain(2) = state3Strain(1) + 0.5*du; 00875 state3Stress(2) = state3Stress(1) + 0.5*df; 00876 } 00877 else if (state3Strain(1) > 0.0) { 00878 // pt 2 should be along a line between 1 and 3 00879 double du = state3Strain(2)-state3Strain(0); 00880 double df = state3Stress(2)-state3Stress(0); 00881 state3Strain(1) = state3Strain(0) + 0.5*du; 00882 state3Stress(1) = state3Stress(0) + 0.5*df; 00883 } 00884 else { 00885 double avgforce = 0.5*(state3Stress(2) + state3Stress(1)); 00886 double dfr = 0.0; 00887 if (avgforce < 0.0){ 00888 dfr = -avgforce/100; 00889 } 00890 else { 00891 dfr = avgforce/100; 00892 } 00893 double slope12 = (state3Stress(1) - state3Stress(0))/(state3Strain(1) - state3Strain(0)); 00894 double slope34 = (state3Stress(3) - state3Stress(2))/(state3Strain(3) - state3Strain(2)); 00895 state3Stress(1) = avgforce - dfr; 00896 state3Stress(2) = avgforce + dfr; 00897 state3Strain(1) = state3Strain(0) + (state3Stress(1) - state3Stress(0))/slope12; 00898 state3Strain(2) = state3Strain(3) - (state3Stress(3) - state3Stress(2))/slope34; 00899 } 00900 } 00901 } 00902 } 00903 else { 00904 // linear unload reload path is expected 00905 double du = state3Strain(3)-state3Strain(0); 00906 double df = state3Stress(3)-state3Stress(0); 00907 state3Strain(1) = state3Strain(0) + 0.33*du; 00908 state3Strain(2) = state3Strain(0) + 0.67*du; 00909 state3Stress(1) = state3Stress(0) + 0.33*df; 00910 state3Stress(2) = state3Stress(0) + 0.67*df; 00911 } 00912 00913 00914 double checkSlope = state3Stress(0)/state3Strain(0); 00915 double slope = 0.0; 00916 00917 // final check 00918 int i = 0; 00919 while (i<3) { 00920 double du = state3Strain(i+1)-state3Strain(i); 00921 double df = state3Stress(i+1)-state3Stress(i); 00922 if (du<0.0 || df<0.0) { 00923 double du = state3Strain(3)-state3Strain(0); 00924 double df = state3Stress(3)-state3Stress(0); 00925 state3Strain(1) = state3Strain(0) + 0.33*du; 00926 state3Strain(2) = state3Strain(0) + 0.67*du; 00927 state3Stress(1) = state3Stress(0) + 0.33*df; 00928 state3Stress(2) = state3Stress(0) + 0.67*df; 00929 slope = df/du; 00930 i = 3; 00931 } 00932 if (slope > 1e-8 && slope < checkSlope) { 00933 state3Strain(1) = 0.0; state3Stress(1) = 0.0; 00934 state3Strain(2) = state3Strain(3)/2; state3Stress(2) = state3Stress(3)/2; 00935 } 00936 i++; 00937 } 00938 00939 00940 } 00941 00942 void Pinching4Material::getState4(Vector& state4Strain,Vector& state4Stress, double kunload) 00943 { 00944 00945 double kmax = (kunload>kElasticPosDamgd) ? kunload:kElasticPosDamgd; 00946 00947 if (state4Strain(0)*state4Strain(3) <0.0){ 00948 // trilinear unload reload path expected 00949 state4Strain(2) = hghTstateStrain*rDispP; 00950 if (uForceP==0.0){ 00951 state4Stress(2) = hghTstateStress*rForceP; 00952 } 00953 else if (rForceP-uForceP > 1e-8) { 00954 state4Stress(2) = hghTstateStress*rForceP; 00955 } 00956 else { 00957 if (TmaxStrainDmnd > envlpPosStrain(3)) { 00958 double st1 = hghTstateStress*uForceP*(1.0+1e-6); 00959 double st2 = envlpPosDamgdStress(4)*(1.0+1e-6); 00960 state4Stress(2) = (st1>st2) ? st1:st2; 00961 } 00962 else { 00963 double st1 = envlpPosDamgdStress(3)*uForceP*(1.0+1e-6); 00964 double st2 = envlpPosDamgdStress(4)*(1.0+1e-6); 00965 state4Stress(2) = (st1>st2) ? st1:st2; 00966 } 00967 } 00968 // if reload stiffness exceeds unload stiffness, reduce reload stiffness to make it equal to unload stiffness 00969 if ((state4Stress(3)-state4Stress(2))/(state4Strain(3)-state4Strain(2)) > kElasticPosDamgd) { 00970 state4Strain(2) = hghTstateStrain - (state4Stress(3)-state4Stress(2))/kElasticPosDamgd; 00971 } 00972 // check that reloading point is not behind point 1 00973 if (state4Strain(2)<state4Strain(0)) { 00974 // path taken to be a straight line between points 1 and 4 00975 double du = state4Strain(3) - state4Strain(0); 00976 double df = state4Stress(3) - state4Stress(0); 00977 state4Strain(1) = state4Strain(0) + 0.33*du; 00978 state4Strain(2) = state4Strain(0) + 0.67*du; 00979 state4Stress(1) = state4Stress(0) + 0.33*df; 00980 state4Stress(2) = state4Stress(0) + 0.67*df; 00981 } 00982 else { 00983 if (TmaxStrainDmnd > envlpPosStrain(3)) { 00984 state4Stress(1) = uForceP*envlpPosDamgdStress(4); 00985 } 00986 else { 00987 state4Stress(1) = uForceP*envlpPosDamgdStress(3); 00988 } 00989 state4Strain(1) = lowTstateStrain + (-lowTstateStress+state4Stress(1))/kunload; 00990 00991 if (state4Strain(1) < state4Strain(0)) { 00992 // point 2 should be along a line between 1 and 3 00993 double du = state4Strain(2) - state4Strain(0); 00994 double df = state4Stress(2) - state4Stress(0); 00995 state4Strain(1) = state4Strain(0) + 0.5*du; 00996 state4Stress(1) = state4Stress(0) + 0.5*df; 00997 } 00998 else if ((state4Stress(2) - state4Stress(1))/(state4Strain(2) - state4Strain(1)) > kmax) { 00999 // linear unload-reload path expected 01000 double du = state4Strain(3) - state4Strain(0); 01001 double df = state4Stress(3) - state4Stress(0); 01002 state4Strain(1) = state4Strain(0) + 0.33*du; 01003 state4Strain(2) = state4Strain(0) + 0.67*du; 01004 state4Stress(1) = state4Stress(0) + 0.33*df; 01005 state4Stress(2) = state4Stress(0) + 0.67*df; 01006 } 01007 else if ((state4Strain(2) < state4Strain(1))||((state4Stress(2)-state4Stress(1))/(state4Strain(2)-state4Strain(1))<0)) { 01008 if (state4Strain(1)>0.0) { 01009 // pt 2 should be along a line between 1 and 3 01010 double du = state4Strain(2)-state4Strain(0); 01011 double df = state4Stress(2)-state4Stress(0); 01012 state4Strain(1) = state4Strain(0) + 0.5*du; 01013 state4Stress(1) = state4Stress(0) + 0.5*df; 01014 } 01015 else if (state4Strain(2) < 0.0) { 01016 // pt 2 should be along a line between 2 and 4 01017 double du = state4Strain(3)-state4Strain(1); 01018 double df = state4Stress(3)-state4Stress(1); 01019 state4Strain(2) = state4Strain(1) + 0.5*du; 01020 state4Stress(2) = state4Stress(1) + 0.5*df; 01021 } 01022 else { 01023 double avgforce = 0.5*(state4Stress(2) + state4Stress(1)); 01024 double dfr = 0.0; 01025 if (avgforce < 0.0){ 01026 dfr = -avgforce/100; 01027 } 01028 else { 01029 dfr = avgforce/100; 01030 } 01031 double slope12 = (state4Stress(1) - state4Stress(0))/(state4Strain(1) - state4Strain(0)); 01032 double slope34 = (state4Stress(3) - state4Stress(2))/(state4Strain(3) - state4Strain(2)); 01033 state4Stress(1) = avgforce - dfr; 01034 state4Stress(2) = avgforce + dfr; 01035 state4Strain(1) = state4Strain(0) + (state4Stress(1) - state4Stress(0))/slope12; 01036 state4Strain(2) = state4Strain(3) - (state4Stress(3) - state4Stress(2))/slope34; 01037 } 01038 } 01039 } 01040 } 01041 else { 01042 // linear unload reload path is expected 01043 double du = state4Strain(3)-state4Strain(0); 01044 double df = state4Stress(3)-state4Stress(0); 01045 state4Strain(1) = state4Strain(0) + 0.33*du; 01046 state4Strain(2) = state4Strain(0) + 0.67*du; 01047 state4Stress(1) = state4Stress(0) + 0.33*df; 01048 state4Stress(2) = state4Stress(0) + 0.67*df; 01049 } 01050 01051 01052 01053 double checkSlope = state4Stress(0)/state4Strain(0); 01054 double slope = 0.0; 01055 01056 // final check 01057 int i = 0; 01058 while (i<3) { 01059 double du = state4Strain(i+1)-state4Strain(i); 01060 double df = state4Stress(i+1)-state4Stress(i); 01061 if (du<0.0 || df<0.0) { 01062 double du = state4Strain(3)-state4Strain(0); 01063 double df = state4Stress(3)-state4Stress(0); 01064 state4Strain(1) = state4Strain(0) + 0.33*du; 01065 state4Strain(2) = state4Strain(0) + 0.67*du; 01066 state4Stress(1) = state4Stress(0) + 0.33*df; 01067 state4Stress(2) = state4Stress(0) + 0.67*df; 01068 slope = df/du; 01069 i = 3; 01070 } 01071 if (slope > 1e-8 && slope < checkSlope) { 01072 state4Strain(1) = 0.0; state4Stress(1) = 0.0; 01073 state4Strain(2) = state4Strain(3)/2; state4Stress(2) = state4Stress(3)/2; 01074 } 01075 01076 i++; 01077 } 01078 } 01079 01080 double Pinching4Material::Envlp3Tangent(Vector s3Strain, Vector s3Stress, double u) 01081 { 01082 double k = 0.0; 01083 int i = 0; 01084 while ((k==0.0||i<=2) && (i<=2)) 01085 { 01086 if (u>= s3Strain(i)) { 01087 k = (s3Stress(i+1)-s3Stress(i))/(s3Strain(i+1)-s3Strain(i)); 01088 } 01089 i++; 01090 } 01091 if (k==0.0) { 01092 if (u<s3Strain(0)) { 01093 i = 0; 01094 } 01095 else { 01096 i = 2; 01097 } 01098 k = (s3Stress(i+1)-s3Stress(i))/(s3Strain(i+1)-s3Strain(i)); 01099 01100 } 01101 return k; 01102 } 01103 01104 double Pinching4Material::Envlp4Tangent(Vector s4Strain, Vector s4Stress, double u) 01105 { 01106 double k = 0.0; 01107 int i = 0; 01108 while ((k==0.0||i<=2) && (i<=2)) 01109 { 01110 if (u>= s4Strain(i)) { 01111 k = (s4Stress(i+1)-s4Stress(i))/(s4Strain(i+1)-s4Strain(i)); 01112 } 01113 i++; 01114 } 01115 if (k==0.0) { 01116 if (u<s4Strain(0)) { 01117 i = 0; 01118 } 01119 else { 01120 i = 2; 01121 } 01122 k = (s4Stress(i+1)-s4Stress(i))/(s4Strain(i+1)-s4Strain(i)); 01123 01124 } 01125 return k; 01126 } 01127 01128 01129 double Pinching4Material::Envlp3Stress(Vector s3Strain, Vector s3Stress, double u) 01130 { 01131 double k = 0.0; 01132 int i = 0; 01133 double f = 0.0; 01134 while ((k==0.0||i<=2) && (i<=2)) 01135 { 01136 if (u>= s3Strain(i)) { 01137 k = (s3Stress(i+1)-s3Stress(i))/(s3Strain(i+1)-s3Strain(i)); 01138 f = s3Stress(i)+(u-s3Strain(i))*k; 01139 } 01140 i++; 01141 } 01142 if (k==0.0) { 01143 if (u<s3Strain(0)) { 01144 i = 0; 01145 } 01146 else { 01147 i = 2; 01148 } 01149 k = (s3Stress(i+1)-s3Stress(i))/(s3Strain(i+1)-s3Strain(i)); 01150 f = s3Stress(i)+(u-s3Strain(i))*k; 01151 } 01152 return f; 01153 } 01154 01155 double Pinching4Material::Envlp4Stress(Vector s4Strain, Vector s4Stress, double u) 01156 { 01157 double k = 0.0; 01158 int i = 0; 01159 double f = 0.0; 01160 while ((k==0.0||i<=2) && (i<=2)) 01161 { 01162 if (u>= s4Strain(i)) { 01163 k = (s4Stress(i+1)-s4Stress(i))/(s4Strain(i+1)-s4Strain(i)); 01164 f = s4Stress(i)+(u-s4Strain(i))*k; 01165 } 01166 i++; 01167 } 01168 if (k==0.0) { 01169 if (u<s4Strain(0)) { 01170 i = 0; 01171 } 01172 else { 01173 i = 2; 01174 } 01175 k = (s4Stress(i+1)-s4Stress(i))/(s4Strain(i+1)-s4Strain(i)); 01176 f = s4Stress(i)+(u-s4Strain(i))*k; 01177 } 01178 return f; 01179 } 01180 01181 void Pinching4Material::updateDmg(double strain, double dstrain) 01182 { 01183 double tes = 0.0; 01184 double umaxAbs = (TmaxStrainDmnd>-TminStrainDmnd) ? TmaxStrainDmnd:-TminStrainDmnd; 01185 double uultAbs = (envlpPosStrain(4)>-envlpNegStrain(4)) ? envlpPosStrain(4):-envlpNegStrain(4); 01186 TnCycle = CnCycle + fabs(dstrain)/(4*umaxAbs); 01187 if ((strain<uultAbs && strain>-uultAbs)&& Tenergy< energyCapacity) 01188 { 01189 TgammaK = gammaK1*pow((umaxAbs/uultAbs),gammaK3); 01190 TgammaD = gammaD1*pow((umaxAbs/uultAbs),gammaD3); 01191 TgammaF = gammaF1*pow((umaxAbs/uultAbs),gammaF3); 01192 01193 if (Tenergy>elasticStrainEnergy && DmgCyc == 0) { 01194 tes = ((Tenergy-elasticStrainEnergy)/energyCapacity); 01195 TgammaK = TgammaK + gammaK2*pow(tes,gammaK4); 01196 TgammaD = TgammaD + gammaD2*pow(tes,gammaD4); 01197 TgammaF = TgammaF + gammaF2*pow(tes,gammaF4); 01198 } else if (DmgCyc == 1) { 01199 TgammaK = TgammaK + gammaK2*pow(TnCycle,gammaK4); 01200 TgammaD = TgammaD + gammaD2*pow(TnCycle,gammaD4); 01201 TgammaF = TgammaF + gammaF2*pow(TnCycle,gammaF4); 01202 } 01203 double kminP = (posEnvlpStress(TmaxStrainDmnd)/TmaxStrainDmnd); 01204 double kminN = (negEnvlpStress(TminStrainDmnd)/TminStrainDmnd); 01205 double kmin = ((kminP/kElasticPos)>(kminN/kElasticNeg)) ? (kminP/kElasticPos):(kminN/kElasticNeg); 01206 double gammaKLimEnv = (0.0>(1.0-kmin)) ? 0.0:(1.0-kmin); 01207 01208 double k1 = (TgammaK<gammaKLimit) ? TgammaK:gammaKLimit; 01209 TgammaK = (k1<gammaKLimEnv) ? k1:gammaKLimEnv; 01210 TgammaD = (TgammaD<gammaDLimit) ? TgammaD:gammaDLimit; 01211 TgammaF = (TgammaF<gammaFLimit) ? TgammaF:gammaFLimit; 01212 } 01213 else if (strain<uultAbs && strain>-uultAbs) { 01214 double kminP = (posEnvlpStress(TmaxStrainDmnd)/TmaxStrainDmnd); 01215 double kminN = (negEnvlpStress(TminStrainDmnd)/TminStrainDmnd); 01216 double kmin = ((kminP/kElasticPos)>(kminN/kElasticNeg)) ? (kminP/kElasticPos):(kminN/kElasticNeg); 01217 double gammaKLimEnv = (0.0>(1.0-kmin)) ? 0.0:(1.0-kmin); 01218 01219 TgammaK = (gammaKLimit<gammaKLimEnv) ? gammaKLimit:gammaKLimEnv; 01220 TgammaD = gammaDLimit; 01221 TgammaF = gammaFLimit; 01222 } 01223 01224 }
Generated on Mon Oct 23 15:05:21 2006 for OpenSees by |
opensees-support @ berkeley.edu
Supported by the National Science Foundation
©2006, UC Regents