4747double SurfaceLoad :: oneOverRoot3 = 1.0 /sqrt(3.0 );
4848double SurfaceLoad :: GsPts[4 ][2 ];
4949
50- Matrix SurfaceLoad::tangentStiffness (SL_NUM_DOF , SL_NUM_DOF );
51- Vector SurfaceLoad::internalForces (SL_NUM_DOF );
52- Vector SurfaceLoad::theVector (SL_NUM_DOF );
50+ Matrix SurfaceLoad::tangentStiffness12 (12 ,12 );
51+ Matrix SurfaceLoad::tangentStiffness24 (24 ,24 );
52+ Vector SurfaceLoad::internalForces12 (12 );
53+ Vector SurfaceLoad::internalForces24 (24 );
5354
5455#include < elementAPI.h>
5556static int num_SurfaceLoad = 0 ;
@@ -100,14 +101,15 @@ OPS_SurfaceLoad(void)
100101SurfaceLoad::SurfaceLoad (int tag, int Nd1, int Nd2, int Nd3, int Nd4, double pressure)
101102 :Element(tag,ELE_TAG_SurfaceLoad),
102103 myExternalNodes(SL_NUM_NODE ),
104+ numDOF(0 ), theMatrix(0 ), theVector(0 ),
103105 g1(SL_NUM_NDF ),
104106 g2(SL_NUM_NDF ),
105107 myNhat(SL_NUM_NDF ),
106108 myNI(SL_NUM_NODE ),
107- dcrd1(SL_NUM_NDF ),
108- dcrd2(SL_NUM_NDF ),
109- dcrd3(SL_NUM_NDF ),
110- dcrd4(SL_NUM_NDF )
109+ dcrd1(3 ),
110+ dcrd2(3 ),
111+ dcrd3(3 ),
112+ dcrd4(3 )
111113{
112114 myExternalNodes (0 ) = Nd1;
113115 myExternalNodes (1 ) = Nd2;
@@ -126,19 +128,26 @@ SurfaceLoad::SurfaceLoad(int tag, int Nd1, int Nd2, int Nd3, int Nd4, double pre
126128 my_pressure = pressure;
127129
128130 mLoadFactor = 1.0 ;
131+
132+ tangentStiffness12.Zero ();
133+ tangentStiffness24.Zero ();
134+
135+ internalForces12.Zero ();
136+ internalForces24.Zero ();
129137}
130138
131139SurfaceLoad::SurfaceLoad ()
132140 :Element(0 ,ELE_TAG_SurfaceLoad),
133141 myExternalNodes(SL_NUM_NODE ),
142+ numDOF(0 ), theMatrix(0 ), theVector(0 ),
134143 g1(SL_NUM_NDF ),
135144 g2(SL_NUM_NDF ),
136145 myNhat(SL_NUM_NDF ),
137146 myNI(SL_NUM_NODE ),
138- dcrd1(SL_NUM_NDF ),
139- dcrd2(SL_NUM_NDF ),
140- dcrd3(SL_NUM_NDF ),
141- dcrd4(SL_NUM_NDF )
147+ dcrd1(3 ),
148+ dcrd2(3 ),
149+ dcrd3(3 ),
150+ dcrd4(3 )
142151{
143152}
144153
@@ -168,7 +177,7 @@ SurfaceLoad::getNodePtrs(void)
168177int
169178SurfaceLoad::getNumDOF (void )
170179{
171- return SL_NUM_DOF ;
180+ return numDOF ;
172181}
173182
174183void
@@ -188,7 +197,24 @@ SurfaceLoad::setDomain(Domain *theDomain)
188197 dcrd2 = theNodes[1 ]->getCrds ();
189198 dcrd3 = theNodes[2 ]->getCrds ();
190199 dcrd4 = theNodes[3 ]->getCrds ();
200+ if (3 != dcrd1.Size () || 3 != dcrd2.Size () || 3 != dcrd3.Size () || 3 != dcrd4.Size ()) {
201+ opserr << " SurfaceLoad::setDomain() - nodes are not defined in three dimensions" << endln;
202+ return ;
203+ }
191204
205+ int ndf1 = theNodes[0 ]->getNumberDOF ();
206+ int ndf2 = theNodes[1 ]->getNumberDOF ();
207+ int ndf3 = theNodes[2 ]->getNumberDOF ();
208+ int ndf4 = theNodes[3 ]->getNumberDOF ();
209+ if (ndf1 != ndf2 || ndf1 != ndf3 || ndf1 != ndf4) {
210+ opserr << " SurfaceLoad::setDomain() - nodes have differing numbers of DOFs" << endln;
211+ return ;
212+ }
213+
214+ numDOF = SL_NUM_NODE *ndf1;
215+ theMatrix = (numDOF == 12 ) ? &tangentStiffness12 : &tangentStiffness24;
216+ theVector = (numDOF == 12 ) ? &internalForces12 : &internalForces24;
217+
192218 // call the base class method
193219 this ->DomainComponent ::setDomain (theDomain);
194220}
@@ -254,14 +280,15 @@ SurfaceLoad::UpdateBase(double Xi, double Eta)
254280const Matrix &
255281SurfaceLoad::getTangentStiff (void )
256282{
257- tangentStiffness. Zero ();
258- return tangentStiffness ;
283+ // theMatrix-> Zero();
284+ return *theMatrix ;
259285}
260286
261287const Matrix &
262288SurfaceLoad::getInitialStiff (void )
263289{
264- return getTangentStiff ();
290+ // theMatrix->Zero();
291+ return *theMatrix;
265292}
266293
267294void
@@ -296,22 +323,23 @@ SurfaceLoad::addInertiaLoadToUnbalance(const Vector &accel)
296323const Vector &
297324SurfaceLoad::getResistingForce ()
298325{
299- internalForces. Zero ();
326+ theVector-> Zero ();
300327
328+ int nodeDOF = (numDOF == 12 ) ? 3 : 6 ;
301329 // loop over Gauss points
302330 for (int i = 0 ; i < 4 ; i++) {
303331 this ->UpdateBase (GsPts[i][0 ],GsPts[i][1 ]);
304332
305333 // loop over nodes
306- for (int j = 0 ; j < 4 ; j++) {
307- // loop over dof
334+ for (int j = 0 ; j < SL_NUM_NODE ; j++) {
335+ // loop over displacement dof
308336 for (int k = 0 ; k < 3 ; k++) {
309- internalForces[j* 3 +k] = internalForces [j*3 +k] - mLoadFactor *my_pressure*myNhat (k)*myNI (j);
337+ (*theVector) [j*nodeDOF +k] -= mLoadFactor *my_pressure*myNhat (k)*myNI (j);
310338 }
311339 }
312340 }
313341
314- return internalForces ;
342+ return *theVector ;
315343}
316344
317345const Vector &
0 commit comments