@@ -192,7 +192,59 @@ Lattice2d :: computeStiffnessMatrix(FloatMatrix &answer, MatResponseMode rMode,
192192 answer .plusProductUnsym (bj , dbj , dV );
193193}
194194
195+ void
196+ Lattice2d :: giveInternalForcesVector (FloatArray & answer ,
197+ TimeStep * tStep , int useUpdatedGpRecord )
198+ {
199+ FloatMatrix b ;
200+ FloatArray u , stress (3 ), strain ;
201+
202+ this -> computeVectorOf (VM_Total , tStep , u );
203+ // subtract initial displacements, if defined
204+ if ( initialDisplacements ) {
205+ u .subtract (* initialDisplacements );
206+ }
195207
208+ answer .clear ();
209+
210+ for ( GaussPoint * gp : * this -> giveDefaultIntegrationRulePtr () ) {
211+ this -> computeBmatrixAt (gp , b );
212+
213+ if ( useUpdatedGpRecord == 1 ) {
214+ auto status = gp -> giveMaterialStatus ();
215+ LatticeMaterialStatus * lmatStat = dynamic_cast < LatticeMaterialStatus * > ( status );
216+ FloatArray stressTemp ;
217+ if ( lmatStat ) {
218+ stressTemp = lmatStat -> giveLatticeStress ();
219+ }
220+ stress .at (1 ) = stressTemp .at (1 );
221+ stress .at (2 ) = stressTemp .at (2 );
222+ stress .at (3 ) = stressTemp .at (6 );
223+ } else {
224+ if ( !this -> isActivated (tStep ) ) {
225+ strain .zero ();
226+ }
227+ strain .beProductOf (b , u );
228+ this -> computeStressVector (stress , strain , gp , tStep );
229+ }
230+
231+ if ( stress .giveSize () == 0 ) {
232+ break ;
233+ }
234+
235+ // compute nodal representation of internal forces using f = B^T*Sigma dV
236+ double dV = this -> computeVolumeAround (gp );
237+ answer .plusProduct (b , stress , dV );
238+
239+ }
240+
241+ // if inactive update state, but no contribution to global system
242+ if ( !this -> isActivated (tStep ) ) {
243+ answer .zero ();
244+ return ;
245+ }
246+ }
247+
196248void Lattice2d :: computeGaussPoints ()
197249// Sets up the array of Gauss Points of the receiver.
198250{
0 commit comments