C++ NDMaterial类代码示例

int
TimoshenkoSection3d::setTrialSectionDeformation (const Vector &deforms)
{
  int res = 0;
  e = deforms;

  kData[0] = 0.0; kData[1] = 0.0; kData[2] = 0.0; kData[3] = 0.0;
  kData[4] = 0.0; kData[5] = 0.0; kData[6] = 0.0; kData[7] = 0.0;
  kData[8] = 0.0; 
  sData[0] = 0.0; sData[1] = 0.0;  sData[2] = 0.0; 

  int loc = 0;

  double d0 = deforms(0);
  double d1 = deforms(1);
  double d2 = deforms(2);

  for (int i = 0; i < numFibers; i++) {
    NDMaterial *theMat = theMaterials[i];
    double y = matData[loc++] - yBar;
    double z = matData[loc++] - zBar;
    double A = matData[loc++];

    // determine material strain and set it
    double strain = d0 + y*d1 + z*d2;

    Vector eps(3);
    eps(0) = strain;
    res = theMat->setTrialStrain(eps);

    const Vector &stress = theMat->getStress();
    const Matrix &tangent = theMat->getTangent();

    double value = tangent(0,0) * A;
    double vas1 = y*value;
    double vas2 = z*value;
    double vas1as2 = vas1*z;

    kData[0] += value;
    kData[1] += vas1;
    kData[2] += vas2;
    
    kData[4] += vas1 * y;
    kData[5] += vas1as2;
    
    kData[8] += vas2 * z; 

    double fs0 = stress(0) * A;

    sData[0] += fs0;
    sData[1] += fs0 * y;
    sData[2] += fs0 * z;
  }

  kData[3] = kData[1];
  kData[6] = kData[2];
  kData[7] = kData[5];

  return res;
}

int
TimoshenkoSection2d::sendSelf(int commitTag, Channel &theChannel)
{
  int res = 0;

  // create an id to send objects tag and numFibers, 
  //     size 3 so no conflict with matData below if just 1 fiber
  static ID data(3);
  data(0) = this->getTag();
  data(1) = numFibers;
  int dbTag = this->getDbTag();
  res += theChannel.sendID(dbTag, commitTag, data);
  if (res < 0) {
    opserr <<  "TimoshenkoSection2d::sendSelf - failed to send ID data\n";
    return res;
  }    

  if (numFibers != 0) {
    
    // create an id containingg classTag and dbTag for each material & send it
    ID materialData(2*numFibers);
    for (int i=0; i<numFibers; i++) {
      NDMaterial *theMat = theMaterials[i];
      materialData(2*i) = theMat->getClassTag();
      int matDbTag = theMat->getDbTag();
      if (matDbTag == 0) {
	matDbTag = theChannel.getDbTag();
	if (matDbTag != 0)
	  theMat->setDbTag(matDbTag);
      }
      materialData(2*i+1) = matDbTag;
    }    
    
    res += theChannel.sendID(dbTag, commitTag, materialData);
    if (res < 0) {
      opserr <<  "TimoshenkoSection2d::sendSelf - failed to send material data\n";
      return res;
    }    

    // send the fiber data, i.e. area and loc
    Vector fiberData(matData, 3*numFibers);
    res += theChannel.sendVector(dbTag, commitTag, fiberData);
    if (res < 0) {
      opserr <<  "TimoshenkoSection2d::sendSelf - failed to send material data\n";
      return res;
    }    

    // now invoke send(0 on all the materials
    for (int j=0; j<numFibers; j++)
      theMaterials[j]->sendSelf(commitTag, theChannel);

  }

  return res;
}

// constructors:
WSection2d::WSection2d(int tag, NDMaterial &theMat,
		       double D, double Tw, double Bf, double Tf,
		       int Nfdw, int Nftf, double shape, double flag):
  SectionForceDeformation(tag, SEC_TAG_WSection2d),
  theFibers(0), yFibers(0), AFibers(0), e(6),
  d(D), tw(Tw), bf(Bf), tf(Tf), nfdw(Nfdw), nftf(Nftf), shapeFactor(shape)
{
  int numFibers = nfdw + 2*nftf;
  
  theFibers = new NDMaterial*[numFibers];
  yFibers   = new double[numFibers];
  AFibers   = new double[numFibers];
  
  for (int i = 0; i < numFibers; i++) {
    theFibers[i] = flag ? theMat.getCopy("BeamFiber") :
      theMat.getCopy("TimoshenkoFiber");
    if (theFibers[i] == 0)
      opserr << "WSection2d::WSection2d -- failed to get copy of beam fiber" << endln;
  }
  
  double dw = d-2*tf;
  
  double a_f = bf*tf/nftf;
  double a_w = dw*tw/nfdw;
  
  int loc = 0;
  
  double yIncr  = tf/nftf;
  double yStart = 0.5*d - 0.5*yIncr;
  
  for (loc = 0; loc < nftf; loc++) {
    AFibers[loc] = AFibers[numFibers-loc-1] = a_f;
    yFibers[loc] = yStart - yIncr*loc;
    yFibers[numFibers-loc-1] = -yFibers[loc];
  }
  
  yIncr  = dw/nfdw;
  yStart = 0.5*dw - 0.5*yIncr;
  
  int count = 0;
  
  for ( ; loc < numFibers-nftf; loc++, count++) {
    AFibers[loc] = a_w;
    yFibers[loc] = yStart - yIncr*count;
  }
  
  code(0) = SECTION_RESPONSE_P;
  code(1) = SECTION_RESPONSE_MZ;
  code(2) = SECTION_RESPONSE_MY;
  code(3) = SECTION_RESPONSE_VY;
  code(4) = SECTION_RESPONSE_VZ;
  code(5) = SECTION_RESPONSE_T;
}

int
TimoshenkoSection3d::revertToStart(void)
{
  // revert the fibers to start    
  int err = 0;


  kData[0] = 0.0; kData[1] = 0.0; kData[2] = 0.0; kData[3] = 0.0;
  kData[4] = 0.0; kData[5] = 0.0; kData[6] = 0.0; kData[7] = 0.0;
  kData[8] = 0.0; 
  sData[0] = 0.0; sData[1] = 0.0;  sData[2] = 0.0; 

  int loc = 0;

  for (int i = 0; i < numFibers; i++) {
    NDMaterial *theMat = theMaterials[i];
    double y = matData[loc++] - yBar;
    double z = matData[loc++] - zBar;
    double A = matData[loc++];

    // invoke revertToStart on the material
    err += theMat->revertToStart();

    const Matrix &tangent = theMat->getTangent();
    const Vector &stress = theMat->getStress();

    double value = tangent(0,0) * A;
    double vas1 = y*value;
    double vas2 = z*value;
    double vas1as2 = vas1*z;

    kData[0] += value;
    kData[1] += vas1;
    kData[2] += vas2;
    
    kData[4] += vas1 * y;
    kData[5] += vas1as2;
    
    kData[8] += vas2 * z; 

    double fs0 = stress(0) * A;
    sData[0] += fs0;
    sData[1] += fs0 * y;
    sData[2] += fs0 * z;
  }

  kData[3] = kData[1];
  kData[6] = kData[2];
  kData[7] = kData[5];

  return err;
}

ZeroLengthND::ZeroLengthND(int tag, int dim, int Nd1, int Nd2, 
	       const Vector& x, const Vector& yprime, 
		   NDMaterial &theNDmat) : 
Element(tag, ELE_TAG_ZeroLengthND),
connectedExternalNodes(2),
dimension(dim), numDOF(0), 
transformation(3,3), A(0), v(0), e(0.0), K(0), P(0),
end1Ptr(0), end2Ptr(0), theNDMaterial(0), the1DMaterial(0), order(0)
{
	// Obtain copy of Nd material model
	theNDMaterial = theNDmat.getCopy();
	
	if (theNDMaterial == 0) {
		opserr << "ZeroLengthND::zeroLengthND-- failed to get copy of NDMaterial\n";
		exit(-1);
	}
	// Get the material order
	order = theNDMaterial->getOrder();

	// Check material order
	if (order < 2 || order > 3) {
		opserr << "ZeroLengthND::  -- NDMaterial not of order 2 or 3\n";
		exit(-1);
	}

	// Set up the transformation matrix of direction cosines
	this->setUp(Nd1, Nd2, x, yprime);
}

//full constructor
PlateFromPlaneStressMaterialThermal::PlateFromPlaneStressMaterialThermal(    
				   int tag, NDMaterial &ndMat, double g ) :
NDMaterial( tag, ND_TAG_PlateFromPlaneStressMaterialThermal ),
strain(5),gmod(g)
{
  theMat = ndMat.getCopy("PlaneStress") ;
}

//full constructor
ConstantPressureVolumeQuad :: ConstantPressureVolumeQuad( 
                            int tag, 
                  	    int node1,
			    int node2,
			    int node3,
			    int node4,
			    NDMaterial &theMaterial ) :
Element( tag, ELE_TAG_ConstantPressureVolumeQuad ),
connectedExternalNodes(4), load(0)
{
  connectedExternalNodes(0) = node1 ;
  connectedExternalNodes(1) = node2 ;
  connectedExternalNodes(2) = node3 ;
  connectedExternalNodes(3) = node4 ;

  int i ;
  for ( i = 0 ;  i < 4; i++ ) {

      materialPointers[i] = theMaterial.getCopy("AxiSymmetric2D") ;

      if (materialPointers[i] == 0) {
	opserr << "ConstantPressureVolumeQuad::constructor - failed to get a material of type: AxiSymmetric2D\n";
	exit(-1);
      } //end if
      
  } //end for i 

}

GenericSectionNd::GenericSectionNd(int tag, NDMaterial &m, const ID &mCode)
:SectionForceDeformation(tag,SEC_TAG_GenericNd),
otherDbTag(0), theModel(0), code(0)
{
    theModel = m.getCopy();

    if (theModel == 0) {
		g3ErrorHandler->fatal("%s -- failed to get copy of material model",
			"GenericSectionNd::GenericSectionNd");
    }

	order = theModel->getOrder();

	code = new ID(mCode);

	if (code == 0) {
		g3ErrorHandler->fatal("%s -- failed to allocate section ID",
			"GenericSectionNd::GenericSectionNd");
	}

    if (order != code->Size()) {
		g3ErrorHandler->warning("%s -- code size does not match order of material model",
			"GenericSectionNd::GenericSectionNd");
    }
}

NineFourNodeQuadUP::NineFourNodeQuadUP(int tag, 
	int nd1, int nd2, int nd3, int nd4,int nd5, int nd6, int nd7, int nd8,int nd9,
	NDMaterial &m, const char *type, double t, double bulk, double r,
		  double p1, double p2, double b1, double b2)
:Element (tag, ELE_TAG_Nine_Four_Node_QuadUP), 
  theMaterial(0), connectedExternalNodes(9), 
  Ki(0), Q(22), thickness(t), kc(bulk), rho(r)
{
    this->shapeFunction(wu, nintu, nenu, 0);
/*	for( int L = 0; L < nintu; L++) {
		for( int j = 0; j < nenu; j++) {
		printf("%5d %5d %15.6e %15.6e %15.6e\n", L+1, j+1,
			shlu[0][j][L],shlu[1][j][L],shlu[2][j][L]);
		}
	}
	exit(-1);
*/
    this->shapeFunction(wp, nintp, nenp, 1);
    this->shapeFunction(wp, nintp, nenu, 2);

	// Body forces
	b[0] = b1;
	b[1] = b2;
	// Permeabilities
    perm[0] = p1;
    perm[1] = p2;

    // Allocate arrays of pointers to NDMaterials
    theMaterial = new NDMaterial *[nintu];
    
    if (theMaterial == 0) {
      opserr << "NineFourNodeQuadUP::NineFourNodeQuadUP - failed allocate material model pointer\n";
      exit(-1);
    }

    for (int i = 0; i < nintu; i++) {
      
      // Get copies of the material model for each integration point
      theMaterial[i] = m.getCopy(type);
      
      // Check allocation
      if (theMaterial[i] == 0) {
	     opserr << "NineFourNodeQuadUP::NineFourNodeQuadUP -- failed to get a copy of material model\n";
	     exit(-1);
      }
    }

    // Set connected external node IDs
    connectedExternalNodes(0) = nd1;
    connectedExternalNodes(1) = nd2;
    connectedExternalNodes(2) = nd3;
    connectedExternalNodes(3) = nd4;
    connectedExternalNodes(4) = nd5;
    connectedExternalNodes(5) = nd6;
    connectedExternalNodes(6) = nd7;
    connectedExternalNodes(7) = nd8;
    connectedExternalNodes(8) = nd9;
}

FourNodeQuadUP::FourNodeQuadUP(int tag, int nd1, int nd2, int nd3, int nd4,
                               NDMaterial &m, const char *type, double t, double bulk,
                               double r, double p1, double p2, double b1, double b2, double p)
    :Element (tag, ELE_TAG_FourNodeQuadUP),
     theMaterial(0), connectedExternalNodes(4),
     nd1Ptr(0), nd2Ptr(0), nd3Ptr(0), nd4Ptr(0), Ki(0),
     Q(12), pressureLoad(12), applyLoad(0), thickness(t), kc(bulk), rho(r), pressure(p),
     end1InitDisp(0),end2InitDisp(0),end3InitDisp(0),end4InitDisp(0)
{
    pts[0][0] = -0.5773502691896258;
    pts[0][1] = -0.5773502691896258;
    pts[1][0] =  0.5773502691896258;
    pts[1][1] = -0.5773502691896258;
    pts[2][0] =  0.5773502691896258;
    pts[2][1] =  0.5773502691896258;
    pts[3][0] = -0.5773502691896258;
    pts[3][1] =  0.5773502691896258;

    wts[0] = 1.0;
    wts[1] = 1.0;
    wts[2] = 1.0;
    wts[3] = 1.0;

    // Body forces
    b[0] = b1;
    b[1] = b2;
    // Permeabilities
    perm[0] = p1;
    perm[1] = p2;

    // Allocate arrays of pointers to NDMaterials
    theMaterial = new NDMaterial *[4];

    if (theMaterial == 0) {
        opserr << "FourNodeQuadUP::FourNodeQuadUP - failed allocate material model pointer\n";
        exit(-1);
    }

    for (int i = 0; i < 4; i++) {

        // Get copies of the material model for each integration point
        theMaterial[i] = m.getCopy(type);

        // Check allocation
        if (theMaterial[i] == 0) {
            opserr << "FourNodeQuadUP::FourNodeQuadUP -- failed to get a copy of material model\n";
            exit(-1);
        }
    }

    // Set connected external node IDs
    connectedExternalNodes(0) = nd1;
    connectedExternalNodes(1) = nd2;
    connectedExternalNodes(2) = nd3;
    connectedExternalNodes(3) = nd4;
}

PlaneStressFiberMaterial::PlaneStressFiberMaterial(int tag, NDMaterial &theMat)
: NDMaterial(tag, ND_TAG_PlaneStressFiberMaterial),
Tstrain22(0.0), Cstrain22(0.0), twoDtgLastCommit(3,3), theMaterial(0), strain(2)
{
  // Get a copy of the material
  theMaterial = theMat.getCopy();
  
  if (theMaterial == 0) {
    opserr << "PlaneStressFiberMaterial::PlaneStressFiberMaterial -- failed to get copy of material\n";
    exit(-1);
  }
}

Tri31::Tri31(int tag, int nd1, int nd2, int nd3,
	     NDMaterial &m, const char *type, double t,
	     double p, double r, double b1, double b2)
:Element (tag, ELE_TAG_Tri31), 
  theMaterial(0), connectedExternalNodes(3), 
 Q(6), pressureLoad(6), thickness(t), pressure(p), rho(r), Ki(0)
{
	pts[0][0] = 0.333333333333333;
	pts[0][1] = 0.333333333333333;

	wts[0] = 0.5;

	if (strcmp(type,"PlaneStrain") != 0 && strcmp(type,"PlaneStress") != 0
	    && strcmp(type,"PlaneStrain2D") != 0 && strcmp(type,"PlaneStress2D") != 0) {
			opserr << "Tri31::Tri31 -- improper material type: " << type << "for Tri31\n";
	        exit(-1);
	}

	// Body forces
	b[0] = b1;
	b[1] = b2;

	numgp = 1;
	numnodes = 3;

    // Allocate arrays of pointers to NDMaterials
    theMaterial = new NDMaterial *[numgp];
    
    if (theMaterial == 0) {
      opserr << "Tri31::Tri31 - failed allocate material model pointer\n";
      exit(-1);
    }

	int i;
    for (i = 0; i < numgp; i++) {
		// Get copies of the material model for each integration point
        theMaterial[i] = m.getCopy(type);
			
        // Check allocation
        if (theMaterial[i] == 0) {
	   	    opserr << "Tri31::Tri31 -- failed to get a copy of material model\n";
	        exit(-1);
        }
    }

    // Set connected external node IDs
    connectedExternalNodes(0) = nd1;
    connectedExternalNodes(1) = nd2;
    connectedExternalNodes(2) = nd3;
    
    for (i=0; i<numnodes; i++) theNodes[i] = 0;
}

// full constructor
SSPquad::SSPquad(int tag, int Nd1, int Nd2, int Nd3, int Nd4, NDMaterial &theMat,
                 const char *type, double thick, double b1, double b2)
    :Element(tag,ELE_TAG_SSPquad),
     theMaterial(0),
     mExternalNodes(SSPQ_NUM_NODE),
     mTangentStiffness(SSPQ_NUM_DOF,SSPQ_NUM_DOF),
     mInternalForces(SSPQ_NUM_DOF),
     Q(SSPQ_NUM_DOF),
     mMass(SSPQ_NUM_DOF,SSPQ_NUM_DOF),
     mNodeCrd(2,4),
     mStrain(3),
     mStress(3),
     Mmem(3,SSPQ_NUM_DOF),
     Kmem(SSPQ_NUM_DOF,SSPQ_NUM_DOF),
     Kstab(SSPQ_NUM_DOF,SSPQ_NUM_DOF),
     mThickness(thick),
     applyLoad(0)
{
    mExternalNodes(0) = Nd1;
    mExternalNodes(1) = Nd2;
    mExternalNodes(2) = Nd3;
    mExternalNodes(3) = Nd4;

    mThickness = thick;

    b[0] = b1;
    b[1] = b2;

    appliedB[0] = 0.0;
    appliedB[1] = 0.0;

    // get copy of the material object
    NDMaterial *theMatCopy = theMat.getCopy(type);
    if (theMatCopy != 0) {
        theMaterial = (NDMaterial *)theMatCopy;
    } else {
        opserr << "SSPquad::SSPquad - failed to get copy of material model\n";;
    }

    // check material
    if (theMaterial == 0) {
        opserr << "SSPquad::SSPquad - failed to allocate material model pointer\n";
        exit(-1);
    }

    // check the type
    if (strcmp(type,"PlaneStrain") != 0 && strcmp(type,"PlaneStress") != 0) {
        opserr << "SSPquad::SSPquad - improper material type: " << type << "for SSPquad\n";
        exit(-1);
    }
}

//full constructor
MembranePlateFiberSection::MembranePlateFiberSection(    
				   int tag, 
                                   double thickness, 
                                   NDMaterial &Afiber ) :
SectionForceDeformation( tag, SEC_TAG_MembranePlateFiberSection ),
strainResultant(8)
{
  this->h  = thickness ;

  int i ;
  for ( i = 0; i < numFibers; i++ )
      theFibers[i] = Afiber.getCopy( "PlateFiber" ) ;

}

BeamFiberMaterial::BeamFiberMaterial(int tag, NDMaterial &theMat)
: NDMaterial(tag, ND_TAG_BeamFiberMaterial),
Tstrain22(0.0), Tstrain33(0.0), Tgamma23(0.0),
Cstrain22(0.0), Cstrain33(0.0), Cgamma23(0.0),
theMaterial(0), strain(3)

{
  // Get a copy of the material
  theMaterial = theMat.getCopy("ThreeDimensional");
  
  if (theMaterial == 0) {
    opserr << "BeamFiberMaterial::BeamFiberMaterial -- failed to get copy of material\n";
    exit(-1);
  }
}