but it seems something wrong.
the node 3 and 4 have the X commitDiap 0.0088, while Y commitDiap 0.0184(this value seems right).
i think X must <<Y
Thank you!
Code: Select all
#include <stdlib.h>
#include <OPS_Globals.h>
#include <StandardStream.h>
#include <ArrayOfTaggedObjects.h>
#include <Domain.h>
#include <Node.h>
#include <ElasticMaterial.h>
#include <SP_Constraint.h>
#include <LoadPattern.h>
#include <LinearSeries.h>
#include <NodalLoad.h>
#include <StaticAnalysis.h>
#include <AnalysisModel.h>
#include <Linear.h>
#include <PenaltyConstraintHandler.h>
#include <DOF_Numberer.h>
#include <RCM.h>
#include <LoadControl.h>
#include <BandSPDLinSOE.h>
#include <BandSPDLinLapackSolver.h>
#include <PlainNumberer.h>
#include <NewtonRaphson.h>
#include <CTestNormDispIncr.h>
#include <concrete01.h>
#include <steel01.h>
#include <D:\projects\_OpenSees_source_\SRC\material\section\repres\patch\QuadPatch.h>
#include <FiberSectionRepr.h>
#include <StraightReinfLayer.h>
#include <LinearCrdTransf2d.h>
#include <FiberSection2d.h>
#include <ForceBeamColumn2d.h>
#include <LobattoBeamIntegration.h>
#include "cell.h"
#include <ReinfBar.h>
#include <UniaxialFiber2d.h>
#include <ElasticBeam2d.h>
#include <TransformationConstraintHandler.h>
#include <BandGenLinLapackSolver.h>
#include <BandGenLinSOE.h>
StandardStream sserr;
OPS_Stream *opserrPtr = &sserr;
double ops_Dt = 0;
Domain *ops_TheActiveDomain = 0;
Element *ops_TheActiveElement = 0;
SectionForceDeformation * buildSection(int secTag,FiberSectionRepr *fiberSectionRepr,UniaxialMaterial **theMaterials);
// main routine
int main(int argc, char **argv)
{
Domain *theDomain = new Domain();
Node *node1=new Node(1,3,0,0);
Node *node2=new Node(2,3,360,0);
Node *node3=new Node(3,3,0,144);
Node *node4=new Node(4,3,360,144);
theDomain->addNode(node1);
theDomain->addNode(node2);
theDomain->addNode(node3);
theDomain->addNode(node4);
SP_Constraint *sp1 = new SP_Constraint(1, 1, 0, 0.0);
SP_Constraint *sp2 = new SP_Constraint(2, 1, 1, 0.0);
SP_Constraint *sp3 = new SP_Constraint(3, 1, 2, 0.0);
SP_Constraint *sp4 = new SP_Constraint(4, 2, 0, 0.0);
SP_Constraint *sp5 = new SP_Constraint(5, 2, 1, 0.0);
SP_Constraint *sp6 = new SP_Constraint(6, 2, 2, 0.0);
theDomain->addSP_Constraint(sp1);
theDomain->addSP_Constraint(sp2);
theDomain->addSP_Constraint(sp3);
theDomain->addSP_Constraint(sp4);
theDomain->addSP_Constraint(sp5);
theDomain->addSP_Constraint(sp6);
UniaxialMaterial *theMaterial1 = new Concrete01(1,-6.0,-0.004,-5.0,-0.014);
UniaxialMaterial *theMaterial2 = new Concrete01(2,-5.0,-0.002,0.0,-0.006);
UniaxialMaterial *theMaterial3 = new Steel01(3,60,30000,0.01);
UniaxialMaterial *theMaterials[3]={theMaterial1,theMaterial2,theMaterial3};
Matrix theCord(4,2);
theCord(0,0)=-10.5;
theCord(0,1)=-6;
theCord(2,0)=10.5;
theCord(2,1)=6;
theCord(1,0) = theCord(2,0);
theCord(1,1) = theCord(0,1);
theCord(3,0) = theCord(0,0);
theCord(3,1) = theCord(2,1);
QuadPatch thePatch1(1,10,1,theCord);
theCord(0,0)=-12;
theCord(0,1)=6;
theCord(2,0)=12;
theCord(2,1)=7.5;
theCord(1,0) = theCord(2,0);
theCord(1,1) = theCord(0,1);
theCord(3,0) = theCord(0,0);
theCord(3,1) = theCord(2,1);
QuadPatch thePatch2(2,10,1,theCord);
theCord(0,0)=-12;
theCord(0,1)=-7.5;
theCord(2,0)=12;
theCord(2,1)=-6;
theCord(1,0) = theCord(2,0);
theCord(1,1) = theCord(0,1);
theCord(3,0) = theCord(0,0);
theCord(3,1) = theCord(2,1);
QuadPatch thePatch3(2,10,1,theCord);
theCord(0,0)=-12;
theCord(0,1)=-6;
theCord(2,0)=-10.5;
theCord(2,1)=6;
theCord(1,0) = theCord(2,0);
theCord(1,1) = theCord(0,1);
theCord(3,0) = theCord(0,0);
theCord(3,1) = theCord(2,1);
QuadPatch thePatch4(2,2,1,theCord);
theCord(0,0)=10.5;
theCord(0,1)=-6;
theCord(2,0)=12;
theCord(2,1)=6;
theCord(1,0) = theCord(2,0);
theCord(1,1) = theCord(0,1);
theCord(3,0) = theCord(0,0);
theCord(3,1) = theCord(2,1);
QuadPatch thePatch5(2,2,1,theCord);
Vector startPt(2);
Vector endPt(2);
startPt(0) = 10.5;
startPt(1) = 6;
endPt(0) = 10.5;
endPt(1) = -6;
StraightReinfLayer reinfLayer1(3,3,0.6,startPt, endPt);
startPt(0) = 0;
startPt(1) = 6;
endPt(0) = 0;
endPt(1) = -6;
StraightReinfLayer reinfLayer2(3,2,0.6,startPt, endPt);
startPt(0) = -10.5;
startPt(1) = 6;
endPt(0) = 10.5;
endPt(1) = -6;
StraightReinfLayer reinfLayer3(3,3,0.6,startPt, endPt);
FiberSectionRepr *theFiberSection=new FiberSectionRepr(1,30,30);
theFiberSection->addPatch(thePatch1);
theFiberSection->addPatch(thePatch2);
theFiberSection->addPatch(thePatch3);
theFiberSection->addPatch(thePatch4);
theFiberSection->addPatch(thePatch5);
theFiberSection->addReinfLayer(reinfLayer1);
theFiberSection->addReinfLayer(reinfLayer2);
theFiberSection->addReinfLayer(reinfLayer3);
FiberSection2d *theSection=(FiberSection2d *)buildSection(1,theFiberSection,theMaterials);
Vector jntOffsetI(2), jntOffsetJ(2);
jntOffsetI.Zero();
jntOffsetJ.Zero();
CrdTransf2d *crdTransf2d;
crdTransf2d = new LinearCrdTransf2d(1, jntOffsetI, jntOffsetJ);
FiberSection2d **sections = new FiberSection2d* [5];
for(int i=0;i<5;i++)
sections[i]=theSection;
LobattoBeamIntegration beamIntegr1;
Element *element1 = new ForceBeamColumn2d(1, 1, 3, 5,
(SectionForceDeformation **)sections,beamIntegr1, *crdTransf2d,
0, 10, 1e-8);
theDomain->addElement(element1);
Element *element2 = new ForceBeamColumn2d(2, 2, 4, 5,
(SectionForceDeformation **)sections,beamIntegr1, *crdTransf2d,
0, 10, 1e-8);
theDomain->addElement(element2);
CrdTransf2d *crdTransf2d2;
crdTransf2d2 = new LinearCrdTransf2d(2, jntOffsetI, jntOffsetJ);
Element *element3=new ElasticBeam2d(3,3,4,360,4030,8640,*crdTransf2d2);
TimeSeries *theSeries = new LinearSeries();
LoadPattern *theLoadPattern = new LoadPattern(1);
theLoadPattern->setTimeSeries(theSeries);
theDomain->addLoadPattern(theLoadPattern);
Vector theLoadValues(3);
theLoadValues(0) = 0;
theLoadValues(1) = -180;
theLoadValues(2) = 0;
NodalLoad *theLoad = new NodalLoad(1, 3, theLoadValues);
theDomain->addNodalLoad(theLoad, 1);
NodalLoad *theLoad1 = new NodalLoad(2, 4, theLoadValues);
theDomain->addNodalLoad(theLoad1, 1);
theDomain->initialize();
opserr<<*theDomain<<endln;
AnalysisModel *theModel = new AnalysisModel();
EquiSolnAlgo *theSolnAlgo = new NewtonRaphson();
StaticIntegrator *theIntegrator = new LoadControl(0.1, 1, 0.1, 0.1);
ConstraintHandler *theHandler = new TransformationConstraintHandler();
RCM *theRCM = new RCM();
DOF_Numberer *theNumberer = new DOF_Numberer(*theRCM);
BandGenLinSolver *theSolver = new BandGenLinLapackSolver();
LinearSOE *theSOE = new BandGenLinSOE(*theSolver);
CTestNormDispIncr *theTest = new CTestNormDispIncr(1e-12,10,3);
StaticAnalysis theAnalysis(*theDomain,
*theHandler,
*theNumberer,
*theModel,
*theSolnAlgo,
*theSOE,
*theIntegrator,theTest);
int numSteps = 10;
theAnalysis.analyze(numSteps);
opserr<<*theDomain;
exit(0);
}
SectionForceDeformation * buildSection(int secTag,FiberSectionRepr *fiberSectionRepr,UniaxialMaterial **theMaterials)
{
int i, j, k;
int numFibers;
int numPatches;
Patch **patch;
int numReinfLayers;
ReinfLayer **reinfLayer;
numPatches = fiberSectionRepr->getNumPatches();
patch = fiberSectionRepr->getPatches();
numReinfLayers = fiberSectionRepr->getNumReinfLayers();
reinfLayer = fiberSectionRepr->getReinfLayers();
int numSectionRepresFibers = fiberSectionRepr->getNumFibers();
Fiber **sectionRepresFibers = fiberSectionRepr->getFibers();
numFibers = numSectionRepresFibers;
for (i = 0; i < numPatches; i++)
numFibers += patch[i]->getNumCells();
for (i = 0; i < numReinfLayers; i++)
numFibers += reinfLayer[i]->getNumReinfBars();
static Vector fiberPosition(2);
int matTag;
ID fibersMaterial(numFibers-numSectionRepresFibers);
Matrix fibersPosition(2,numFibers-numSectionRepresFibers);
Vector fibersArea(numFibers-numSectionRepresFibers);
int numCells;
Cell **cell;
k = 0;
for (i = 0; i < numPatches; i++)
{
numCells = patch[i]->getNumCells();
matTag = patch[i]->getMaterialID();
cell = patch[i]->getCells();
for (j = 0; j < numCells; j++)
{
fibersMaterial(k) = matTag;
fibersArea(k) = cell[j]->getArea();
fiberPosition = cell[j]->getCentroidPosition();
fibersPosition(0,k) = fiberPosition(0);
fibersPosition(1,k) = fiberPosition(1);
k++;
}
for (j = 0; j < numCells; j++)
delete cell[j];
delete [] cell;
}
ReinfBar *reinfBar;
int numReinfBars;
for (i = 0; i < numReinfLayers; i++)
{
numReinfBars = reinfLayer[i]->getNumReinfBars();
reinfBar = reinfLayer[i]->getReinfBars();
matTag = reinfLayer[i]->getMaterialID();
for (j = 0; j < numReinfBars; j++)
{
fibersMaterial(k) = matTag;
fibersArea(k) = reinfBar[j].getArea();
fiberPosition = reinfBar[j].getPosition();
fibersPosition(0,k) = fiberPosition(0);
fibersPosition(1,k) = fiberPosition(1);
k++;
}
delete [] reinfBar;
}
UniaxialMaterial *material;
Fiber **fiber = new Fiber *[numFibers];
for (i=0; i<numSectionRepresFibers; i++)
fiber[i] = sectionRepresFibers[i];
k = 0;
for (i = numSectionRepresFibers; i < numFibers; i++)
{
material = theMaterials[fibersMaterial(k)-1];
fiber[i] = new UniaxialFiber2d(k, *material, fibersArea(k), fibersPosition(0,k));
k++;
}
SectionForceDeformation *section = new FiberSection2d(secTag, numFibers, fiber);
for (i = 0; i < numFibers; i++)
delete fiber[i];
return section;
}