diff --git a/SRC/Makefile b/SRC/Makefile index 7567a9c9a..eea6c00be 100644 --- a/SRC/Makefile +++ b/SRC/Makefile @@ -1065,6 +1065,7 @@ SECTION_LIBS = $(FE)/material/section/SectionForceDeformation.o \ $(FE)/material/section/integration/WideFlangeSectionIntegration.o \ $(FE)/material/section/integration/RCTBeamSectionIntegration.o \ $(FE)/material/section/integration/TubeSectionIntegration.o \ + $(FE)/material/section/integration/HSSSectionIntegration.o \ $(FE)/material/section/repres/patch/Patch.o \ $(FE)/material/section/repres/patch/QuadPatch.o \ $(FE)/material/section/repres/patch/CircPatch.o \ diff --git a/SRC/classTags.h b/SRC/classTags.h index 3021d506b..97650667c 100644 --- a/SRC/classTags.h +++ b/SRC/classTags.h @@ -340,6 +340,7 @@ #define SECTION_INTEGRATION_TAG_RCCIRCULAR 5 #define SECTION_INTEGRATION_TAG_RCTUNNEL 6 #define SECTION_INTEGRATION_TAG_Tube 7 +#define SECTION_INTEGRATION_TAG_HSS 8 #define ND_TAG_WrapperNDMaterial 9 #define ND_TAG_ElasticIsotropic 10 diff --git a/SRC/interpreter/OpenSeesSectionCommands.cpp b/SRC/interpreter/OpenSeesSectionCommands.cpp index 94f3e45d6..7a40909f7 100644 --- a/SRC/interpreter/OpenSeesSectionCommands.cpp +++ b/SRC/interpreter/OpenSeesSectionCommands.cpp @@ -47,6 +47,7 @@ #include #include #include +#include #include #include #include @@ -91,6 +92,7 @@ void* OPS_Bidirectional(); void* OPS_Elliptical2(); void* OPS_Isolator2spring(); void* OPS_FiberSection2dThermal(); +void* OPS_HSSSection(); namespace { static FiberSection2d* theActiveFiberSection2d = 0; @@ -1086,6 +1088,7 @@ namespace { functionMap.insert(std::make_pair("ElasticWarpingShear", &OPS_ElasticWarpingShearSection2d)); functionMap.insert(std::make_pair("ElasticTube", &OPS_ElasticTubeSection3d)); functionMap.insert(std::make_pair("Tube", &OPS_TubeSection)); + functionMap.insert(std::make_pair("HSS", &OPS_HSSSection)); functionMap.insert(std::make_pair("WFSection2d", &OPS_WFSection2d)); functionMap.insert(std::make_pair("WSection2d", &OPS_WFSection2d)); functionMap.insert(std::make_pair("RCSection2d", &OPS_RCSection2d)); diff --git a/SRC/material/nD/J2BeamFiber2d.cpp b/SRC/material/nD/J2BeamFiber2d.cpp index 1ff76a80e..94cfe2056 100644 --- a/SRC/material/nD/J2BeamFiber2d.cpp +++ b/SRC/material/nD/J2BeamFiber2d.cpp @@ -700,7 +700,7 @@ J2BeamFiber2d::setParameter(const char **argv, int argc, param.setValue(nu); return param.addObject(2, this); } - else if (strcmp(argv[0],"sigmaY") == 0 || strcmp(argv[0],"fy") == 0) { + else if (strcmp(argv[0],"sigmaY") == 0 || strcmp(argv[0],"fy") == 0 || strcmp(argv[0],"Fy") == 0) { param.setValue(sigmaY); return param.addObject(5, this); } diff --git a/SRC/material/nD/J2BeamFiber3d.cpp b/SRC/material/nD/J2BeamFiber3d.cpp index e555857ee..3c32ce39a 100644 --- a/SRC/material/nD/J2BeamFiber3d.cpp +++ b/SRC/material/nD/J2BeamFiber3d.cpp @@ -759,7 +759,7 @@ J2BeamFiber3d::setParameter(const char **argv, int argc, param.setValue(nu); return param.addObject(2, this); } - else if (strcmp(argv[0],"sigmaY") == 0 || strcmp(argv[0],"fy") == 0) { + else if (strcmp(argv[0],"sigmaY") == 0 || strcmp(argv[0],"fy") == 0 || strcmp(argv[0],"Fy") == 0) { param.setValue(sigmaY); return param.addObject(5, this); } diff --git a/SRC/material/section/integration/CMakeLists.txt b/SRC/material/section/integration/CMakeLists.txt index d4d786865..545bf2613 100644 --- a/SRC/material/section/integration/CMakeLists.txt +++ b/SRC/material/section/integration/CMakeLists.txt @@ -13,6 +13,7 @@ target_sources(OPS_Material RCCircularSectionIntegration.cpp TubeSectionIntegration.cpp RCTunnelSectionIntegration.cpp + HSSSectionIntegration.cpp PUBLIC SectionIntegration.h WideFlangeSectionIntegration.h @@ -21,6 +22,7 @@ target_sources(OPS_Material RCCircularSectionIntegration.h TubeSectionIntegration.h RCTunnelSectionIntegration.h + HSSSectionIntegration.h ) target_include_directories(OPS_Material PUBLIC ${CMAKE_CURRENT_LIST_DIR}) diff --git a/SRC/material/section/integration/HSSSectionIntegration.cpp b/SRC/material/section/integration/HSSSectionIntegration.cpp new file mode 100644 index 000000000..56f0aebb4 --- /dev/null +++ b/SRC/material/section/integration/HSSSectionIntegration.cpp @@ -0,0 +1,544 @@ +/* ****************************************************************** ** +** OpenSees - Open System for Earthquake Engineering Simulation ** +** Pacific Earthquake Engineering Research Center ** +** ** +** ** +** (C) Copyright 1999, The Regents of the University of California ** +** All Rights Reserved. ** +** ** +** Commercial use of this program without express permission of the ** +** University of California, Berkeley, is strictly prohibited. See ** +** file 'COPYRIGHT' in main directory for information on usage and ** +** redistribution, and for a DISCLAIMER OF ALL WARRANTIES. ** +** ** +** Developed by: ** +** Frank McKenna (fmckenna@ce.berkeley.edu) ** +** Gregory L. Fenves (fenves@ce.berkeley.edu) ** +** Filip C. Filippou (filippou@ce.berkeley.edu) ** +** ** +** ****************************************************************** */ + +// $Revision: 1.4 $ +// $Date: 2010-09-13 21:31:26 $ +// $Source: /usr/local/cvs/OpenSees/SRC/material/section/integration/HSSSectionIntegration.cpp,v $ + +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include + +void* OPS_HSSSection() +{ + if (OPS_GetNumRemainingInputArgs() < 7) { + opserr << "WARNING insufficient arguments\n"; + opserr << "Want: section HSS tag? matTag? h? b? t? nfh? nfb? <-nd> <-shape shape?> <-GJ GJ?> <-torsion tag?>" << endln; + return 0; + } + + int ndm = OPS_GetNDM(); + + int tag, matTag; + double h, b, t; + int nfh, nfb; + int nft = 1; + + SectionForceDeformation* theSection = 0; + + int numdata = 1; + if (OPS_GetIntInput(&numdata, &tag) < 0) { + opserr << "WARNING invalid section HSS tag" << endln; + return 0; + } + + if (OPS_GetIntInput(&numdata, &matTag) < 0) { + opserr << "WARNING invalid section HSS matTag" << endln; + return 0; + } + + if (OPS_GetDoubleInput(&numdata, &h) < 0) { + opserr << "WARNING invalid h" << endln; + opserr << "HSS section: " << tag << endln; + return 0; + } + + if (OPS_GetDoubleInput(&numdata, &b) < 0) { + opserr << "WARNING invalid b" << endln; + opserr << "HSS section: " << tag << endln; + return 0; + } + + if (OPS_GetDoubleInput(&numdata, &t) < 0) { + opserr << "WARNING invalid t" << endln; + opserr << "HSS section: " << tag << endln; + return 0; + } + + if (OPS_GetIntInput(&numdata, &nfh) < 0) { + opserr << "WARNING invalid nfh" << endln; + opserr << "HSS section: " << tag << endln; + return 0; + } + + if (OPS_GetIntInput(&numdata, &nfb) < 0) { + opserr << "WARNING invalid nfb" << endln; + opserr << "HSS section: " << tag << endln; + return 0; + } + + HSSSectionIntegration hsssect(h, b, t, nfh, nfb, nft); + + int numFibers = hsssect.getNumFibers(); + bool isND = false; double shape = 1.0; + UniaxialMaterial *torsion = 0; + bool deleteTorsion = false; + while (OPS_GetNumRemainingInputArgs() > 0) { + const char* flag = OPS_GetString(); + // read <-nd> + if (strcmp(flag,"-nd") == 0) + isND = true; + // read <-shape shape> + if (strcmp(flag,"-shape") == 0 && OPS_GetNumRemainingInputArgs() > 0) { + if (OPS_GetDoubleInput(&numdata, &shape) < 0) { + opserr << "WARNING invalid shape" << endln; + opserr << "HSS section: " << tag << endln; + return 0; + } + isND = true; + } + // read <-GJ GJ> + if (strcmp(flag,"-GJ") == 0 && OPS_GetNumRemainingInputArgs() > 0) { + double GJ; + if (OPS_GetDoubleInput(&numdata, &GJ) < 0) { + opserr << "WARNING: failed to read GJ\n"; + return 0; + } + torsion = new ElasticMaterial(0,GJ); + deleteTorsion = true; + } + // read <-torsion tag> + if (strcmp(flag,"-torsion") == 0 && OPS_GetNumRemainingInputArgs() > 0) { + int torsionTag; + if (OPS_GetIntInput(&numdata, &torsionTag) < 0) { + opserr << "WARNING: failed to read torsion\n"; + return 0; + } + torsion = OPS_getUniaxialMaterial(torsionTag); + } + } + + if (isND) { + NDMaterial *theSteel = OPS_getNDMaterial(matTag); + if (theSteel == 0) { + opserr << "WARNING ND material does not exist\n"; + opserr << "material: " << matTag; + opserr << "\nHSS section: " << tag << endln; + return 0; + } + + NDMaterial **theMats = new NDMaterial *[numFibers]; + hsssect.arrangeFibers(theMats, theSteel); + + if (ndm == 2) + theSection = new NDFiberSection2d(tag, numFibers, theMats, hsssect); + if (ndm == 3) + theSection = new NDFiberSection3d(tag, numFibers, theMats, hsssect, shape); + + delete [] theMats; + } + + else { + UniaxialMaterial *theSteel = OPS_getUniaxialMaterial(matTag); + if (theSteel == 0) { + opserr << "WARNING uniaxial material does not exist\n"; + opserr << "material: " << matTag; + opserr << "\nHSS section: " << tag << endln; + return 0; + } + + if (torsion == 0) { + opserr << "WARNING torsion not speified for FiberSection\n"; + opserr << "\nHSS section: " << tag << endln; + return 0; + } + + UniaxialMaterial **theMats = new UniaxialMaterial *[numFibers]; + hsssect.arrangeFibers(theMats, theSteel); + + if (ndm == 2) + theSection = new FiberSection2d(tag, numFibers, theMats, hsssect); + if (ndm == 3) + theSection = new FiberSection3d(tag, numFibers, theMats, hsssect, *torsion); + + if (deleteTorsion) + delete torsion; + + delete [] theMats; + } + + return theSection; +} + +HSSSectionIntegration::HSSSectionIntegration(double H, double B, double T, + int NFH, int NFB, int NFT): + SectionIntegration(SECTION_INTEGRATION_TAG_HSS), + h(H), b(B), t(T), Nfh(NFH), Nfb(NFB), Nft(1), parameterID(0) +{ + +} + +HSSSectionIntegration::HSSSectionIntegration(): + SectionIntegration(SECTION_INTEGRATION_TAG_HSS), + h(0.0), b(0.0), t(0.0), Nfh(0), Nfb(0), Nft(0), parameterID(0) +{ + +} + +HSSSectionIntegration::~HSSSectionIntegration() +{ + +} + +int +HSSSectionIntegration::getNumFibers(FiberType type) +{ + return 2*Nfh + 2*Nfb + 4*Nft*Nft; +} + +int +HSSSectionIntegration::arrangeFibers(UniaxialMaterial **theMaterials, + UniaxialMaterial *theSteel) +{ + int numFibers = this->getNumFibers(); + + for (int i = 0; i < numFibers; i++) + theMaterials[i] = theSteel; + + return 0; +} + +int +HSSSectionIntegration::arrangeFibers(NDMaterial **theMaterials, + NDMaterial *theSteel) +{ + int numFibers = this->getNumFibers(); + + for (int i = 0; i < numFibers; i++) + theMaterials[i] = theSteel; + + return 0; +} + +void +HSSSectionIntegration::getFiberLocations(int nFibers, double *yi, double *zi) +{ + double hw = h - 2*t; + double bw = b - 2*t; + + int i, loc; + + double yIncr = hw/Nfh; + double yStart = 0.5 * (hw-yIncr); + + for (loc = 0, i = 0; loc < Nfh; loc++, i++) { + yi[loc] = yStart - yIncr*i; + yi[loc+Nfh] = yi[loc]; + } + if (zi != 0) { + for (loc = 0; loc < Nfh; loc++) { + zi[loc] = 0.5*(bw + t); + zi[loc+Nfh] = -zi[loc]; + } + } + + for (loc = 2*Nfh, i = 0; i < Nfb; loc++, i++) { + yi[loc] = 0.5*(hw + t); + yi[loc+Nfb] = -yi[loc]; + } + if (zi != 0) { + double zIncr = bw/Nfb; + double zStart = 0.5 * (bw-zIncr); + for (loc = 2*Nfh, i = 0; i < Nfb; loc++, i++) { + zi[loc] = zStart - zIncr*i; + zi[loc+Nfb] = zi[loc]; + } + } + + // Four corners + loc = 2*Nfh + 2*Nfb; + double yy = 0.5*(hw+t); + yi[loc++] = yy; + yi[loc++] = -yy; + yi[loc++] = -yy; + yi[loc++] = yy; + if (zi != 0) { + loc = 2*Nfh + 2*Nfb; + double zz = 0.5*(bw+t); + zi[loc++] = zz; + zi[loc++] = zz; + zi[loc++] = -zz; + zi[loc++] = -zz; + } + + return; +} + +void +HSSSectionIntegration::getFiberWeights(int nFibers, double *wt) +{ + double hw = h - 2*t; + double bw = b - 2*t; + + // Assuming Nft = 1 + double a_h = hw*t/(Nfh); + double a_b = bw*t/(Nfb); + + int i, loc; + + for (loc = 0; loc < Nfh; loc++) { + wt[loc] = a_h; + wt[loc+Nfh] = a_h; + } + + for (loc = 2*Nfh, i = 0; i < Nfb; i++, loc++) { + wt[loc] = a_b; + wt[loc+Nfb] = a_b; + } + + // Four corners + double a_t = t*t; + loc = 2*Nfh + 2*Nfb; + for (i = 0; i < 4; i++, loc++) { + wt[loc] = a_t; + } + + return; +} + +SectionIntegration* +HSSSectionIntegration::getCopy(void) +{ + return new HSSSectionIntegration(h, b, t, Nfh, Nfb, Nft); +} + +int +HSSSectionIntegration::setParameter(const char **argv, int argc, + Parameter ¶m) +{ + if (argc < 1) + return -1; + + if (strcmp(argv[0],"h") == 0) { + param.setValue(h); + return param.addObject(1, this); + } + if (strcmp(argv[0],"b") == 0) { + param.setValue(b); + return param.addObject(2, this); + } + if (strcmp(argv[0],"t") == 0) { + param.setValue(t); + return param.addObject(3, this); + } + + return -1; +} + +int +HSSSectionIntegration::updateParameter(int parameterID, + Information &info) +{ + switch (parameterID) { + case 1: + h = info.theDouble; + return 0; + case 2: + b = info.theDouble; + return 0; + case 3: + t = info.theDouble; + return 0; + default: + return -1; + } +} + +int +HSSSectionIntegration::activateParameter(int paramID) +{ + parameterID = paramID; + + return 0; +} + +void +HSSSectionIntegration::getLocationsDeriv(int nFibers, double *dyidh, double *dzidh) +{ + /* + //double dw = d-2*tf; + + double dddh = 0.0; + double ddwdh = 0.0; + //double dtwdh = 0.0; + //double dbfdh = 0.0; + double dtfdh = 0.0; + + if (parameterID == 1) { // d + dddh = 1.0; + ddwdh = 1.0; + } + //if (parameterID == 2) // tw + // dtwdh = 1.0; + //if (parameterID == 3) // bf + // dbfdh = 1.0; + if (parameterID == 4) { // tf + dtfdh = 1.0; + ddwdh = -2.0; + } + + //double yIncr = tf/Nftf; + //double yStart = 0.5*d - 0.5*yIncr; + double dyIncrdh = dtfdh/Nftf; + double dyStartdh = 0.5 * (dddh-dyIncrdh); + + int loc; + + for (loc = 0; loc < Nftf; loc++) { + //yi[loc] = yStart - yIncr*loc; + //yi[nFibers-loc-1] = -yi[loc]; + dyidh[loc] = dyStartdh - dyIncrdh*loc; + dyidh[nFibers-loc-1] = -dyidh[loc]; + } + + //yIncr = dw/Nfdw; + //yStart = 0.5*dw - 0.5*yIncr; + dyIncrdh = ddwdh/Nfdw; + dyStartdh = 0.5 * (ddwdh-dyIncrdh); + + for (int count = 0; loc < nFibers-Nftf; loc++, count++) { + //yi[loc] = yStart - yIncr*count; + dyidh[loc] = dyStartdh - dyIncrdh*count; + } + + if (dzidh != 0) { + for (int i = 0; i < nFibers; i++) + dzidh[i] = 0.0; + } + */ + + return; +} + +void +HSSSectionIntegration::getWeightsDeriv(int nFibers, double *dwtdh) +{ + /* + double dw = d-2*tf; + + double ddwdh = 0.0; + double dtwdh = 0.0; + double dbfdh = 0.0; + double dtfdh = 0.0; + + if (parameterID == 1) // d + ddwdh = 1.0; + if (parameterID == 2) // tw + dtwdh = 1.0; + if (parameterID == 3) // bf + dbfdh = 1.0; + if (parameterID == 4) { // tf + dtfdh = 1.0; + ddwdh = -2.0; + } + + double dAfdh = (bf*dtfdh + dbfdh*tf) / Nftf; + double dAwdh = (dw*dtwdh + ddwdh*tw) / Nfdw; + + int loc = 0; + + for (loc = 0; loc < Nftf; loc++) { + dwtdh[loc] = dAfdh; + dwtdh[nFibers-loc-1] = dAfdh; + } + + for ( ; loc < nFibers-Nftf; loc++) + dwtdh[loc] = dAwdh; + + //for (int i = 0; i < nFibers; i++) + // opserr << dwtsdh[i] << ' '; + //opserr << endln; + */ + + return; +} + +void +HSSSectionIntegration::Print(OPS_Stream &s, int flag) +{ + s << "HSS" << endln; + s << " h = " << h << endln; + s << " b = " << b << endln; + s << " t = " << t << endln; + s << " Nfh = " << Nfh << endln; + s << " Nfb = " << Nfb << endln; + s << " Nft = " << Nft << endln; + + return; +} + +int +HSSSectionIntegration::sendSelf(int cTag, Channel &theChannel) +{ + static Vector data(6); + + data(0) = h; + data(1) = b; + data(2) = t; + data(3) = Nfh; + data(4) = Nfb; + data(5) = Nft; + + int dbTag = this->getDbTag(); + + if (theChannel.sendVector(dbTag, cTag, data) < 0) { + opserr << "HSSSectionIntegration::sendSelf() - failed to send Vector data\n"; + return -1; + } + + return 0; +} + +int +HSSSectionIntegration::recvSelf(int cTag, Channel &theChannel, + FEM_ObjectBroker &theBroker) +{ + static Vector data(6); + + int dbTag = this->getDbTag(); + + if (theChannel.recvVector(dbTag, cTag, data) < 0) { + opserr << "HSSSectionIntegration::recvSelf() - failed to receive Vector data\n"; + return -1; + } + + h = data(0); + b = data(1); + t = data(2); + Nfh = (int)data(3); + Nfb = (int)data(4); + Nft = (int)data(5); + + return 0; +} diff --git a/SRC/material/section/integration/HSSSectionIntegration.h b/SRC/material/section/integration/HSSSectionIntegration.h new file mode 100644 index 000000000..02dd045a5 --- /dev/null +++ b/SRC/material/section/integration/HSSSectionIntegration.h @@ -0,0 +1,77 @@ +/* ****************************************************************** ** +** OpenSees - Open System for Earthquake Engineering Simulation ** +** Pacific Earthquake Engineering Research Center ** +** ** +** ** +** (C) Copyright 1999, The Regents of the University of California ** +** All Rights Reserved. ** +** ** +** Commercial use of this program without express permission of the ** +** University of California, Berkeley, is strictly prohibited. See ** +** file 'COPYRIGHT' in main directory for information on usage and ** +** redistribution, and for a DISCLAIMER OF ALL WARRANTIES. ** +** ** +** Developed by: ** +** Frank McKenna (fmckenna@ce.berkeley.edu) ** +** Gregory L. Fenves (fenves@ce.berkeley.edu) ** +** Filip C. Filippou (filippou@ce.berkeley.edu) ** +** ** +** ****************************************************************** */ + +// $Revision: 1.5 $ +// $Date: 2010-09-13 21:31:26 $ +// $Source: /usr/local/cvs/OpenSees/SRC/material/section/integration/HSSSectionIntegration.h,v $ + +#ifndef HSSSectionIntegration_h +#define HSSSectionIntegration_h + +#include + +class UniaxialMaterial; +class NDMaterial; + +class HSSSectionIntegration : public SectionIntegration +{ + public: + HSSSectionIntegration(double h, double b, double t, + int Nfh, int Nfb, int Nft = 1); + HSSSectionIntegration(); + ~HSSSectionIntegration(); + + int getNumFibers(FiberType type = all); + + void getFiberLocations(int nFibers, double *yi, double *zi = 0); + void getFiberWeights(int nFibers, double *wt); + + SectionIntegration *getCopy(void); + + int sendSelf(int cTag, Channel &theChannel); + int recvSelf(int cTag, Channel &theChannel, FEM_ObjectBroker &theBroker); + + int setParameter(const char **argv, int argc, Parameter ¶m); + int updateParameter(int parameterID, Information &info); + int activateParameter(int parameterID); + + void getLocationsDeriv(int nFibers, double *dyidh, double *dzidh = 0); + void getWeightsDeriv(int nFibers, double *dwtdh); + + void Print(OPS_Stream &s, int flag = 0); + + int arrangeFibers(UniaxialMaterial **theMaterials, + UniaxialMaterial *theSteel); + int arrangeFibers(NDMaterial **theMaterials, + NDMaterial *theSteel); + + private: + double h; + double b; + double t; + + int Nfh; + int Nfb; + int Nft; + + int parameterID; +}; + +#endif diff --git a/SRC/material/section/integration/Makefile b/SRC/material/section/integration/Makefile index c492ccd20..4f3190b63 100644 --- a/SRC/material/section/integration/Makefile +++ b/SRC/material/section/integration/Makefile @@ -3,7 +3,8 @@ include ../../../../Makefile.def OBJS = SectionIntegration.o \ WideFlangeSectionIntegration.o RCSectionIntegration.o \ RCTBeamSectionIntegration.o RCCircularSectionIntegration.o \ - TubeSectionIntegration.o RCTunnelSectionIntegration.o + TubeSectionIntegration.o RCTunnelSectionIntegration.o \ + HSSSectionIntegration.o # Compilation control diff --git a/SRC/material/uniaxial/DamperMaterial.cpp b/SRC/material/uniaxial/DamperMaterial.cpp index 6ecba2cd0..27e20183d 100644 --- a/SRC/material/uniaxial/DamperMaterial.cpp +++ b/SRC/material/uniaxial/DamperMaterial.cpp @@ -88,6 +88,9 @@ DamperMaterial::DamperMaterial(int tag, trialStrain(0.0), trialStrainRate(0.0), theMaterial(0) { theMaterial = theMaterialModel->getCopy(); + + if (theMaterial == 0) + opserr << "DamperMaterial::DamperMaterial -- failed to get copy of material\n"; } @@ -105,7 +108,8 @@ DamperMaterial::DamperMaterial() DamperMaterial::~DamperMaterial() { - delete theMaterial; + if (theMaterial) + delete theMaterial; } @@ -118,9 +122,10 @@ DamperMaterial::setTrialStrain(double strain, double strainRate) trialStrain = strain; trialStrainRate = strainRate; - theMaterial->setTrialStrain(strainRate, 0); - - return 0; + if (theMaterial) + return theMaterial->setTrialStrain(strainRate, 0); + else + return -1; } @@ -139,39 +144,49 @@ DamperMaterial::getStrainRate(void) double DamperMaterial::getStress(void) { - return theMaterial->getStress(); + if (theMaterial) + return theMaterial->getStress(); + else + return 0.0; } - - double DamperMaterial::getTangent(void) { - return 0; + return 0.0; } double DamperMaterial::getInitialTangent(void) { - return 0; + return 0.0; } double DamperMaterial::getDampTangent(void) { - return theMaterial->getTangent(); + if (theMaterial) + return theMaterial->getTangent(); + else + return 0.0; } int DamperMaterial::commitState(void) { - return theMaterial->commitState(); + if (theMaterial) + return theMaterial->commitState(); + else + return -1; } int DamperMaterial::revertToLastCommit(void) { - return theMaterial->revertToLastCommit(); + if (theMaterial) + return theMaterial->revertToLastCommit(); + else + return -1; } @@ -180,9 +195,11 @@ DamperMaterial::revertToStart(void) { trialStrain = 0.0; trialStrainRate = 0.0; - theMaterial->revertToStart(); - - return 0; + + if (theMaterial) + return theMaterial->revertToStart(); + else + return -1; } @@ -190,11 +207,13 @@ DamperMaterial::revertToStart(void) UniaxialMaterial * DamperMaterial::getCopy(void) { - DamperMaterial *theCopy = new - DamperMaterial(this->getTag(), theMaterial); + DamperMaterial *theCopy = 0; + if (theMaterial) { + theCopy = new DamperMaterial(this->getTag(), theMaterial); - theCopy->trialStrain = trialStrain; - theCopy->trialStrainRate = trialStrainRate; + theCopy->trialStrain = trialStrain; + theCopy->trialStrainRate = trialStrainRate; + } return theCopy; } @@ -203,6 +222,11 @@ DamperMaterial::getCopy(void) int DamperMaterial::sendSelf(int cTag, Channel &theChannel) { + if (theMaterial == 0) { + opserr << "DamperMaterial::sendSelf() - theMaterial is null, nothing to send\n"; + return -1; + } + int res = 0; static ID data(3); @@ -253,7 +277,8 @@ DamperMaterial::recvSelf(int cTag, Channel &theChannel, if (theMaterial == 0) { opserr << "FATAL DamperMaterial::recvSelf() "; opserr << " could not get a UniaxialMaterial \n"; - exit(-1); + //exit(-1); + return -1; } theMaterial->setDbTag(dbTag); theMaterial->recvSelf(cTag, theChannel, theBroker); @@ -264,11 +289,11 @@ DamperMaterial::recvSelf(int cTag, Channel &theChannel, void DamperMaterial::Print(OPS_Stream &s, int flag) { - if (flag == OPS_PRINT_PRINTMODEL_MATERIAL) { - s << "DamperMaterial tag: " << this->getTag() << endln; - s << " "; - theMaterial->Print(s, flag); - } + s << "DamperMaterial tag: " << this->getTag() << endln; + if (theMaterial) + s << "\tMaterial: " << theMaterial->getTag() << endln; + else + s << "\tMaterial is NULL" << endln; } diff --git a/SRC/material/uniaxial/HyperbolicGapMaterial.cpp b/SRC/material/uniaxial/HyperbolicGapMaterial.cpp index a7b355316..ea5d50cce 100644 --- a/SRC/material/uniaxial/HyperbolicGapMaterial.cpp +++ b/SRC/material/uniaxial/HyperbolicGapMaterial.cpp @@ -93,12 +93,14 @@ HyperbolicGapMaterial::HyperbolicGapMaterial(int tag, double kmax, double kur, d Kmax(kmax), Kur(kur), Rf(rf), Fult(fult), gap(gap0) { if (gap>=0) { - opserr << "HyperbolicGapMaterial::HyperbolicGapMaterial -- Initial gap size must be negative for compression-only material\n"; - exit(-1); + opserr << "HyperbolicGapMaterial::HyperbolicGapMaterial -- Initial gap size must be negative for compression-only material, setting to negative\n"; + //exit(-1); + gap = -gap; } if (Fult>0) { - opserr << "HyperbolicGapMaterial::HyperbolicGapMaterial -- Fult must be negative for compression-only material\n"; - exit(-1); + opserr << "HyperbolicGapMaterial::HyperbolicGapMaterial -- Fult must be negative for compression-only material, setting to negative\n"; + //exit(-1); + Fult = -Fult; } if (Kmax == 0.0) { opserr << "HyperbolicGapMaterial::HyperbolicGapMaterial -- Kmax is zero, continuing with Kmax = Fult/0.002\n"; diff --git a/SRC/material/uniaxial/InitStrainMaterial.cpp b/SRC/material/uniaxial/InitStrainMaterial.cpp index 6c663793b..2d5b304d7 100644 --- a/SRC/material/uniaxial/InitStrainMaterial.cpp +++ b/SRC/material/uniaxial/InitStrainMaterial.cpp @@ -88,10 +88,11 @@ InitStrainMaterial::InitStrainMaterial(int tag, if (theMaterial == 0) { opserr << "InitStrainMaterial::InitStrainMaterial -- failed to get copy of material\n"; - exit(-1); + //exit(-1); + } else { + theMaterial->setTrialStrain(epsInit); + theMaterial->commitState(); } - theMaterial->setTrialStrain(epsInit); - theMaterial->commitState(); } InitStrainMaterial::InitStrainMaterial() @@ -111,67 +112,95 @@ int InitStrainMaterial::setTrialStrain(double strain, double strainRate) { localStrain = strain; - - return theMaterial->setTrialStrain(strain+epsInit, strainRate); + + if (theMaterial) + return theMaterial->setTrialStrain(strain+epsInit, strainRate); + else + return -1; } double InitStrainMaterial::getStress(void) { - return theMaterial->getStress(); + if (theMaterial) + return theMaterial->getStress(); + else + return 0.0; } double InitStrainMaterial::getTangent(void) { - return theMaterial->getTangent(); + if (theMaterial) + return theMaterial->getTangent(); + else + return 0.0; } double InitStrainMaterial::getDampTangent(void) { - return theMaterial->getDampTangent(); + if (theMaterial) + return theMaterial->getDampTangent(); + else + return 0.0; } double InitStrainMaterial::getStrain(void) { - return theMaterial->getStrain(); + if (theMaterial) + return theMaterial->getStrain(); + else + return 0.0; } double InitStrainMaterial::getStrainRate(void) { - return theMaterial->getStrainRate(); + if (theMaterial) + return theMaterial->getStrainRate(); + else + return 0.0; } int InitStrainMaterial::commitState(void) -{ - return theMaterial->commitState(); +{ + if (theMaterial) + return theMaterial->commitState(); + else + return -1; } int InitStrainMaterial::revertToLastCommit(void) { - return theMaterial->revertToLastCommit(); + if (theMaterial) + return theMaterial->revertToLastCommit(); + else + return -1; } int InitStrainMaterial::revertToStart(void) { int res = 0; - res = theMaterial->revertToStart(); - res += theMaterial->setTrialStrain(epsInit); - res += theMaterial->commitState(); - return res; + if (theMaterial) { + res = theMaterial->revertToStart(); + res += theMaterial->setTrialStrain(epsInit); + res += theMaterial->commitState(); + return res; + } else + return -1; } UniaxialMaterial * InitStrainMaterial::getCopy(void) { - InitStrainMaterial *theCopy = - new InitStrainMaterial(this->getTag(), *theMaterial, epsInit); + InitStrainMaterial *theCopy = 0; + if (theMaterial) + theCopy = new InitStrainMaterial(this->getTag(), *theMaterial, epsInit); return theCopy; } @@ -179,6 +208,11 @@ InitStrainMaterial::getCopy(void) int InitStrainMaterial::sendSelf(int cTag, Channel &theChannel) { + if (theMaterial == 0) { + opserr << "InitStrainMaterial::sendSelf() - theMaterial is null, nothing to send\n"; + return -1; + } + int dbTag = this->getDbTag(); static ID dataID(3); @@ -260,11 +294,17 @@ InitStrainMaterial::Print(OPS_Stream &s, int flag) s << "\t\t\t{"; s << "\"name\": \"" << this->getTag() << "\", "; s << "\"type\": \"InitStrainMaterial\", "; - s << "\"Material\": " << theMaterial->getTag() << ", "; + if (theMaterial) + s << "\"Material\": " << theMaterial->getTag() << ", "; + else + s << "\"Material\": " << "NULL" << ", "; s << "\"initialStrain\": " << epsInit << "}"; } else { s << "InitStrainMaterial tag: " << this->getTag() << endln; - s << "\tMaterial: " << theMaterial->getTag() << endln; + if (theMaterial) + s << "\tMaterial: " << theMaterial->getTag() << endln; + else + s << "\tMaterial is NULL" << endln; s << "\tinitital strain: " << epsInit << endln; } } @@ -278,7 +318,10 @@ InitStrainMaterial::setParameter(const char **argv, int argc, Parameter ¶m) } // Otherwise, pass it on to the wrapped material - return theMaterial->setParameter(argv, argc, param); + if (theMaterial) + return theMaterial->setParameter(argv, argc, param); + else + return -1; } int @@ -286,8 +329,11 @@ InitStrainMaterial::updateParameter(int parameterID, Information &info) { if (parameterID == 1) { this->epsInit = info.theDouble; - theMaterial->setTrialStrain(localStrain+epsInit); - theMaterial->commitState(); + if (theMaterial) { + theMaterial->setTrialStrain(localStrain+epsInit); + theMaterial->commitState(); + } else + return -1; } return 0; @@ -296,18 +342,27 @@ InitStrainMaterial::updateParameter(int parameterID, Information &info) double InitStrainMaterial::getStressSensitivity(int gradIndex, bool conditional) { - return theMaterial->getStressSensitivity(gradIndex, conditional); + if (theMaterial) + return theMaterial->getStressSensitivity(gradIndex, conditional); + else + return 0.0; } double InitStrainMaterial::getInitialTangentSensitivity(int gradIndex) { - return theMaterial->getInitialTangentSensitivity(gradIndex); + if (theMaterial) + return theMaterial->getInitialTangentSensitivity(gradIndex); + else + return 0.0; } int InitStrainMaterial::commitSensitivity(double strainGradient, int gradIndex, int numGrads) { - return theMaterial->commitSensitivity(strainGradient, gradIndex, numGrads); + if (theMaterial) + return theMaterial->commitSensitivity(strainGradient, gradIndex, numGrads); + else + return -1; } diff --git a/SRC/material/uniaxial/MultiplierMaterial.cpp b/SRC/material/uniaxial/MultiplierMaterial.cpp index 3af29b7c4..075a64937 100644 --- a/SRC/material/uniaxial/MultiplierMaterial.cpp +++ b/SRC/material/uniaxial/MultiplierMaterial.cpp @@ -99,7 +99,7 @@ MultiplierMaterial::MultiplierMaterial(int tag, UniaxialMaterial &material, doub if (theMaterial == 0) { opserr << "MultiplierMaterial::MultiplierMaterial -- failed to get copy of material\n"; - exit(-1); + //exit(-1); } } @@ -118,33 +118,48 @@ MultiplierMaterial::~MultiplierMaterial() int MultiplierMaterial::setTrialStrain(double strain, double strainRate) { - return theMaterial->setTrialStrain(strain, strainRate); + if (theMaterial) + return theMaterial->setTrialStrain(strain, strainRate); + else + return -1; } int MultiplierMaterial::setTrialStrain(double strain, double temp, double strainRate) { - return theMaterial->setTrialStrain(strain, temp, strainRate); + if (theMaterial) + return theMaterial->setTrialStrain(strain, temp, strainRate); + else + return -1; } double MultiplierMaterial::getStress(void) { - return multiplier*theMaterial->getStress(); + if (theMaterial) + return multiplier*theMaterial->getStress(); + else + return 0.0; } double MultiplierMaterial::getTangent(void) { - return multiplier*theMaterial->getTangent(); + if (theMaterial) + return multiplier*theMaterial->getTangent(); + else + return 0.0; } double MultiplierMaterial::getDampTangent(void) { - return multiplier*theMaterial->getDampTangent(); + if (theMaterial) + return multiplier*theMaterial->getDampTangent(); + else + return 0.0; } @@ -152,38 +167,54 @@ MultiplierMaterial::getDampTangent(void) double MultiplierMaterial::getStrain(void) { - return theMaterial->getStrain(); + if (theMaterial) + return theMaterial->getStrain(); + else + return 0.0; } double MultiplierMaterial::getStrainRate(void) { - return theMaterial->getStrainRate(); + if (theMaterial) + return theMaterial->getStrainRate(); + else + return 0.0; } int MultiplierMaterial::commitState(void) -{ - return theMaterial->commitState(); +{ + if (theMaterial) + return theMaterial->commitState(); + else + return -1; } int MultiplierMaterial::revertToLastCommit(void) { - return theMaterial->revertToLastCommit(); + if (theMaterial) + return theMaterial->revertToLastCommit(); + else + return -1; } int MultiplierMaterial::revertToStart(void) { - return theMaterial->revertToStart(); + if (theMaterial) + return theMaterial->revertToStart(); + else + return -1; } UniaxialMaterial * MultiplierMaterial::getCopy(void) { - MultiplierMaterial *theCopy = - new MultiplierMaterial(this->getTag(), *theMaterial, multiplier); + MultiplierMaterial *theCopy = 0; + if (theMaterial) + theCopy = new MultiplierMaterial(this->getTag(), *theMaterial, multiplier); return theCopy; } @@ -191,6 +222,11 @@ MultiplierMaterial::getCopy(void) int MultiplierMaterial::sendSelf(int cTag, Channel &theChannel) { + if (theMaterial == 0) { + opserr << "MultiplierMaterial::sendSelf() - theMaterial is null, nothing to send\n"; + return -1; + } + int dbTag = this->getDbTag(); static ID dataID(3); @@ -267,7 +303,10 @@ void MultiplierMaterial::Print(OPS_Stream &s, int flag) { s << "MultiplierMaterial tag: " << this->getTag() << endln; - s << "\tMaterial: " << theMaterial->getTag() << endln; + if (theMaterial) + s << "\tMaterial: " << theMaterial->getTag() << endln; + else + s << "\tMaterial is NULL" << endln; s << "\tMultiplier: " << multiplier << endln; } @@ -278,7 +317,10 @@ MultiplierMaterial::setParameter(const char **argv, int argc, Parameter ¶m) param.setValue(multiplier); return param.addObject(1,this); } - return theMaterial->setParameter(argv, argc, param); + if (theMaterial) + return theMaterial->setParameter(argv, argc, param); + else + return -1; } int @@ -306,6 +348,9 @@ MultiplierMaterial::activateParameter(int paramID) double MultiplierMaterial::getStressSensitivity(int gradIndex, bool conditional) { + if (theMaterial == 0) + return 0.0; + // dsig = dF*sigma + F*dsigma if (parameterID == 1) return theMaterial->getStress(); // dF*sigma where dF=1 @@ -316,12 +361,18 @@ MultiplierMaterial::getStressSensitivity(int gradIndex, bool conditional) double MultiplierMaterial::getStrainSensitivity(int gradIndex) { - return theMaterial->getStrainSensitivity(gradIndex); + if (theMaterial) + return theMaterial->getStrainSensitivity(gradIndex); + else + return 0.0; } double MultiplierMaterial::getInitialTangentSensitivity(int gradIndex) { + if (theMaterial == 0) + return 0.0; + if (parameterID == 1) return theMaterial->getInitialTangent(); else @@ -331,6 +382,9 @@ MultiplierMaterial::getInitialTangentSensitivity(int gradIndex) double MultiplierMaterial::getDampTangentSensitivity(int gradIndex) { + if (theMaterial == 0) + return 0.0; + if (parameterID == 1) return theMaterial->getDampTangent(); else @@ -340,11 +394,17 @@ MultiplierMaterial::getDampTangentSensitivity(int gradIndex) double MultiplierMaterial::getRhoSensitivity(int gradIndex) { - return theMaterial->getRhoSensitivity(gradIndex); + if (theMaterial) + return theMaterial->getRhoSensitivity(gradIndex); + else + return 0.0; } int MultiplierMaterial::commitSensitivity(double strainGradient, int gradIndex, int numGrads) { - return theMaterial->commitSensitivity(strainGradient, gradIndex, numGrads); + if (theMaterial) + return theMaterial->commitSensitivity(strainGradient, gradIndex, numGrads); + else + return -1; } diff --git a/SRC/material/uniaxial/PathIndependentMaterial.cpp b/SRC/material/uniaxial/PathIndependentMaterial.cpp index 50c17b598..342beb83b 100644 --- a/SRC/material/uniaxial/PathIndependentMaterial.cpp +++ b/SRC/material/uniaxial/PathIndependentMaterial.cpp @@ -79,7 +79,7 @@ PathIndependentMaterial::PathIndependentMaterial(int tag, UniaxialMaterial &mate if (theMaterial == 0) { opserr << "PathIndependentMaterial::PathIndependentMaterial -- failed to get copy of material\n"; - exit(-1); + //exit(-1); } } @@ -98,44 +98,65 @@ PathIndependentMaterial::~PathIndependentMaterial() int PathIndependentMaterial::setTrialStrain(double strain, double strainRate) { + if (theMaterial) return theMaterial->setTrialStrain(strain, strainRate); + else + return -1; } double PathIndependentMaterial::getStress(void) { + if (theMaterial) return theMaterial->getStress(); + else + return 0.0; } double PathIndependentMaterial::getTangent(void) { + if (theMaterial) return theMaterial->getTangent(); + else + return 0.0; } double PathIndependentMaterial::getDampTangent(void) { + if (theMaterial) return theMaterial->getDampTangent(); + else + return 0.0; } double PathIndependentMaterial::getInitialTangent(void) { + if (theMaterial) return theMaterial->getInitialTangent(); + else + return 0.0; } double PathIndependentMaterial::getStrain(void) { + if (theMaterial) return theMaterial->getStrain(); + else + return 0.0; } double PathIndependentMaterial::getStrainRate(void) { + if (theMaterial) return theMaterial->getStrainRate(); + else + return 0.0; } int @@ -153,22 +174,31 @@ PathIndependentMaterial::revertToLastCommit(void) int PathIndependentMaterial::revertToStart(void) { + if (theMaterial) return theMaterial->revertToStart(); + else + return -1; } UniaxialMaterial * PathIndependentMaterial::getCopy(void) { - PathIndependentMaterial *theCopy = - new PathIndependentMaterial(this->getTag(), *theMaterial); + PathIndependentMaterial *theCopy = 0; + if (theMaterial) + theCopy = new PathIndependentMaterial(this->getTag(), *theMaterial); - return theCopy; + return theCopy; } int PathIndependentMaterial::sendSelf(int cTag, Channel &theChannel) { + if (theMaterial == 0) { + opserr << "PathIndependentMaterial::sendSelf() - theMaterial is null, nothing to send\n"; + return -1; + } + int res = 0; static ID classTags(3); @@ -235,7 +265,8 @@ PathIndependentMaterial::recvSelf(int cTag, Channel &theChannel, theMaterial = theBroker.getNewUniaxialMaterial(classTags(0)); if (theMaterial == 0) { opserr << "PathIndependentMaterial::recvSelf -- could not get a UniaxialMaterial\n"; - exit(-1); + //exit(-1); + return -1; } } @@ -254,13 +285,19 @@ void PathIndependentMaterial::Print(OPS_Stream &s, int flag) { s << "PathIndependentMaterial tag: " << this->getTag() << endln; - s << "\tmaterial: " << theMaterial->getTag() << endln; + if (theMaterial) + s << "\tMaterial: " << theMaterial->getTag() << endln; + else + s << "\tMaterial is NULL" << endln; } int PathIndependentMaterial::setParameter(const char **argv, int argc, Parameter ¶m) { - return theMaterial->setParameter(argv, argc, param); + if (theMaterial) + return theMaterial->setParameter(argv, argc, param); + else + return -1; } int @@ -272,31 +309,46 @@ PathIndependentMaterial::updateParameter(int parameterID, Information &info) double PathIndependentMaterial::getStressSensitivity(int gradIndex, bool conditional) { - return theMaterial->getStressSensitivity(gradIndex, conditional); + if (theMaterial) + return theMaterial->getStressSensitivity(gradIndex, conditional); + else + return 0.0; } double PathIndependentMaterial::getStrainSensitivity(int gradIndex) { - return theMaterial->getStrainSensitivity(gradIndex); + if (theMaterial) + return theMaterial->getStrainSensitivity(gradIndex); + else + return 0.0; } double PathIndependentMaterial::getInitialTangentSensitivity(int gradIndex) { - return theMaterial->getInitialTangentSensitivity(gradIndex); + if (theMaterial) + return theMaterial->getInitialTangentSensitivity(gradIndex); + else + return 0.0; } double PathIndependentMaterial::getDampTangentSensitivity(int gradIndex) { - return theMaterial->getDampTangentSensitivity(gradIndex); + if (theMaterial) + return theMaterial->getDampTangentSensitivity(gradIndex); + else + return 0.0; } double PathIndependentMaterial::getRhoSensitivity(int gradIndex) { - return theMaterial->getRhoSensitivity(gradIndex); + if (theMaterial) + return theMaterial->getRhoSensitivity(gradIndex); + else + return 0.0; } int diff --git a/SRC/material/uniaxial/PenaltyMaterial.cpp b/SRC/material/uniaxial/PenaltyMaterial.cpp index 3d6b2f1a1..faa9f3501 100644 --- a/SRC/material/uniaxial/PenaltyMaterial.cpp +++ b/SRC/material/uniaxial/PenaltyMaterial.cpp @@ -98,7 +98,7 @@ PenaltyMaterial::PenaltyMaterial(int tag, UniaxialMaterial &material, double mul if (theMaterial == 0) { opserr << "PenaltyMaterial::PenaltyMaterial -- failed to get copy of material\n"; - exit(-1); + //exit(-1); } } @@ -117,72 +117,101 @@ PenaltyMaterial::~PenaltyMaterial() int PenaltyMaterial::setTrialStrain(double strain, double strainRate) { - return theMaterial->setTrialStrain(strain, strainRate); + if (theMaterial) + return theMaterial->setTrialStrain(strain, strainRate); + else + return -1; } int PenaltyMaterial::setTrialStrain(double strain, double temp, double strainRate) { - return theMaterial->setTrialStrain(strain, temp, strainRate); + if (theMaterial) + return theMaterial->setTrialStrain(strain, temp, strainRate); + else + return -1; } double PenaltyMaterial::getStress(void) { - return theMaterial->getStress() + penalty*theMaterial->getStrain(); + if (theMaterial) + return theMaterial->getStress() + penalty*theMaterial->getStrain(); + else + return 0.0; } double PenaltyMaterial::getTangent(void) { - return theMaterial->getTangent() + penalty; + if (theMaterial) + return theMaterial->getTangent() + penalty; + else + return 0.0; } double PenaltyMaterial::getDampTangent(void) { - return theMaterial->getDampTangent(); + if (theMaterial) + return theMaterial->getDampTangent(); + else + return 0.0; } - - double PenaltyMaterial::getStrain(void) { - return theMaterial->getStrain(); + if (theMaterial) + return theMaterial->getStrain(); + else + return 0.0; } double PenaltyMaterial::getStrainRate(void) { - return theMaterial->getStrainRate(); + if (theMaterial) + return theMaterial->getStrainRate(); + else + return 0.0; } int PenaltyMaterial::commitState(void) -{ - return theMaterial->commitState(); +{ + if (theMaterial) + return theMaterial->commitState(); + else + return -1; } int PenaltyMaterial::revertToLastCommit(void) { - return theMaterial->revertToLastCommit(); + if (theMaterial) + return theMaterial->revertToLastCommit(); + else + return -1; } int PenaltyMaterial::revertToStart(void) { - return theMaterial->revertToStart(); + if (theMaterial) + return theMaterial->revertToStart(); + else + return -1; } UniaxialMaterial * PenaltyMaterial::getCopy(void) { - PenaltyMaterial *theCopy = - new PenaltyMaterial(this->getTag(), *theMaterial, penalty); + PenaltyMaterial *theCopy = 0; + if (theMaterial) + theCopy = new PenaltyMaterial(this->getTag(), *theMaterial, penalty); return theCopy; } @@ -190,6 +219,11 @@ PenaltyMaterial::getCopy(void) int PenaltyMaterial::sendSelf(int cTag, Channel &theChannel) { + if (theMaterial == 0) { + opserr << "PenaltyMaterial::sendSelf() - theMaterial is null, nothing to send\n"; + return -1; + } + int dbTag = this->getDbTag(); static ID dataID(3); @@ -266,7 +300,10 @@ void PenaltyMaterial::Print(OPS_Stream &s, int flag) { s << "PenaltyMaterial tag: " << this->getTag() << endln; - s << "\tMaterial: " << theMaterial->getTag() << endln; + if (theMaterial) + s << "\tMaterial: " << theMaterial->getTag() << endln; + else + s << "\tMaterial is NULL" << endln; s << "\tPenalty: " << penalty << endln; } @@ -277,7 +314,10 @@ PenaltyMaterial::setParameter(const char **argv, int argc, Parameter ¶m) param.setValue(penalty); return param.addObject(1,this); } - return theMaterial->setParameter(argv, argc, param); + if (theMaterial) + return theMaterial->setParameter(argv, argc, param); + else + return -1; } int @@ -305,6 +345,9 @@ PenaltyMaterial::activateParameter(int paramID) double PenaltyMaterial::getStressSensitivity(int gradIndex, bool conditional) { + if (theMaterial == 0) + return 0.0; + // dsig = dsigma + dalpha*strain < + alpha*dstrain> if (parameterID == 1) return theMaterial->getStrain(); // dalpha*strain where dalpha=1 @@ -315,7 +358,10 @@ PenaltyMaterial::getStressSensitivity(int gradIndex, bool conditional) double PenaltyMaterial::getStrainSensitivity(int gradIndex) { - return theMaterial->getStrainSensitivity(gradIndex); + if (theMaterial) + return theMaterial->getStrainSensitivity(gradIndex); + else + return 0.0; } double @@ -323,24 +369,37 @@ PenaltyMaterial::getInitialTangentSensitivity(int gradIndex) { if (parameterID == 1) return 1.0; - else - return theMaterial->getInitialTangentSensitivity(gradIndex); + else { + if (theMaterial) + return theMaterial->getInitialTangentSensitivity(gradIndex); + else + return 0.0; + } } double PenaltyMaterial::getDampTangentSensitivity(int gradIndex) { - theMaterial->getDampTangentSensitivity(gradIndex); + if (theMaterial) + return theMaterial->getDampTangentSensitivity(gradIndex); + else + return 0.0; } double PenaltyMaterial::getRhoSensitivity(int gradIndex) { - return theMaterial->getRhoSensitivity(gradIndex); + if (theMaterial) + return theMaterial->getRhoSensitivity(gradIndex); + else + return 0.0; } int PenaltyMaterial::commitSensitivity(double strainGradient, int gradIndex, int numGrads) { - return theMaterial->commitSensitivity(strainGradient, gradIndex, numGrads); + if (theMaterial) + return theMaterial->commitSensitivity(strainGradient, gradIndex, numGrads); + else + return -1; } diff --git a/SRC/material/uniaxial/SimpleFractureMaterial.cpp b/SRC/material/uniaxial/SimpleFractureMaterial.cpp index 659442bc8..828fc6b19 100644 --- a/SRC/material/uniaxial/SimpleFractureMaterial.cpp +++ b/SRC/material/uniaxial/SimpleFractureMaterial.cpp @@ -83,17 +83,19 @@ SimpleFractureMaterial::SimpleFractureMaterial(int tag, UniaxialMaterial &materi { theMaterial = material.getCopy(); - Cstress = theMaterial->getStress(); - Ctangent = theMaterial->getTangent(); - Cstrain = theMaterial->getStrain(); - Tstress = Cstress; - Ttangent = Ctangent; - Tstrain = Cstrain; - if (theMaterial == 0) { opserr << "SimpleFractureMaterial::SimpleFractureMaterial -- failed to get copy of material\n"; - exit(-1); + Cstress = 0.0; + Ctangent = 0.0; + Cstrain = 0.0; + } else { + Cstress = theMaterial->getStress(); + Ctangent = theMaterial->getTangent(); + Cstrain = theMaterial->getStrain(); } + Tstress = Cstress; + Ttangent = Ctangent; + Tstrain = Cstrain; } SimpleFractureMaterial::SimpleFractureMaterial() @@ -125,6 +127,9 @@ SimpleFractureMaterial::setTrialStrain(double strain, double strainRate) int SimpleFractureMaterial::setTrialStrain(double strain, double temp, double strainRate) { + if (theMaterial == 0) + return -1; + // determine based on last converged step Tfailed = Cfailed; TstartCompStrain = CstartCompStrain; @@ -208,7 +213,10 @@ SimpleFractureMaterial::getTangent(void) double SimpleFractureMaterial::getDampTangent(void) { - return theMaterial->getDampTangent(); + if (theMaterial) + return theMaterial->getDampTangent(); + else + return 0.0; } @@ -221,12 +229,18 @@ SimpleFractureMaterial::getStrain(void) double SimpleFractureMaterial::getStrainRate(void) { - return theMaterial->getStrainRate(); + if (theMaterial) + return theMaterial->getStrainRate(); + else + return 0.0; } int SimpleFractureMaterial::commitState(void) -{ +{ + if (theMaterial == 0) + return -1; + Cfailed = Tfailed; Cstress = Tstress; Ctangent = Ttangent; @@ -239,6 +253,9 @@ SimpleFractureMaterial::commitState(void) int SimpleFractureMaterial::revertToLastCommit(void) { + if (theMaterial == 0) + return -1; + Tfailed = Cfailed; Tstress = Cstress; Ttangent = Ctangent; @@ -251,6 +268,9 @@ SimpleFractureMaterial::revertToLastCommit(void) int SimpleFractureMaterial::revertToStart(void) { + if (theMaterial == 0) + return -1; + Tfailed = false; Cstrain = 0; theMaterial->revertToStart(); @@ -262,8 +282,9 @@ SimpleFractureMaterial::revertToStart(void) UniaxialMaterial * SimpleFractureMaterial::getCopy(void) { - SimpleFractureMaterial *theCopy = - new SimpleFractureMaterial(this->getTag(), *theMaterial, maxStrain); + SimpleFractureMaterial *theCopy = 0; + if (theMaterial) + theCopy = new SimpleFractureMaterial(this->getTag(), *theMaterial, maxStrain); return theCopy; } @@ -271,6 +292,11 @@ SimpleFractureMaterial::getCopy(void) int SimpleFractureMaterial::sendSelf(int cTag, Channel &theChannel) { + if (theMaterial == 0) { + opserr << "SimpleFractureMaterial::sendSelf() - theMaterial is null, nothing to send\n"; + return -1; + } + int dbTag = this->getDbTag(); static ID dataID(3); @@ -290,11 +316,7 @@ SimpleFractureMaterial::sendSelf(int cTag, Channel &theChannel) static Vector dataVec(6); dataVec(0) = maxStrain; - if (Cfailed == true) - dataVec(1) = 1.0; - else - dataVec(1) = 0.0; - + dataVec(1) = Cfailed ? 1.0 : 0.0; dataVec(2) = Cstress; dataVec(3) = Cstrain; dataVec(4) = Ctangent; @@ -347,13 +369,8 @@ SimpleFractureMaterial::recvSelf(int cTag, Channel &theChannel, } maxStrain = dataVec(0); - - if (dataVec(2) == 1.0) - Cfailed = true; - else - Cfailed = false; - - Cstress = dataVec(2); + Cfailed = dataVec(1) == 1.0 ? true : false; + Cstress = dataVec(2); Cstrain = dataVec(3); Ctangent = dataVec(4); CstartCompStrain = dataVec(5); @@ -371,6 +388,9 @@ void SimpleFractureMaterial::Print(OPS_Stream &s, int flag) { s << "SimpleFractureMaterial tag: " << this->getTag() << endln; - s << "\tMaterial: " << theMaterial->getTag() << endln; + if (theMaterial) + s << "\tMaterial: " << theMaterial->getTag() << endln; + else + s << "\tMaterial is NULL" << endln; s << "\tMax strain: " << maxStrain << endln; } diff --git a/SRC/material/uniaxial/SteelBRB.cpp b/SRC/material/uniaxial/SteelBRB.cpp index 65acd9a46..632bc3267 100755 --- a/SRC/material/uniaxial/SteelBRB.cpp +++ b/SRC/material/uniaxial/SteelBRB.cpp @@ -482,7 +482,8 @@ double SteelBRB::Newton_BRB(double CStress, double beta, double CPlastStrain, d if (fabs(F)>Tol){ opserr<< "Fatal error: SteelBRB::Newton_BRB does not converge ===============\n"; - exit(-1); + //exit(-1); + return 0.0; } return x0; diff --git a/SRC/material/uniaxial/Trilinwpd.cpp b/SRC/material/uniaxial/Trilinwpd.cpp index aeec3a4c0..529aac772 100644 --- a/SRC/material/uniaxial/Trilinwpd.cpp +++ b/SRC/material/uniaxial/Trilinwpd.cpp @@ -652,8 +652,8 @@ trilinwpd::recvSelf(int commitTag, Channel &theChannel, Cstrain = data(25); Ttangent = data(26); pt = data(27); - pb = data(29); - itype=data(30); + pb = data(28); + itype=data(29); // set the trial values TrotMax = CrotMax; TrotMin = CrotMin; diff --git a/Win64/proj/material/material.vcxproj b/Win64/proj/material/material.vcxproj index dcaf0a307..6add27679 100644 --- a/Win64/proj/material/material.vcxproj +++ b/Win64/proj/material/material.vcxproj @@ -239,6 +239,7 @@ + diff --git a/Win64/proj/material/material.vcxproj.filters b/Win64/proj/material/material.vcxproj.filters index 17327347c..0f2d6e0ad 100644 --- a/Win64/proj/material/material.vcxproj.filters +++ b/Win64/proj/material/material.vcxproj.filters @@ -1382,6 +1382,9 @@ uniaxial\py_tz_qz + + section\integration +