EL_NLEeq.cppGo to the documentation of this file.00001 /* 00002 //================================================================================ 00003 # COPYRIGHT (C): :-)) # 00004 # PROJECT: Object Oriented Finite Element Program # 00005 # PURPOSE: General platform for elaso-plastic constitutive model # 00006 # implementation # 00007 # # 00008 # CLASS: EvolutionLaw_NL_Eeq (on plastic equivalent strain) # 00009 # # 00010 # VERSION: # 00011 # LANGUAGE: C++.ver >= 2.0 ( Borland C++ ver=3.00, SUN C++ ver=2.1 ) # 00012 # TARGET OS: DOS || UNIX || . . . # 00013 # DESIGNER(S): Boris Jeremic, Zhaohui Yang # 00014 # PROGRAMMER(S): Boris Jeremic, Zhaohui Yang # 00015 # # 00016 # # 00017 # DATE: 09-02-2000 # 00018 # UPDATE HISTORY: # 00019 # # 00020 # # 00021 # # 00022 # SHORT EXPLANATION: This is a nonlinear evolution law for the evoltion of a # 00023 # scalar variable k which depends on plastic equi. strain # 00024 # i.e. dk = f(e_eq)*de_eq # 00025 //================================================================================ 00026 */ 00027 00028 #ifndef EL_NLEeq_CPP 00029 #define EL_NLEeq_CPP 00030 00031 #include "EL_NLEeq.h" 00032 #include <basics.h> 00033 00034 //================================================================================ 00035 // Copy constructor 00036 //================================================================================ 00037 00038 EvolutionLaw_NL_Eeq::EvolutionLaw_NL_Eeq(const EvolutionLaw_NL_Eeq &NLE ) { 00039 00040 this->eeqEtaPeak = NLE.geteeqEtaPeak(); 00041 this->etaResidual = NLE.getetaResidual(); 00042 this->etaStart = NLE.getetaStart(); 00043 this->etaPeak = NLE.getetaPeak(); 00044 this->e = NLE.gete(); 00045 this->d = NLE.getd(); 00046 } 00047 00048 00049 //================================================================================ 00050 // Create a clone of itself 00051 //================================================================================ 00052 EvolutionLaw_S* EvolutionLaw_NL_Eeq::newObj() { 00053 00054 EvolutionLaw_S *newEL = new EvolutionLaw_NL_Eeq( *this ); 00055 00056 return newEL; 00057 00058 } 00059 00064 // 00065 //void EvolutionLaw_NL_Eeq::InitVars(EPState *EPS) { 00066 // 00067 // // set initial E_Young corresponding to current stress state 00068 // //double p_atm = 100.0; //Kpa atmospheric pressure 00069 // //double p = EPS->getStress().p_hydrostatic(); 00070 // //double E = EPS->getEo() * pow( (p/p_atm), geta()); 00071 // EPS->setE( EPS->getEo() ); 00072 // 00073 //} 00074 00075 00076 //================================================================================ 00077 // Set initial value of D once the current stress hit the yield surface 00078 // for NL model only 00079 // 00080 // 00081 //================================================================================ 00082 // 00083 //void EvolutionLaw_NL_Eeq::setInitD(EPState *EPS) { 00084 // 00085 //} 00086 00090 // 00091 //void EvolutionLaw_NL_Eeq::UpdateVar( EPState *EPS, int WhichOne) { 00092 // 00093 // //========================================================================= 00094 // // Updating corresponding internal var by a nonliear function f 00095 // 00096 // // Getting e_eq 00097 // straintensor pstrain = EPS->getPlasticStrain(); 00098 // double e_eq = pstrain.equivalent(); 00099 // double e_eq2 = e_eq * e_eq; 00100 // //cout << "e_eq = " << e_eq << endlnn; 00101 // 00102 // double upper1 = getd()*getetaResidual()*e_eq2; 00103 // double upper2 = (gete()*getetaPeak()+2.0*getd()*geteeqEtaPeak()*(-getetaResidual()+getetaPeak()))*e_eq; 00104 // double upper3 = getetaStart()*getd()*pow(geteeqEtaPeak(),2.0)*(-getetaResidual()+getetaPeak())/(getetaPeak()-getetaStart()); 00105 // double upper = upper1 + upper2 + upper3; 00106 // 00107 // double lower = getd()*e_eq2 + gete()*e_eq +getd()*pow(geteeqEtaPeak(), 2.0)*(-getetaResidual()+getetaPeak())/(getetaPeak()-getetaStart()); 00108 // 00109 // double new_S = upper/lower; 00110 // 00111 // EPS->setScalarVar(WhichOne, new_S); 00112 // 00113 //} 00114 00115 00116 00117 //================================================================================ 00118 // Evaluating h_s ( for the evaluation of Kp ) 00119 //================================================================================ 00120 //double EvolutionLaw_NL_Eeq::h( EPState *EPS, double norm_dQods_dev ) { 00121 00122 double EvolutionLaw_NL_Eeq::h_s( EPState *EPS, PotentialSurface *PS) 00123 { 00124 double h; 00125 00126 //========================================================================= 00127 // Getting de_eq/dLambda 00128 //double de_eq = EPS->getdPlasticStrain().equivalent(); 00129 stresstensor dQods = PS->dQods( EPS ); 00130 tensor dQods_dev = dQods.deviator(); 00131 00132 //Evaluate the norm of the deviator of dQods 00133 //temp1 = dQods("ij")*dQods("ij"); 00134 tensor temp1 = dQods_dev("ij")*dQods_dev("ij"); 00135 double norm_dQods_dev = pow( temp1.trace(), 0.5 ); 00136 double de_eqodL = pow( 2.0 / 3.0, 0.5 ) * norm_dQods_dev; 00137 00138 //Evaluating dSodeeq 00139 double dSodeeq; 00140 straintensor pstrain = EPS->getPlasticStrain(); //bug! should be total plastic strain's deviatoric part???? 00141 straintensor pstrain_dev = pstrain.deviator(); 00142 00143 //Checking equivalent() 00144 //Using deviatoric plastic strain to evolve the scalar internal var! 00145 //tensor temptx = pstrain("ij") * pstrain("ij"); 00146 //double tempdx = temptx.trace(); 00147 //double e_eq = pow( 2.0 * tempdx / 3.0, 0.5 ); 00149 00150 tensor temptx = pstrain_dev("ij") * pstrain_dev("ij"); 00151 double tempdx = temptx.trace(); 00152 double e_eq = pow( 2.0 * tempdx / 3.0, 0.5 ); 00153 00154 00155 //double e_eq = pstrain.equivalent(); ????????????why this wouldn't work? 00156 double e_eq2 = e_eq * e_eq; 00157 00158 double upper1 = -getd()*(getetaPeak() - getetaResidual())*(getetaPeak() - getetaStart())*( e_eq- geteeqEtaPeak()); 00159 double upper21 = gete()*(getetaPeak() - getetaStart())* (e_eq + geteeqEtaPeak()); 00160 double upper22 = 2.0*getd()* geteeqEtaPeak()* (-(getetaStart()*e_eq)-getetaResidual()* geteeqEtaPeak() + getetaPeak()* (e_eq+geteeqEtaPeak())); 00161 double upper = upper1 * (upper21 + upper22); 00162 00163 double lower1 = gete()*(getetaPeak()-getetaStart())*e_eq; 00164 double lower2 = getd()*(getetaPeak()*(e_eq2+pow(geteeqEtaPeak(), 2.0) )-getetaStart()*e_eq2-getetaResidual()*pow(geteeqEtaPeak(), 2.0)); 00165 double lower = pow(lower1 + lower2, 2.0); 00166 00167 dSodeeq = upper / lower; 00168 00169 h = dSodeeq * de_eqodL; 00170 00171 // Drucker-Prager's evolution law 00172 // Get the current stress's I1 00173 //double I1 = EPS->getStress().Iinvariant1(); 00174 00175 //Von Mises 00176 //double k = EPS->getScalarVar(1); 00177 //double Kp = -2.0 * k * geta() * temp; 00178 00179 return h; 00180 00181 } 00182 00183 00184 //================================================================================ 00185 // Print vars defined in NLinear Evolution Law 00186 //================================================================================ 00187 void EvolutionLaw_NL_Eeq::print() 00188 { 00189 opserr << (*this); 00190 } 00191 00192 00193 //================================================================================ 00194 double EvolutionLaw_NL_Eeq::geteeqEtaPeak() const 00195 { 00196 return eeqEtaPeak; 00197 } 00198 00199 //================================================================================ 00200 double EvolutionLaw_NL_Eeq::getetaResidual() const 00201 { 00202 return etaResidual; 00203 } 00204 00205 //================================================================================ 00206 double EvolutionLaw_NL_Eeq::getetaStart() const 00207 { 00208 return etaStart; 00209 } 00210 00211 //================================================================================ 00212 double EvolutionLaw_NL_Eeq::getetaPeak() const 00213 { 00214 return etaPeak; 00215 } 00216 00217 //================================================================================ 00218 double EvolutionLaw_NL_Eeq::gete() const 00219 { 00220 return e; 00221 } 00222 00223 //================================================================================ 00224 double EvolutionLaw_NL_Eeq::getd() const 00225 { 00226 return d; 00227 } 00228 00229 //================================================================================ 00230 OPS_Stream& operator<< (OPS_Stream& os, const EvolutionLaw_NL_Eeq & NLEL) 00231 { 00232 // os.unsetf( ios::scientific ); 00233 os.precision(5); 00234 00235 //os.width(10); 00236 os << endln << "Nonlinear Evolution(plastic equivalent strain ) Law's parameters:" << endln; 00237 os << "eeqEtaPeak = " << NLEL.geteeqEtaPeak() << "; "; 00238 os << "etaResidual = " << NLEL.getetaResidual() << "; " << endln; 00239 os << "etaStart = " << NLEL.getetaStart() << "; "; 00240 os << "etaPeak = " << NLEL.getetaPeak() << "; "; 00241 os << "e = " << NLEL.gete() << "; "; 00242 os << "d = " << NLEL.getd() << "; " << endln; 00243 //os.width(10); 00244 00245 return os; 00246 } 00247 00248 #endif 00249 |