MBDyn-1.7.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups
strforce.cc
Go to the documentation of this file.
1 /* $Header: /var/cvs/mbdyn/mbdyn/mbdyn-1.0/mbdyn/struct/strforce.cc,v 1.47 2017/07/03 22:22:09 masarati Exp $ */
2 /*
3  * MBDyn (C) is a multibody analysis code.
4  * http://www.mbdyn.org
5  *
6  * Copyright (C) 1996-2017
7  *
8  * Pierangelo Masarati <masarati@aero.polimi.it>
9  * Paolo Mantegazza <mantegazza@aero.polimi.it>
10  *
11  * Dipartimento di Ingegneria Aerospaziale - Politecnico di Milano
12  * via La Masa, 34 - 20156 Milano, Italy
13  * http://www.aero.polimi.it
14  *
15  * Changing this copyright notice is forbidden.
16  *
17  * This program is free software; you can redistribute it and/or modify
18  * it under the terms of the GNU General Public License as published by
19  * the Free Software Foundation (version 2 of the License).
20  *
21  *
22  * This program is distributed in the hope that it will be useful,
23  * but WITHOUT ANY WARRANTY; without even the implied warranty of
24  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25  * GNU General Public License for more details.
26  *
27  * You should have received a copy of the GNU General Public License
28  * along with this program; if not, write to the Free Software
29  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
30  */
31 
32 /* Forze */
33 
34 #include "mbconfig.h" /* This goes first in every *.c,*.cc file */
35 
36 #include <cfloat>
37 
38 #include "dataman.h"
39 #include "strforce.h"
40 #include "strforce_impl.h"
41 #include "tpldrive_impl.h"
42 
43 /* AbsoluteDispForce - begin */
44 
45 /* Costruttore non banale */
46 
47 AbsoluteDispForce::AbsoluteDispForce(unsigned int uL,
48  const StructDispNode* pN,
49  const TplDriveCaller<Vec3>* pDC,
50  flag fOut)
51 : Elem(uL, fOut),
52 Force(uL, fOut),
53 f(pDC),
54 pNode(pN)
55 #ifdef USE_NETCDF
56 ,
57 Var_F(0)
58 #endif // USE_NETCDF
59 {
60  NO_OP;
61 }
62 
63 
64 AbsoluteDispForce::~AbsoluteDispForce(void)
65 {
66  NO_OP;
67 }
68 
69 void
70 AbsoluteDispForce::WorkSpaceDim(integer* piNumRows, integer* piNumCols) const
71 {
72  *piNumRows = 3;
73  *piNumCols = 1;
74 }
75 
76 void
77 AbsoluteDispForce::InitialWorkSpaceDim(
78  integer* piNumRows,
79  integer* piNumCols) const
80 {
81  *piNumRows = 3;
82  *piNumCols = 1;
83 }
84 
85 /* Contributo al file di restart */
86 std::ostream&
87 AbsoluteDispForce::Restart(std::ostream& out) const
88 {
89  Force::Restart(out) << ", absolute displacement, "
90  << pNode->GetLabel() << ", ",
91  f.pGetDriveCaller()->Restart(out) << ';' << std::endl;
92  return out;
93 }
94 
95 /* Assembla il residuo */
96 SubVectorHandler&
97 AbsoluteDispForce::AssRes(SubVectorHandler& WorkVec,
98  doublereal /* dCoef */ ,
99  const VectorHandler& /* XCurr */ ,
100  const VectorHandler& /* XPrimeCurr */ )
101 {
102  DEBUGCOUT("Entering AbsoluteDispForce::AssRes()" << std::endl);
103 
104  integer iNumRows;
105  integer iNumCols;
106  WorkSpaceDim(&iNumRows, &iNumCols);
107  WorkVec.ResizeReset(iNumRows);
108 
109  /* Indici delle incognite del nodo */
110  integer iFirstMomentumIndex = pNode->iGetFirstMomentumIndex();
111  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
112  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex + iCnt);
113  }
114 
115  WorkVec.Add(1, f.Get());
116 
117  return WorkVec;
118 }
119 
120 /* Inverse Dynamics*/
121 SubVectorHandler&
122 AbsoluteDispForce::AssRes(SubVectorHandler& WorkVec,
123  const VectorHandler& XCurr,
124  const VectorHandler& XPrimeCurr,
125  const VectorHandler& /* XPrimePrimeCurr */ ,
126  InverseDynamics::Order iOrder)
127 
128 {
129  DEBUGCOUT("Entering AbsoluteDispForce::AssRes()" << std::endl);
130 
131  return AssRes(WorkVec, 1., XCurr, XPrimeCurr);
132 }
133 
134 void
135 AbsoluteDispForce::OutputPrepare(OutputHandler& OH)
136 {
137  if (bToBeOutput()) {
138 #ifdef USE_NETCDF
139  if (OH.UseNetCDF(OutputHandler::FORCES)) {
140  ASSERT(OH.IsOpen(OutputHandler::NETCDF));
141 
142  std::ostringstream os;
143  os << "elem.force." << GetLabel();
144  (void)OH.CreateVar(os.str(), "absolute displacement");
145 
146  // joint sub-data
147  os << '.';
148  Var_F = OH.CreateVar<Vec3>(os.str() + "F", "N",
149  "force components (x, y, z)");
150  }
151 #endif // USE_NETCDF
152  }
153 }
154 
155 void
156 AbsoluteDispForce::Output(OutputHandler& OH) const
157 {
158  if (bToBeOutput()) {
159 #ifdef USE_NETCDF
160  if (OH.UseNetCDF(OutputHandler::FORCES)) {
161  Var_F->put_rec(f.Get().pGetVec(), OH.GetCurrentStep());
162  }
163 #endif // USE_NETCDF
164 
165  if (OH.UseText(OutputHandler::FORCES)) {
166  OH.Forces()
167  << GetLabel()
168  << " " << pNode->GetLabel()
169  << " " << f.Get()
170  << " " << pNode->GetXCurr()
171  << std::endl;
172  }
173  }
174 }
175 
176 
177 /* Contributo al residuo durante l'assemblaggio iniziale */
178 SubVectorHandler&
179 AbsoluteDispForce::InitialAssRes(SubVectorHandler& WorkVec,
180  const VectorHandler& /* XCurr */ )
181 {
182  DEBUGCOUT("Entering AbsoluteDispForce::InitialAssRes()" << std::endl);
183 
184  integer iNumRows;
185  integer iNumCols;
186  InitialWorkSpaceDim(&iNumRows, &iNumCols);
187  WorkVec.ResizeReset(iNumRows);
188 
189  /* Indici delle incognite del nodo */
190  integer iFirstPositionIndex = pNode->iGetFirstPositionIndex();
191  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
192  WorkVec.PutRowIndex(iCnt, iFirstPositionIndex + iCnt);
193  }
194 
195  WorkVec.Add(1, f.Get());
196 
197  return WorkVec;
198 }
199 
200 void
201 AbsoluteDispForce::GetConnectedNodes(std::vector<const Node *>& connectedNodes) const
202 {
203  connectedNodes.resize(1);
204  connectedNodes[0] = pNode;
205 }
206 
207 /* AbsoluteDispForce - end */
208 
209 
210 /* AbsoluteInternalDispForce - begin */
211 
212 /* Costruttore non banale */
213 
214 AbsoluteInternalDispForce::AbsoluteInternalDispForce(unsigned int uL,
215  const StructDispNode* pN1,
216  const StructDispNode* pN2,
217  const TplDriveCaller<Vec3>* pDC,
218  flag fOut)
219 : Elem(uL, fOut),
220 Force(uL, fOut),
221 f(pDC),
222 pNode1(pN1),
223 pNode2(pN2)
224 #ifdef USE_NETCDF
225 ,
226 Var_F(0)
227 #endif // USE_NETCDF
228 {
229  NO_OP;
230 }
231 
232 
233 AbsoluteInternalDispForce::~AbsoluteInternalDispForce(void)
234 
235 {
236  NO_OP;
237 }
238 
239 void
240 AbsoluteInternalDispForce::WorkSpaceDim(integer* piNumRows, integer* piNumCols) const
241 {
242  *piNumRows = 6;
243  *piNumCols = 1;
244 }
245 
246 void
247 AbsoluteInternalDispForce::InitialWorkSpaceDim(
248  integer* piNumRows,
249  integer* piNumCols) const
250 {
251  *piNumRows = 6;
252  *piNumCols = 1;
253 }
254 
255 /* Contributo al file di restart */
256 std::ostream&
257 AbsoluteInternalDispForce::Restart(std::ostream& out) const
258 {
259  Force::Restart(out) << ", absolute internal displacement, "
260  << pNode1->GetLabel() << ", "
261  << pNode2->GetLabel() << ", ",
262  f.pGetDriveCaller()->Restart(out) << ';' << std::endl;
263  return out;
264 }
265 
266 /* Assembla il residuo */
267 SubVectorHandler&
268 AbsoluteInternalDispForce::AssRes(SubVectorHandler& WorkVec,
269  doublereal /* dCoef */ ,
270  const VectorHandler& /* XCurr */ ,
271  const VectorHandler& /* XPrimeCurr */ )
272 {
273  DEBUGCOUT("Entering AbsoluteDispForce::AssRes()" << std::endl);
274 
275  integer iNumRows;
276  integer iNumCols;
277  WorkSpaceDim(&iNumRows, &iNumCols);
278  WorkVec.ResizeReset(iNumRows);
279 
280  /* Indici delle incognite del nodo */
281  integer iFirstMomentumIndex1 = pNode1->iGetFirstMomentumIndex();
282  integer iFirstMomentumIndex2 = pNode2->iGetFirstMomentumIndex();
283  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
284  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex1 + iCnt);
285  WorkVec.PutRowIndex(3 + iCnt, iFirstMomentumIndex2 + iCnt);
286  }
287 
288  Vec3 fTmp(f.Get());
289  WorkVec.Add(1, fTmp);
290  WorkVec.Sub(3 + 1, fTmp);
291 
292  return WorkVec;
293 }
294 
295 /* Inverse Dynamics*/
296 SubVectorHandler&
297 AbsoluteInternalDispForce::AssRes(SubVectorHandler& WorkVec,
298  const VectorHandler& XCurr,
299  const VectorHandler& XPrimeCurr,
300  const VectorHandler& /* XPrimePrimeCurr */ ,
301  InverseDynamics::Order iOrder)
302 {
303  DEBUGCOUT("Entering AbsoluteDispForce::AssRes()" << std::endl);
304 
305  return AssRes(WorkVec, 1., XCurr, XPrimeCurr);
306 }
307 
308 void
309 AbsoluteInternalDispForce::OutputPrepare(OutputHandler& OH)
310 {
311  if (bToBeOutput()) {
312 #ifdef USE_NETCDF
313  if (OH.UseNetCDF(OutputHandler::FORCES)) {
314  ASSERT(OH.IsOpen(OutputHandler::NETCDF));
315 
316  std::ostringstream os;
317  os << "elem.force." << GetLabel();
318  (void)OH.CreateVar(os.str(), "internal absolute displacement");
319 
320  // joint sub-data
321  os << '.';
322  Var_F = OH.CreateVar<Vec3>(os.str() + "F", "N",
323  "force components (x, y, z)");
324  }
325 #endif // USE_NETCDF
326  }
327 }
328 
329 void
330 AbsoluteInternalDispForce::Output(OutputHandler& OH) const
331 {
332  if (bToBeOutput()) {
333 #ifdef USE_NETCDF
334  if (OH.UseNetCDF(OutputHandler::FORCES)) {
335  Var_F->put_rec(f.Get().pGetVec(), OH.GetCurrentStep());
336  }
337 #endif // USE_NETCDF
338 
339  if (OH.UseText(OutputHandler::FORCES)) {
340  Vec3 fTmp(f.Get());
341  OH.Forces()
342  << GetLabel()
343  << " " << pNode1->GetLabel()
344  << " " << fTmp
345  << " " << pNode1->GetXCurr()
346  << " " << pNode2->GetLabel()
347  << " " << -fTmp
348  << " " << pNode2->GetXCurr()
349  << std::endl;
350  }
351 
352  /* TODO: NetCDF */
353  }
354 }
355 
356 
357 /* Contributo al residuo durante l'assemblaggio iniziale */
358 SubVectorHandler&
359 AbsoluteInternalDispForce::InitialAssRes(SubVectorHandler& WorkVec,
360  const VectorHandler& /* XCurr */ )
361 {
362  DEBUGCOUT("Entering AbsoluteDispForce::InitialAssRes()" << std::endl);
363 
364  integer iNumRows;
365  integer iNumCols;
366  InitialWorkSpaceDim(&iNumRows, &iNumCols);
367  WorkVec.ResizeReset(iNumRows);
368 
369  /* Indici delle incognite del nodo */
370  integer iFirstPositionIndex1 = pNode1->iGetFirstPositionIndex();
371  integer iFirstPositionIndex2 = pNode2->iGetFirstPositionIndex();
372  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
373  WorkVec.PutRowIndex(iCnt, iFirstPositionIndex1 + iCnt);
374  WorkVec.PutRowIndex(3 + iCnt, iFirstPositionIndex2 + iCnt);
375  }
376 
377  Vec3 fTmp(f.Get());
378  WorkVec.Add(1, fTmp);
379  WorkVec.Sub(3 + 1, fTmp);
380 
381  return WorkVec;
382 }
383 
384 void
385 AbsoluteInternalDispForce::GetConnectedNodes(std::vector<const Node *>& connectedNodes) const
386 {
387  connectedNodes.resize(2);
388  connectedNodes[0] = pNode1;
389  connectedNodes[1] = pNode2;
390 }
391 
392 /* AbsoluteDispForce - end */
393 
394 
395 /* StructuralForce - begin */
396 
397 /* Costruttore */
398 StructuralForce::StructuralForce(unsigned int uL,
399  const StructNode* pN,
400  const TplDriveCaller<Vec3>* pDC,
401  flag fOut)
402 : Elem(uL, fOut),
403 Force(uL, fOut),
404 f(pDC),
405 pNode(pN)
406 #ifdef USE_NETCDF
407 ,
408 Var_F(0)
409 #endif // USE_NETCDF
410 {
411  ASSERT(pNode != NULL);
412  ASSERT(pNode->GetNodeType() == Node::STRUCTURAL);
413  ASSERT(pDC != NULL);
414 }
415 
416 
417 StructuralForce::~StructuralForce(void)
418 {
419  NO_OP;
420 }
421 
422 void
423 StructuralForce::GetConnectedNodes(std::vector<const Node *>& connectedNodes) const {
424  connectedNodes.resize(1);
425  connectedNodes[0] = pNode;
426 }
427 
428 /* StructuralForce - end */
429 
430 
431 /* AbsoluteForce - begin */
432 
433 /* Costruttore non banale */
434 
435 AbsoluteForce::AbsoluteForce(unsigned int uL,
436  const StructNode* pN,
437  const TplDriveCaller<Vec3>* pDC,
438  const Vec3& TmpArm,
439  flag fOut)
440 : Elem(uL, fOut),
441 StructuralForce(uL, pN, pDC, fOut),
442 Arm(TmpArm)
443 #ifdef USE_NETCDF
444 ,
445 Var_A(0)
446 #endif // USE_NETCDF
447 {
448  NO_OP;
449 }
450 
451 
452 AbsoluteForce::~AbsoluteForce(void)
453 {
454  NO_OP;
455 }
456 
457 void
458 AbsoluteForce::WorkSpaceDim(integer* piNumRows, integer* piNumCols) const
459 {
460  *piNumRows = 6;
461  *piNumCols = 3;
462 }
463 
464 void
465 AbsoluteForce::InitialWorkSpaceDim(
466  integer* piNumRows,
467  integer* piNumCols) const
468 {
469  *piNumRows = 12;
470  *piNumCols = 6;
471 }
472 
473 /* Contributo al file di restart */
474 std::ostream&
475 AbsoluteForce::Restart(std::ostream& out) const
476 {
477  Force::Restart(out) << ", absolute, "
478  << pNode->GetLabel() << ", position, reference, node, ",
479  Arm.Write(out, ", ") << ", ";
480  return f.pGetDriveCaller()->Restart(out) << ';' << std::endl;
481 }
482 
483 VariableSubMatrixHandler&
484 AbsoluteForce::AssJac(VariableSubMatrixHandler& WorkMat,
485  doublereal dCoef,
486  const VectorHandler& /* XCurr */ ,
487  const VectorHandler& /* XPrimeCurr */ )
488 {
489  DEBUGCOUT("Entering AbsoluteForce::AssJac()" << std::endl);
490 
491  FullSubMatrixHandler& WM = WorkMat.SetFull();
492 
493  /* Dimensiona e resetta la matrice di lavoro */
494  WM.ResizeReset(3, 3);
495 
496  integer iFirstPositionIndex = pNode->iGetFirstPositionIndex() + 3;
497  integer iFirstMomentumIndex = pNode->iGetFirstMomentumIndex() + 3;
498  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
499  WM.PutRowIndex(iCnt, iFirstMomentumIndex + iCnt);
500  WM.PutColIndex(iCnt, iFirstPositionIndex + iCnt);
501  }
502 
503  /* Dati */
504  Vec3 TmpArm(pNode->GetRRef()*Arm);
505 
506  /* | F/\ | | F |
507  * | | Delta_g = | |
508  * | (d/\F)/\ | | d/\F |
509  */
510 
511  WM.Sub(1, 1, Mat3x3(MatCrossCross, f.Get(), TmpArm*dCoef));
512 
513  return WorkMat;
514 }
515 
516 
517 /* Assembla il residuo */
518 SubVectorHandler&
519 AbsoluteForce::AssRes(SubVectorHandler& WorkVec,
520  doublereal /* dCoef */ ,
521  const VectorHandler& /* XCurr */ ,
522  const VectorHandler& /* XPrimeCurr */ )
523 {
524  DEBUGCOUT("Entering AbsoluteForce::AssRes()" << std::endl);
525 
526  integer iNumRows;
527  integer iNumCols;
528  WorkSpaceDim(&iNumRows, &iNumCols);
529  WorkVec.ResizeReset(iNumRows);
530 
531  /* Indici delle incognite del nodo */
532  integer iFirstMomentumIndex = pNode->iGetFirstMomentumIndex();
533  for (integer iCnt = 1; iCnt <= 6; iCnt++) {
534  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex + iCnt);
535  }
536 
537  const Mat3x3& R(pNode->GetRCurr());
538  Vec3 F(f.Get());
539  Vec3 M((R*Arm).Cross(F));
540 
541  WorkVec.Add(1, F);
542  WorkVec.Add(4, M);
543 
544  return WorkVec;
545 }
546 
547 /* Inverse Dynamics*/
548 SubVectorHandler&
549 AbsoluteForce::AssRes(SubVectorHandler& WorkVec,
550  const VectorHandler& /* XCurr */ ,
551  const VectorHandler& /* XPrimeCurr */ ,
552  const VectorHandler& /* XPrimePrimeCurr */ ,
553  InverseDynamics::Order iOrder)
554 {
555  DEBUGCOUT("Entering AbsoluteForce::AssRes()" << std::endl);
556 
557  integer iNumRows;
558  integer iNumCols;
559  WorkSpaceDim(&iNumRows, &iNumCols);
560  WorkVec.ResizeReset(iNumRows);
561 
562  /* Indici delle incognite del nodo */
563  integer iFirstMomentumIndex = pNode->iGetFirstPositionIndex();
564  for (integer iCnt = 1; iCnt <= 6; iCnt++) {
565  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex + iCnt);
566  }
567 
568  const Mat3x3& R(pNode->GetRCurr());
569  Vec3 F(f.Get());
570  Vec3 M((R*Arm).Cross(F));
571 
572  WorkVec.Add(1, F);
573  WorkVec.Add(4, M);
574 
575  return WorkVec;
576 }
577 
578 void
579 AbsoluteForce::OutputPrepare(OutputHandler& OH)
580 {
581  if (bToBeOutput()) {
582 #ifdef USE_NETCDF
583  if (OH.UseNetCDF(OutputHandler::FORCES)) {
584  ASSERT(OH.IsOpen(OutputHandler::NETCDF));
585 
586  std::ostringstream os;
587  os << "elem.force." << GetLabel();
588  (void)OH.CreateVar(os.str(), "absolute");
589 
590  // joint sub-data
591  os << '.';
592  Var_F = OH.CreateVar<Vec3>(os.str() + "F", "N",
593  "global force components (x, y, z)");
594 
595  Var_A = OH.CreateVar<Vec3>(os.str() + "Arm", "m",
596  "arm in global frame (x, y, z)");
597  }
598 #endif // USE_NETCDF
599  }
600 }
601 
602 void
603 AbsoluteForce::Output(OutputHandler& OH) const
604 {
605  if (bToBeOutput()) {
606 #ifdef USE_NETCDF
607  if (OH.UseNetCDF(OutputHandler::FORCES)) {
608  Var_F->put_rec(f.Get().pGetVec(), OH.GetCurrentStep());
609  Var_A->put_rec((pNode->GetXCurr() + pNode->GetRCurr()*Arm).pGetVec(), OH.GetCurrentStep());
610  }
611 #endif // USE_NETCDF
612 
613  if (OH.UseText(OutputHandler::FORCES)) {
614  OH.Forces()
615  << GetLabel()
616  << " " << pNode->GetLabel()
617  << " " << f.Get()
618  << " " << pNode->GetXCurr() + pNode->GetRCurr()*Arm
619  << std::endl;
620  }
621 
622  /* TODO: NetCDF */
623  }
624 }
625 
626 
627 /* Contributo allo jacobiano durante l'assemblaggio iniziale */
628 VariableSubMatrixHandler&
629 AbsoluteForce::InitialAssJac(VariableSubMatrixHandler& WorkMat,
630  const VectorHandler& /* XCurr */ )
631 {
632  DEBUGCOUT("Entering AbsoluteForce::InitialAssJac()" << std::endl);
633 
634  FullSubMatrixHandler& WM = WorkMat.SetFull();
635 
636  /* Dimensiona e resetta la matrice di lavoro */
637  WM.ResizeReset(6, 6);
638 
639  integer iFirstPositionIndex = pNode->iGetFirstPositionIndex() + 3;
640  integer iFirstVelocityIndex = iFirstPositionIndex + 6;
641  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
642  WM.PutRowIndex(iCnt, iFirstPositionIndex + iCnt);
643  WM.PutRowIndex(3+iCnt, iFirstVelocityIndex + iCnt);
644  WM.PutColIndex(iCnt, iFirstPositionIndex + iCnt);
645  WM.PutColIndex(3+iCnt, iFirstVelocityIndex + iCnt);
646  }
647 
648  /* Dati */
649  Vec3 TmpArm(pNode->GetRRef()*Arm);
650  Vec3 TmpDir = f.Get();
651  const Vec3& Omega(pNode->GetWRef());
652 
653  /* | F/\ | | F |
654  * | | Delta_g = | |
655  * | (d/\F)/\ | | d/\F |
656  */
657 
658  Mat3x3 MTmp(MatCrossCross, TmpDir, TmpArm);
659 
660  WM.Sub(1, 1, MTmp);
661  WM.Sub(4, 1, Mat3x3(MatCrossCross, TmpDir, Omega)*Mat3x3(MatCross, TmpArm));
662  WM.Sub(4, 4, MTmp);
663 
664  return WorkMat;
665 }
666 
667 
668 /* Contributo al residuo durante l'assemblaggio iniziale */
669 SubVectorHandler&
670 AbsoluteForce::InitialAssRes(SubVectorHandler& WorkVec,
671  const VectorHandler& /* XCurr */ )
672 {
673  DEBUGCOUT("Entering AbsoluteForce::InitialAssRes()" << std::endl);
674 
675  integer iNumRows;
676  integer iNumCols;
677  InitialWorkSpaceDim(&iNumRows, &iNumCols);
678  WorkVec.ResizeReset(iNumRows);
679 
680  /* Indici delle incognite del nodo */
681  integer iFirstPositionIndex = pNode->iGetFirstPositionIndex();
682  integer iFirstVelocityIndex = iFirstPositionIndex + 6;
683  for (integer iCnt = 1; iCnt <= 6; iCnt++) {
684  WorkVec.PutRowIndex(iCnt, iFirstPositionIndex + iCnt);
685  WorkVec.PutRowIndex(6+iCnt, iFirstVelocityIndex + iCnt);
686  }
687 
688  /* Dati */
689  const Mat3x3& R(pNode->GetRCurr());
690  Vec3 TmpDir(f.Get());
691  Vec3 TmpArm(R*Arm);
692  const Vec3& Omega(pNode->GetWCurr());
693 
694  WorkVec.Add(1, TmpDir);
695  WorkVec.Add(4, TmpArm.Cross(TmpDir));
696  /* In 7 non c'e' nulla */
697  WorkVec.Add(10, (Omega.Cross(TmpArm)).Cross(TmpDir));
698 
699  return WorkVec;
700 }
701 
702 /* AbsoluteForce - end */
703 
704 
705 /* FollowerForce - begin */
706 
707 /* Costruttore non banale */
708 
709 FollowerForce::FollowerForce(unsigned int uL, const StructNode* pN,
710  const TplDriveCaller<Vec3>* pDC,
711  const Vec3& TmpArm,
712  flag fOut)
713 : Elem(uL, fOut),
714 StructuralForce(uL, pN, pDC, fOut),
715 Arm(TmpArm)
716 #ifdef USE_NETCDF
717 ,
718 Var_A(0)
719 #endif // USE_NETCDF
720 {
721  NO_OP;
722 }
723 
724 
725 FollowerForce::~FollowerForce(void)
726 {
727  NO_OP;
728 }
729 
730 
731 void
732 FollowerForce::WorkSpaceDim(integer* piNumRows, integer* piNumCols) const
733 {
734  *piNumRows = 6;
735  *piNumCols = 3;
736 }
737 
738 
739 void
740 FollowerForce::InitialWorkSpaceDim(
741  integer* piNumRows,
742  integer* piNumCols) const
743 {
744  *piNumRows = 12;
745  *piNumCols = 6;
746 }
747 
748 
749 /* Contributo al file di restart */
750 std::ostream&
751 FollowerForce::Restart(std::ostream& out) const
752 {
753  Force::Restart(out) << ", follower, "
754  << pNode->GetLabel()
755  << ", position, reference, node, ",
756  Arm.Write(out, ", ") << ", ";
757  return f.pGetDriveCaller()->Restart(out) << ';' << std::endl;
758 }
759 
760 
761 VariableSubMatrixHandler&
762 FollowerForce::AssJac(VariableSubMatrixHandler& WorkMat,
763  doublereal dCoef,
764  const VectorHandler& /* XCurr */ ,
765  const VectorHandler& /* XPrimeCurr */ )
766 {
767  DEBUGCOUT("Entering FollowerForce::AssJac()" << std::endl);
768 
769  FullSubMatrixHandler& WM = WorkMat.SetFull();
770 
771  /* Dimensiona e resetta la matrice di lavoro */
772  integer iNumRows = 0;
773  integer iNumCols = 0;
774  WorkSpaceDim(&iNumRows, &iNumCols);
775  WM.ResizeReset(iNumRows, iNumCols);
776 
777  integer iFirstRotationIndex = pNode->iGetFirstPositionIndex() + 3;
778  integer iFirstMomentumIndex = pNode->iGetFirstMomentumIndex();
779  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
780  WM.PutRowIndex(iCnt, iFirstMomentumIndex + iCnt); /* forza */
781  WM.PutRowIndex(3 + iCnt, iFirstMomentumIndex + 3 + iCnt); /* coppia */
782  WM.PutColIndex(iCnt, iFirstRotationIndex + iCnt); /* rotazione */
783  }
784 
785  /* Dati */
786  const Mat3x3& R(pNode->GetRRef());
787  Vec3 TmpDir(R*(f.Get()*dCoef));
788  Vec3 TmpArm(R*Arm);
789 
790  /* | F/\ | | F |
791  * | | Delta_g = | |
792  * | (d/\F)/\ | | d/\F |
793  */
794 
795  WM.Add(1, 1, Mat3x3(MatCross, TmpDir));
796  WM.Add(4, 1, Mat3x3(MatCross, TmpArm.Cross(TmpDir)));
797 
798  return WorkMat;
799 }
800 
801 
802 /* Assembla il residuo */
803 SubVectorHandler&
804 FollowerForce::AssRes(SubVectorHandler& WorkVec,
805  doublereal /* dCoef */ ,
806  const VectorHandler& /* XCurr */ ,
807  const VectorHandler& /* XPrimeCurr */ )
808 {
809  DEBUGCOUT("Entering FollowerForce::AssRes()" << std::endl);
810 
811  integer iNumRows;
812  integer iNumCols;
813  WorkSpaceDim(&iNumRows, &iNumCols);
814  WorkVec.ResizeReset(iNumRows);
815 
816  /* Indici delle incognite del nodo */
817  integer iFirstMomentumIndex = pNode->iGetFirstMomentumIndex();
818  for (integer iCnt = 1; iCnt <= 6; iCnt++) {
819  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex + iCnt);
820  }
821 
822  /* Dati */
823  const Mat3x3& R(pNode->GetRCurr());
824  Vec3 TmpDir = f.Get();
825  Vec3 F(R*TmpDir);
826  Vec3 M(R*Arm.Cross(TmpDir));
827 
828  WorkVec.Add(1, F);
829  WorkVec.Add(4, M);
830 
831  return WorkVec;
832 }
833 
834 
835 /* Inverse Dynamics*/
836 SubVectorHandler&
837 FollowerForce::AssRes(SubVectorHandler& WorkVec,
838  const VectorHandler& /* XCurr */ ,
839  const VectorHandler& /* XPrimeCurr */ ,
840  const VectorHandler& /* XPrimePrimeCurr */ ,
841  InverseDynamics::Order iOrder)
842 {
843  DEBUGCOUT("Entering FollowerForce::AssRes()" << std::endl);
844 
845  integer iNumRows;
846  integer iNumCols;
847  WorkSpaceDim(&iNumRows, &iNumCols);
848  WorkVec.ResizeReset(iNumRows);
849 
850  /* Indici delle incognite del nodo */
851  integer iFirstMomentumIndex = pNode->iGetFirstPositionIndex();
852  for (integer iCnt = 1; iCnt <= 6; iCnt++) {
853  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex+iCnt);
854  }
855 
856  /* Dati */
857  const Mat3x3& R(pNode->GetRCurr());
858  Vec3 TmpDir = f.Get();
859  Vec3 F(R*TmpDir);
860  Vec3 M(R*Arm.Cross(TmpDir));
861 
862  WorkVec.Add(1, F);
863  WorkVec.Add(4, M);
864 
865  return WorkVec;
866 }
867 
868 void
869 FollowerForce::OutputPrepare(OutputHandler& OH)
870 {
871  if (bToBeOutput()) {
872 #ifdef USE_NETCDF
873  if (OH.UseNetCDF(OutputHandler::FORCES)) {
874  ASSERT(OH.IsOpen(OutputHandler::NETCDF));
875 
876  std::ostringstream os;
877  os << "elem.force." << GetLabel();
878  (void)OH.CreateVar(os.str(), "follower");
879 
880  // joint sub-data
881  if (fToBeOutput() & StructuralForce::OUTPUT_REL) {
882  os << '.';
883  Var_F = OH.CreateVar<Vec3>(os.str() + "f", "N",
884  "local force components (x, y, z)");
885 
886  } else {
887  os << '.';
888  Var_F = OH.CreateVar<Vec3>(os.str() + "F", "N",
889  "global force components (x, y, z)");
890  }
891 
892  Var_A = OH.CreateVar<Vec3>(os.str() + "Arm", "m",
893  "arm in global frame (x, y, z)");
894  }
895 #endif // USE_NETCDF
896  }
897 }
898 
899 void
900 FollowerForce::Output(OutputHandler& OH) const
901 {
902  if (bToBeOutput()) {
903 #ifdef USE_NETCDF
904  if (OH.UseNetCDF(OutputHandler::FORCES)) {
905  if (fToBeOutput() & OUTPUT_REL) {
906  Var_F->put_rec((f.Get()).pGetVec(), OH.GetCurrentStep());
907  } else {
908  Var_F->put_rec((pNode->GetRCurr()*f.Get()).pGetVec(), OH.GetCurrentStep());
909  }
910  Var_A->put_rec((pNode->GetXCurr() + pNode->GetRCurr()*Arm).pGetVec(), OH.GetCurrentStep());
911  }
912 #endif // USE_NETCDF
913 
914  if (OH.UseText(OutputHandler::FORCES)) {
915  OH.Forces()
916  << GetLabel()
917  << " " << pNode->GetLabel();
918 
919  if (fToBeOutput() & OUTPUT_REL) {
920  OH.Forces()
921  << " " << f.Get();
922 
923  } else {
924  OH.Forces()
925  << " " << pNode->GetRCurr()*f.Get();
926  }
927 
928  OH.Forces()
929  << " " << pNode->GetXCurr() + pNode->GetRCurr()*Arm
930  << std::endl;
931  }
932  }
933 }
934 
935 
936 /* Contributo allo jacobiano durante l'assemblaggio iniziale */
937 VariableSubMatrixHandler&
938 FollowerForce::InitialAssJac(VariableSubMatrixHandler& WorkMat,
939  const VectorHandler& /* XCurr */ )
940 {
941  DEBUGCOUT("Entering FollowerForce::InitialAssJac()" << std::endl);
942 
943  FullSubMatrixHandler& WM = WorkMat.SetFull();
944 
945  /* Dimensiona e resetta la matrice di lavoro */
946  WM.ResizeReset(12, 6);
947 
948  integer iFirstPositionIndex = pNode->iGetFirstPositionIndex();
949  integer iFirstVelocityIndex = iFirstPositionIndex + 6;
950  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
951  WM.PutRowIndex(iCnt, iFirstPositionIndex + iCnt);
952  WM.PutRowIndex(3 + iCnt, iFirstPositionIndex + 3 + iCnt);
953  WM.PutRowIndex(6 + iCnt, iFirstVelocityIndex + iCnt);
954  WM.PutRowIndex(9 + iCnt, iFirstVelocityIndex + 3 + iCnt);
955  WM.PutColIndex(iCnt, iFirstPositionIndex + 3 + iCnt);
956  WM.PutColIndex(3 + iCnt, iFirstVelocityIndex + 3 + iCnt);
957  }
958 
959  /* Dati */
960  const Mat3x3& R(pNode->GetRRef());
961  Vec3 TmpArm(R*Arm);
962  Vec3 TmpDir = R*f.Get();
963  const Vec3& Omega(pNode->GetWRef());
964 
965  /* | F/\ | | F |
966  * | | Delta_g = | |
967  * | (d/\F)/\ | | d/\F |
968  */
969 
970  WM.Add(1, 1, Mat3x3(MatCross, TmpDir));
971  WM.Add(4, 1, Mat3x3(MatCross, TmpArm.Cross(TmpDir)));
972  WM.Add(7, 1, Mat3x3(MatCrossCross, Omega, TmpDir));
973  WM.Add(7, 4, Mat3x3(MatCross, TmpDir));
974  WM.Add(10, 1, Mat3x3(MatCrossCross, Omega, TmpArm.Cross(TmpDir)));
975  WM.Add(10, 4, Mat3x3(MatCross, TmpArm.Cross(TmpDir)));
976 
977  return WorkMat;
978 }
979 
980 
981 /* Contributo al residuo durante l'assemblaggio iniziale */
982 SubVectorHandler&
983 FollowerForce::InitialAssRes(SubVectorHandler& WorkVec,
984  const VectorHandler& /* XCurr */ )
985 {
986  DEBUGCOUT("Entering FollowerForce::InitialAssRes()" << std::endl);
987 
988  integer iNumRows;
989  integer iNumCols;
990  InitialWorkSpaceDim(&iNumRows, &iNumCols);
991  WorkVec.ResizeReset(iNumRows);
992 
993  /* Indici delle incognite del nodo */
994  integer iFirstPositionIndex = pNode->iGetFirstPositionIndex();
995  integer iFirstVelocityIndex = iFirstPositionIndex + 6;
996  for (integer iCnt = 1; iCnt <= 6; iCnt++) {
997  WorkVec.PutRowIndex(iCnt, iFirstPositionIndex + iCnt);
998  WorkVec.PutRowIndex(6 + iCnt, iFirstVelocityIndex + iCnt);
999  }
1000 
1001  /* Dati */
1002  const Mat3x3& R(pNode->GetRCurr());
1003  Vec3 TmpDir(R*f.Get());
1004  Vec3 TmpArm(R*Arm);
1005  const Vec3& Omega(pNode->GetWCurr());
1006 
1007  WorkVec.Add(1, TmpDir);
1008  WorkVec.Add(4, TmpArm.Cross(TmpDir));
1009  WorkVec.Add(7, Omega.Cross(TmpDir));
1010  WorkVec.Add(10, (Omega.Cross(TmpArm)).Cross(TmpDir)
1011  + TmpArm.Cross(Omega.Cross(TmpDir)));
1012 
1013  return WorkVec;
1014 }
1015 
1016 /* FollowerForce - end */
1017 
1018 
1019 /* AbsoluteCouple - begin */
1020 
1021 /* Costruttore non banale */
1022 
1023 AbsoluteCouple::AbsoluteCouple(unsigned int uL, const StructNode* pN,
1024  const TplDriveCaller<Vec3>* pDC,
1025  flag fOut)
1026 : Elem(uL, fOut),
1027 StructuralForce(uL, pN, pDC, fOut)
1028 {
1029  NO_OP;
1030 }
1031 
1032 
1033 AbsoluteCouple::~AbsoluteCouple(void)
1034 {
1035  NO_OP;
1036 }
1037 
1038 
1039 void
1040 AbsoluteCouple::WorkSpaceDim(integer* piNumRows, integer* piNumCols) const
1041 {
1042  *piNumRows = 3;
1043  *piNumCols = 1;
1044 }
1045 
1046 
1047 void
1048 AbsoluteCouple::InitialWorkSpaceDim(
1049  integer* piNumRows,
1050  integer* piNumCols) const
1051 {
1052  *piNumRows = 3;
1053  *piNumCols = 1;
1054 }
1055 
1056 
1057 /* Contributo al file di restart */
1058 std::ostream&
1059 AbsoluteCouple::Restart(std::ostream& out) const
1060 {
1061  out << " couple: " << GetLabel() << ", absolute, "
1062  << pNode->GetLabel() << ", ";
1063  return f.pGetDriveCaller()->Restart(out) << ';' << std::endl;
1064 }
1065 
1066 
1067 /* Assembla il residuo */
1068 SubVectorHandler&
1069 AbsoluteCouple::AssRes(SubVectorHandler& WorkVec,
1070  doublereal /* dCoef */ ,
1071  const VectorHandler& /* XCurr */ ,
1072  const VectorHandler& /* XPrimeCurr */ )
1073 {
1074  DEBUGCOUT("Entering AbsoluteCouple::AssRes()" << std::endl);
1075 
1076  integer iNumRows;
1077  integer iNumCols;
1078  WorkSpaceDim(&iNumRows, &iNumCols);
1079  WorkVec.ResizeReset(iNumRows);
1080 
1081  /* Indici delle incognite del nodo */
1082  integer iFirstMomentumIndex = pNode->iGetFirstMomentumIndex()+3;
1083  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
1084  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex+iCnt);
1085  }
1086 
1087  WorkVec.Add(1, f.Get());
1088 
1089  return WorkVec;
1090 }
1091 
1092 /* Inverse Dynamics*/
1093 SubVectorHandler&
1094 AbsoluteCouple::AssRes(SubVectorHandler& WorkVec,
1095  const VectorHandler& /* XCurr */ ,
1096  const VectorHandler& /* XPrimeCurr */ ,
1097  const VectorHandler& /* XPrimePrimeCurr */ ,
1098  InverseDynamics::Order iOrder)
1099 {
1100  DEBUGCOUT("Entering AbsoluteCouple::AssRes()" << std::endl);
1101 
1102  integer iNumRows;
1103  integer iNumCols;
1104  WorkSpaceDim(&iNumRows, &iNumCols);
1105  WorkVec.ResizeReset(iNumRows);
1106 
1107  /* Indici delle incognite del nodo */
1108  integer iFirstMomentumIndex = pNode->iGetFirstPositionIndex()+3;
1109  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
1110  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex+iCnt);
1111  }
1112 
1113  WorkVec.Add(1, f.Get());
1114 
1115  return WorkVec;
1116 }
1117 
1118 void
1119 AbsoluteCouple::OutputPrepare(OutputHandler& OH)
1120 {
1121  if (bToBeOutput()) {
1122 #ifdef USE_NETCDF
1123  if (OH.UseNetCDF(OutputHandler::FORCES)) {
1124  ASSERT(OH.IsOpen(OutputHandler::NETCDF));
1125 
1126  std::ostringstream os;
1127  os << "elem.couple." << GetLabel();
1128  (void)OH.CreateVar(os.str(), "absolute");
1129 
1130  // joint sub-data
1131  os << '.';
1132  Var_F = OH.CreateVar<Vec3>(os.str() + "M", "Nm",
1133  "global couple components (x, y, z)");
1134  }
1135 #endif // USE_NETCDF
1136  }
1137 }
1138 
1139 
1140 void
1141 AbsoluteCouple::Output(OutputHandler& OH) const
1142 {
1143  if (bToBeOutput()) {
1144 #ifdef USE_NETCDF
1145  if (OH.UseNetCDF(OutputHandler::FORCES)) {
1146  Var_F->put_rec(f.Get().pGetVec(), OH.GetCurrentStep());
1147  }
1148 #endif // USE_NETCDF
1149 
1150  if (OH.UseText(OutputHandler::FORCES)) {
1151  OH.Forces()
1152  << GetLabel()
1153  << " " << pNode->GetLabel()
1154  << " " << f.Get()
1155  << std::endl;
1156  }
1157 
1158  /* TODO: NetCDF */
1159  }
1160 }
1161 
1162 
1163 /* Contributo al residuo durante l'assemblaggio iniziale */
1164 SubVectorHandler&
1165 AbsoluteCouple::InitialAssRes(SubVectorHandler& WorkVec,
1166  const VectorHandler& /* XCurr */ )
1167 {
1168  DEBUGCOUT("Entering AbsoluteCouple::InitialAssRes()" << std::endl);
1169 
1170  WorkVec.Resize(3);
1171 
1172  /* Indici delle incognite del nodo */
1173  integer iFirstPositionIndex = pNode->iGetFirstPositionIndex() + 3;
1174  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
1175  WorkVec.PutRowIndex(iCnt, iFirstPositionIndex + iCnt);
1176  }
1177 
1178  WorkVec.Add(1, f.Get());
1179 
1180  return WorkVec;
1181 }
1182 
1183 /* AbsoluteCouple - end */
1184 
1185 
1186 /* FollowerCouple - begin */
1187 
1188 FollowerCouple::FollowerCouple(unsigned int uL, const StructNode* pN,
1189  const TplDriveCaller<Vec3>* pDC,
1190  flag fOut)
1191 : Elem(uL, fOut),
1192 StructuralForce(uL, pN, pDC, fOut)
1193 {
1194  NO_OP;
1195 }
1196 
1197 
1198 FollowerCouple::~FollowerCouple(void)
1199 {
1200  NO_OP;
1201 }
1202 
1203 
1204 void
1205 FollowerCouple::WorkSpaceDim(integer* piNumRows, integer* piNumCols) const
1206 {
1207  *piNumRows = 3;
1208  *piNumCols = 3;
1209 }
1210 
1211 
1212 void
1213 FollowerCouple::InitialWorkSpaceDim(
1214  integer* piNumRows,
1215  integer* piNumCols) const
1216 {
1217  *piNumRows = 6;
1218  *piNumCols = 6;
1219 }
1220 
1221 
1222 /* Contributo al file di restart */
1223 std::ostream&
1224 FollowerCouple::Restart(std::ostream& out) const
1225 {
1226  out << " couple: " << GetLabel() << ", follower, "
1227  << pNode->GetLabel() << ", ";
1228  return f.pGetDriveCaller()->Restart(out) << ';' << std::endl;
1229 }
1230 
1231 
1232 VariableSubMatrixHandler&
1233 FollowerCouple::AssJac(VariableSubMatrixHandler& WorkMat,
1234  doublereal dCoef,
1235  const VectorHandler& /* XCurr */ ,
1236  const VectorHandler& /* XPrimeCurr */ )
1237 {
1238  DEBUGCOUT("Entering FollowerCouple::AssJac()" << std::endl);
1239 
1240  FullSubMatrixHandler& WM = WorkMat.SetFull();
1241 
1242  /* Dimensiona e resetta la matrice di lavoro */
1243  integer iNumRows = 0;
1244  integer iNumCols = 0;
1245  WorkSpaceDim(&iNumRows, &iNumCols);
1246  WM.ResizeReset(iNumRows, iNumCols);
1247 
1248  integer iFirstRotationIndex = pNode->iGetFirstPositionIndex() + 3;
1249  integer iFirstMomentumIndex = pNode->iGetFirstMomentumIndex() + 3;
1250  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
1251  WM.PutRowIndex(iCnt, iFirstMomentumIndex + iCnt); /* coppia */
1252  WM.PutColIndex(iCnt, iFirstRotationIndex + iCnt); /* rotazione */
1253  }
1254 
1255  /* Dati */
1256  const Mat3x3& R(pNode->GetRRef());
1257  Mat3x3 MWedge(MatCross, R*(f.Get()*dCoef));
1258 
1259  /* | M /\| Delta_g = | M | */
1260 
1261  WM.Add(1, 1, MWedge);
1262 
1263  return WorkMat;
1264 }
1265 
1266 
1267 /* Assembla il residuo */
1268 SubVectorHandler&
1269 FollowerCouple::AssRes(SubVectorHandler& WorkVec,
1270  doublereal /* dCoef */ ,
1271  const VectorHandler& /* XCurr */ ,
1272  const VectorHandler& /* XPrimeCurr */ )
1273 {
1274  DEBUGCOUT("Entering FollowerCouple::AssRes()" << std::endl);
1275 
1276  integer iNumRows;
1277  integer iNumCols;
1278  WorkSpaceDim(&iNumRows, &iNumCols);
1279  WorkVec.ResizeReset(iNumRows);
1280 
1281  /* Indici delle incognite del nodo */
1282  integer iFirstMomentumIndex = pNode->iGetFirstMomentumIndex() + 3;
1283  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
1284  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex + iCnt);
1285  }
1286 
1287  /* Dati */
1288  const Mat3x3& R(pNode->GetRCurr());
1289  WorkVec.Add(1, R*f.Get());
1290 
1291  return WorkVec;
1292 }
1293 
1294 
1295 /* Inverse Dynamics*/
1296 SubVectorHandler&
1297 FollowerCouple::AssRes(SubVectorHandler& WorkVec,
1298  const VectorHandler& /* XCurr */ ,
1299  const VectorHandler& /* XPrimeCurr */ ,
1300  const VectorHandler& /* XPrimePrimeCurr */ ,
1301  InverseDynamics::Order iOrder)
1302 {
1303  DEBUGCOUT("Entering FollowerCouple::AssRes()" << std::endl);
1304 
1305  integer iNumRows;
1306  integer iNumCols;
1307  WorkSpaceDim(&iNumRows, &iNumCols);
1308  WorkVec.ResizeReset(iNumRows);
1309 
1310  /* Indici delle incognite del nodo */
1311  integer iFirstMomentumIndex = pNode->iGetFirstPositionIndex()+3;
1312  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
1313  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex+iCnt);
1314  }
1315 
1316  /* Dati */
1317  Mat3x3 R(pNode->GetRCurr());
1318  WorkVec.Add(1, R*f.Get());
1319 
1320  return WorkVec;
1321 }
1322 
1323 void
1324 FollowerCouple::OutputPrepare(OutputHandler& OH)
1325 {
1326  if (bToBeOutput()) {
1327 #ifdef USE_NETCDF
1328  if (OH.UseNetCDF(OutputHandler::FORCES)) {
1329  ASSERT(OH.IsOpen(OutputHandler::NETCDF));
1330 
1331  std::ostringstream os;
1332  os << "elem.couple." << GetLabel();
1333  (void)OH.CreateVar(os.str(), "follower");
1334 
1335  // joint sub-data
1336  os << '.';
1337 
1338  if (fToBeOutput() & StructuralForce::OUTPUT_REL) {
1339  Var_F = OH.CreateVar<Vec3>(os.str() + "m", "Nm",
1340  "local couple components (x, y, z)");
1341 
1342  } else {
1343  Var_F = OH.CreateVar<Vec3>(os.str() + "M", "Nm",
1344  "global couple components (x, y, z)");
1345  }
1346  }
1347 #endif // USE_NETCDF
1348  }
1349 }
1350 
1351 
1352 void
1353 FollowerCouple::Output(OutputHandler& OH) const
1354 {
1355  if (bToBeOutput()) {
1356 #ifdef USE_NETCDF
1357  if (OH.UseNetCDF(OutputHandler::FORCES)) {
1358  if (fToBeOutput() & OUTPUT_REL) {
1359  Var_F->put_rec(f.Get().pGetVec(), OH.GetCurrentStep());
1360  } else {
1361  Var_F->put_rec((pNode->GetRCurr()*f.Get()).pGetVec(), OH.GetCurrentStep());
1362  }
1363  }
1364 #endif // USE_NETCDF
1365 
1366  if (OH.UseText(OutputHandler::FORCES)) {
1367  OH.Forces()
1368  << GetLabel()
1369  << " " << pNode->GetLabel();
1370  if (fToBeOutput() & OUTPUT_REL) {
1371  OH.Forces()
1372  << " " << f.Get()
1373  << std::endl;
1374 
1375  } else {
1376  OH.Forces()
1377  << " " << pNode->GetRCurr()*f.Get()
1378  << std::endl;
1379  }
1380  }
1381  }
1382 }
1383 
1384 
1385 /* Contributo allo jacobiano durante l'assemblaggio iniziale */
1386 VariableSubMatrixHandler&
1387 FollowerCouple::InitialAssJac(VariableSubMatrixHandler& WorkMat,
1388  const VectorHandler& /* XCurr */ )
1389 {
1390  DEBUGCOUT("Entering FollowerCouple::InitialAssJac()" << std::endl);
1391 
1392  FullSubMatrixHandler& WM = WorkMat.SetFull();
1393 
1394  /* Dimensiona e resetta la matrice di lavoro */
1395  WM.ResizeReset(6, 6);
1396 
1397  integer iFirstPositionIndex = pNode->iGetFirstPositionIndex() + 3;
1398  integer iFirstVelocityIndex = iFirstPositionIndex + 6;
1399  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
1400  WM.PutRowIndex(iCnt, iFirstPositionIndex + iCnt);
1401  WM.PutRowIndex(3 + iCnt, iFirstVelocityIndex + iCnt);
1402  WM.PutColIndex(iCnt, iFirstPositionIndex + iCnt);
1403  WM.PutColIndex(3 + iCnt, iFirstVelocityIndex + iCnt);
1404  }
1405 
1406  /* Dati */
1407  Vec3 TmpDir(pNode->GetRRef()*f.Get());
1408  const Vec3& Omega(pNode->GetWRef());
1409 
1410  /* | F/\ | | F |
1411  * | | Delta_g = | |
1412  * | (d/\F)/\ | | d/\F |
1413  */
1414 
1415  WM.Add(1, 1, Mat3x3(MatCross, TmpDir));
1416  WM.Add(4, 1, Mat3x3(MatCrossCross, Omega, TmpDir));
1417  WM.Add(4, 4, Mat3x3(MatCross, TmpDir));
1418 
1419  return WorkMat;
1420 }
1421 
1422 
1423 /* Contributo al residuo durante l'assemblaggio iniziale */
1424 SubVectorHandler&
1425 FollowerCouple::InitialAssRes(SubVectorHandler& WorkVec,
1426  const VectorHandler& /* XCurr */ )
1427 {
1428  DEBUGCOUT("Entering FollowerCouple::InitialAssRes()" << std::endl);
1429 
1430  WorkVec.ResizeReset(6);
1431 
1432  /* Indici delle incognite del nodo */
1433  integer iFirstPositionIndex = pNode->iGetFirstPositionIndex() + 3;
1434  integer iFirstVelocityIndex = iFirstPositionIndex + 6;
1435  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
1436  WorkVec.PutRowIndex(iCnt, iFirstPositionIndex + iCnt);
1437  WorkVec.PutRowIndex(6 + iCnt, iFirstVelocityIndex + iCnt);
1438  }
1439 
1440  /* Dati */
1441  const Mat3x3& R(pNode->GetRCurr());
1442  Vec3 TmpDir(R*f.Get());
1443  const Vec3& Omega(pNode->GetWCurr());
1444 
1445  WorkVec.Add(1, TmpDir);
1446  WorkVec.Add(4, Omega.Cross(TmpDir));
1447 
1448  return WorkVec;
1449 }
1450 
1451 /* FollowerCouple - end */
1452 
1453 
1454 /* StructuralInternalForce - begin */
1455 
1456 /* Costruttore */
1457 StructuralInternalForce::StructuralInternalForce(unsigned int uL,
1458  const StructNode* pN1, const StructNode* pN2,
1459  const TplDriveCaller<Vec3>* pDC,
1460  flag fOut)
1461 : Elem(uL, fOut),
1462 Force(uL, fOut),
1463 f(pDC),
1464 pNode1(pN1), pNode2(pN2)
1465 #ifdef USE_NETCDF
1466 ,
1467 Var_F(0)
1468 #endif // USE_NETCDF
1469 {
1470  ASSERT(pNode1 != NULL);
1471  ASSERT(pNode1->GetNodeType() == Node::STRUCTURAL);
1472  ASSERT(pNode2 != NULL);
1473  ASSERT(pNode2->GetNodeType() == Node::STRUCTURAL);
1474  ASSERT(pDC != NULL);
1475 }
1476 
1477 
1478 StructuralInternalForce::~StructuralInternalForce(void)
1479 {
1480  NO_OP;
1481 }
1482 
1483 void
1484 StructuralInternalForce::GetConnectedNodes(
1485  std::vector<const Node *>& connectedNodes) const {
1486  connectedNodes.resize(2);
1487  connectedNodes[0] = pNode1;
1488  connectedNodes[1] = pNode2;
1489 }
1490 
1491 /* StructuralInternalForce - end */
1492 
1493 
1494 /* AbsoluteInternalForce - begin */
1495 
1496 /* Costruttore non banale */
1497 
1498 AbsoluteInternalForce::AbsoluteInternalForce(unsigned int uL,
1499  const StructNode* pN1, const StructNode* pN2,
1500  const TplDriveCaller<Vec3>* pDC,
1501  const Vec3& TmpArm1, const Vec3& TmpArm2,
1502  flag fOut)
1503 : Elem(uL, fOut),
1504 StructuralInternalForce(uL, pN1, pN2, pDC, fOut),
1505 Arm1(TmpArm1), Arm2(TmpArm2)
1506 #ifdef USE_NETCDF
1507 ,
1508 Var_A1(0),
1509 Var_A2(0)
1510 #endif // USE_NETCDF
1511 {
1512  NO_OP;
1513 }
1514 
1515 
1516 AbsoluteInternalForce::~AbsoluteInternalForce(void)
1517 {
1518  NO_OP;
1519 }
1520 
1521 
1522 void
1523 AbsoluteInternalForce::WorkSpaceDim(
1524  integer* piNumRows,
1525  integer* piNumCols) const
1526 {
1527  *piNumRows = 12;
1528  *piNumCols = 6;
1529 }
1530 
1531 
1532 void
1533 AbsoluteInternalForce::InitialWorkSpaceDim(
1534  integer* piNumRows,
1535  integer* piNumCols) const
1536 {
1537  *piNumRows = 24;
1538  *piNumCols = 12;
1539 }
1540 
1541 
1542 /* Contributo al file di restart */
1543 std::ostream&
1544 AbsoluteInternalForce::Restart(std::ostream& out) const
1545 {
1546  Force::Restart(out) << ", absolute internal, "
1547  << pNode1->GetLabel() << ", position, reference, node, ",
1548  Arm1.Write(out, ", ") << ", "
1549  << pNode2->GetLabel() << ", position, reference, node, ",
1550  Arm2.Write(out, ", ") << ", ";
1551  return f.pGetDriveCaller()->Restart(out) << ';' << std::endl;
1552 }
1553 
1554 
1555 VariableSubMatrixHandler&
1556 AbsoluteInternalForce::AssJac(VariableSubMatrixHandler& WorkMat,
1557  doublereal dCoef,
1558  const VectorHandler& /* XCurr */ ,
1559  const VectorHandler& /* XPrimeCurr */ )
1560 {
1561  DEBUGCOUT("Entering AbsoluteInternalForce::AssJac()" << std::endl);
1562 
1563  FullSubMatrixHandler& WM = WorkMat.SetFull();
1564 
1565  /* Dimensiona e resetta la matrice di lavoro */
1566  WM.ResizeReset(6, 6);
1567 
1568  integer iFirstPositionIndex1 = pNode1->iGetFirstPositionIndex() + 3;
1569  integer iFirstMomentumIndex1 = pNode1->iGetFirstMomentumIndex() + 3;
1570 
1571  integer iFirstPositionIndex2 = pNode2->iGetFirstPositionIndex() + 3;
1572  integer iFirstMomentumIndex2 = pNode2->iGetFirstMomentumIndex() + 3;
1573 
1574  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
1575  WM.PutRowIndex(iCnt, iFirstMomentumIndex1 + iCnt);
1576  WM.PutColIndex(iCnt, iFirstPositionIndex1 + iCnt);
1577 
1578  WM.PutRowIndex(3 + iCnt, iFirstMomentumIndex2 + iCnt);
1579  WM.PutColIndex(3 + iCnt, iFirstPositionIndex2 + iCnt);
1580  }
1581 
1582  /* Dati */
1583  Vec3 TmpArm1(pNode1->GetRRef()*Arm1);
1584  Vec3 TmpArm2(pNode2->GetRRef()*Arm2);
1585  Vec3 TmpDir = f.Get()*dCoef;
1586 
1587  /* | F/\ | | F |
1588  * | | Delta_g = | |
1589  * | (d/\F)/\ | | d/\F |
1590  */
1591 
1592  WM.Sub(1, 1, Mat3x3(MatCrossCross, TmpDir, TmpArm1));
1593  WM.Add(4, 4, Mat3x3(MatCrossCross, TmpDir, TmpArm2));
1594 
1595  return WorkMat;
1596 }
1597 
1598 
1599 /* Assembla il residuo */
1600 SubVectorHandler&
1601 AbsoluteInternalForce::AssRes(SubVectorHandler& WorkVec,
1602  doublereal /* dCoef */ ,
1603  const VectorHandler& /* XCurr */ ,
1604  const VectorHandler& /* XPrimeCurr */ )
1605 {
1606  DEBUGCOUT("Entering AbsoluteInternalForce::AssRes()" << std::endl);
1607 
1608  integer iNumRows;
1609  integer iNumCols;
1610  WorkSpaceDim(&iNumRows, &iNumCols);
1611  WorkVec.ResizeReset(iNumRows);
1612 
1613  /* Indici delle incognite del nodo */
1614  integer iFirstMomentumIndex1 = pNode1->iGetFirstMomentumIndex();
1615  integer iFirstMomentumIndex2 = pNode2->iGetFirstMomentumIndex();
1616  for (integer iCnt = 1; iCnt <= 6; iCnt++) {
1617  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex1 + iCnt);
1618  WorkVec.PutRowIndex(6 + iCnt, iFirstMomentumIndex2 + iCnt);
1619  }
1620 
1621  /* Dati */
1622  Vec3 F(f.Get());
1623  Vec3 M1((pNode1->GetRCurr()*Arm1).Cross(F));
1624  Vec3 M2(F.Cross(pNode2->GetRCurr()*Arm1)); /* - x2 /\ F */
1625 
1626  WorkVec.Add(1, F);
1627  WorkVec.Add(4, M1);
1628  WorkVec.Sub(7, F);
1629  WorkVec.Sub(10, M2);
1630 
1631  return WorkVec;
1632 }
1633 
1634 /* Inverse Dynamics*/
1635 SubVectorHandler&
1636 AbsoluteInternalForce::AssRes(SubVectorHandler& WorkVec,
1637  const VectorHandler& /* XCurr */ ,
1638  const VectorHandler& /* XPrimeCurr */ ,
1639  const VectorHandler& /* XPrimePrimeCurr */ ,
1640  InverseDynamics::Order iOrder)
1641 {
1642  DEBUGCOUT("Entering AbsoluteInternalForce::AssRes()" << std::endl);
1643 
1644  integer iNumRows;
1645  integer iNumCols;
1646  WorkSpaceDim(&iNumRows, &iNumCols);
1647  WorkVec.ResizeReset(iNumRows);
1648 
1649  /* Indici delle incognite del nodo */
1650  integer iFirstMomentumIndex1 = pNode1->iGetFirstPositionIndex();
1651  integer iFirstMomentumIndex2 = pNode2->iGetFirstPositionIndex();
1652  for (integer iCnt = 1; iCnt <= 6; iCnt++) {
1653  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex1 + iCnt);
1654  WorkVec.PutRowIndex(6 + iCnt, iFirstMomentumIndex2 + iCnt);
1655  }
1656 
1657  /* Dati */
1658  Vec3 F(f.Get());
1659  Vec3 M1((pNode1->GetRCurr()*Arm1).Cross(F));
1660  Vec3 M2(F.Cross(pNode2->GetRCurr()*Arm1)); /* - x2 /\ F */
1661 
1662  WorkVec.Add(1, F);
1663  WorkVec.Add(4, M1);
1664  WorkVec.Sub(7, F);
1665  WorkVec.Sub(10, M2);
1666 
1667  return WorkVec;
1668 }
1669 
1670 void
1671 AbsoluteInternalForce::OutputPrepare(OutputHandler& OH)
1672 {
1673  if (bToBeOutput()) {
1674 #ifdef USE_NETCDF
1675  if (OH.UseNetCDF(OutputHandler::FORCES)) {
1676  ASSERT(OH.IsOpen(OutputHandler::NETCDF));
1677 
1678  std::ostringstream os;
1679  os << "elem.force." << GetLabel();
1680  (void)OH.CreateVar(os.str(), "internal absolute");
1681 
1682  // joint sub-data
1683  os << '.';
1684  Var_F = OH.CreateVar<Vec3>(os.str() + "F", "N",
1685  "global force components (x, y, z)");
1686 
1687  Var_A1 = OH.CreateVar<Vec3>(os.str() + "Arm1", "m",
1688  "node 1 arm in global frame (x, y, z)");
1689 
1690  Var_A2 = OH.CreateVar<Vec3>(os.str() + "Arm2", "m",
1691  "node 2 arm in global frame (x, y, z)");
1692  }
1693 #endif // USE_NETCDF
1694  }
1695 }
1696 
1697 void
1698 AbsoluteInternalForce::Output(OutputHandler& OH) const
1699 {
1700  if (bToBeOutput()) {
1701 #ifdef USE_NETCDF
1702  if (OH.UseNetCDF(OutputHandler::FORCES)) {
1703  Var_F->put_rec(f.Get().pGetVec(), OH.GetCurrentStep());
1704  Var_A1->put_rec((pNode1->GetXCurr() + pNode1->GetRCurr()*Arm1).pGetVec(), OH.GetCurrentStep());
1705  Var_A2->put_rec((pNode2->GetXCurr() + pNode2->GetRCurr()*Arm2).pGetVec(), OH.GetCurrentStep());
1706  }
1707 #endif // USE_NETCDF
1708 
1709  if (OH.UseText(OutputHandler::FORCES)) {
1710  Vec3 F(f.Get());
1711  OH.Forces()
1712  << GetLabel()
1713  << " " << pNode1->GetLabel()
1714  << " " << F
1715  << " " << pNode1->GetXCurr() + pNode1->GetRCurr()*Arm1
1716  << " " << pNode2->GetLabel()
1717  << " " << -F
1718  << " " << pNode2->GetXCurr() + pNode2->GetRCurr()*Arm2
1719  << std::endl;
1720  }
1721 
1722  /* TODO: NetCDF */
1723  }
1724 }
1725 
1726 
1727 /* Contributo allo jacobiano durante l'assemblaggio iniziale */
1728 VariableSubMatrixHandler&
1729 AbsoluteInternalForce::InitialAssJac(VariableSubMatrixHandler& WorkMat,
1730  const VectorHandler& /* XCurr */ )
1731 {
1732  DEBUGCOUT("Entering AbsoluteInternalForce::InitialAssJac()"
1733  << std::endl);
1734 
1735  FullSubMatrixHandler& WM = WorkMat.SetFull();
1736 
1737  /* Dimensiona e resetta la matrice di lavoro */
1738  WM.ResizeReset(12, 12);
1739 
1740  integer iFirstPositionIndex1 = pNode1->iGetFirstPositionIndex() + 3;
1741  integer iFirstVelocityIndex1 = iFirstPositionIndex1 + 6;
1742 
1743  integer iFirstPositionIndex2 = pNode2->iGetFirstPositionIndex() + 3;
1744  integer iFirstVelocityIndex2 = iFirstPositionIndex2 + 6;
1745 
1746  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
1747  WM.PutRowIndex(iCnt, iFirstPositionIndex1 + iCnt);
1748  WM.PutRowIndex(3 + iCnt, iFirstVelocityIndex1 + iCnt);
1749 
1750  WM.PutColIndex(iCnt, iFirstPositionIndex1 + iCnt);
1751  WM.PutColIndex(3 + iCnt, iFirstVelocityIndex1 + iCnt);
1752 
1753  WM.PutRowIndex(6 + iCnt, iFirstPositionIndex2 + iCnt);
1754  WM.PutRowIndex(9 + iCnt, iFirstVelocityIndex2 + iCnt);
1755 
1756  WM.PutColIndex(6 + iCnt, iFirstPositionIndex2 + iCnt);
1757  WM.PutColIndex(9 + iCnt, iFirstVelocityIndex2 + iCnt);
1758  }
1759 
1760  /* Dati */
1761  Vec3 TmpArm1(pNode1->GetRRef()*Arm1);
1762  Vec3 TmpArm2(pNode2->GetRRef()*Arm2);
1763  Vec3 TmpDir = f.Get();
1764  Vec3 Omega1(pNode1->GetWRef());
1765  Vec3 Omega2(pNode2->GetWRef());
1766 
1767  /* | F/\ | | F |
1768  * | | Delta_g = | |
1769  * | (d/\F)/\ | | d/\F |
1770  */
1771 
1772  Mat3x3 MTmp(MatCrossCross, TmpDir, TmpArm1);
1773  WM.Sub(1, 1, MTmp);
1774  WM.Sub(4, 1, Mat3x3(MatCrossCross, TmpDir, Omega1)*Mat3x3(MatCross, TmpArm1));
1775  WM.Sub(4, 4, MTmp);
1776 
1777  MTmp = Mat3x3(MatCrossCross, TmpDir, TmpArm2);
1778  WM.Add(7, 7, MTmp);
1779  WM.Add(10, 7, Mat3x3(MatCrossCross, TmpDir, Omega2)*Mat3x3(MatCross, TmpArm2));
1780  WM.Add(10, 10, MTmp);
1781 
1782  return WorkMat;
1783 }
1784 
1785 
1786 /* Contributo al residuo durante l'assemblaggio iniziale */
1787 SubVectorHandler&
1788 AbsoluteInternalForce::InitialAssRes(SubVectorHandler& WorkVec,
1789  const VectorHandler& /* XCurr */ )
1790 {
1791  DEBUGCOUT("Entering AbsoluteInternalForce::InitialAssRes()"
1792  << std::endl);
1793 
1794  integer iNumRows;
1795  integer iNumCols;
1796  InitialWorkSpaceDim(&iNumRows, &iNumCols);
1797  WorkVec.ResizeReset(iNumRows);
1798 
1799  /* Indici delle incognite del nodo */
1800  integer iFirstPositionIndex1 = pNode1->iGetFirstPositionIndex();
1801  integer iFirstVelocityIndex1 = iFirstPositionIndex1 + 6;
1802 
1803  integer iFirstPositionIndex2 = pNode2->iGetFirstPositionIndex();
1804  integer iFirstVelocityIndex2 = iFirstPositionIndex2 + 6;
1805 
1806  for (integer iCnt = 1; iCnt <= 6; iCnt++) {
1807  WorkVec.PutRowIndex(iCnt, iFirstPositionIndex1 + iCnt);
1808  WorkVec.PutRowIndex(6 + iCnt, iFirstVelocityIndex1 + iCnt);
1809 
1810  WorkVec.PutRowIndex(12 + iCnt, iFirstPositionIndex2 + iCnt);
1811  WorkVec.PutRowIndex(18 + iCnt, iFirstVelocityIndex2 + iCnt);
1812  }
1813 
1814  /* Dati */
1815  Vec3 TmpDir(f.Get());
1816  Vec3 TmpArm1(pNode1->GetRCurr()*Arm1);
1817  Vec3 TmpArm2(pNode2->GetRCurr()*Arm2);
1818  const Vec3& Omega1(pNode1->GetWCurr());
1819  const Vec3& Omega2(pNode2->GetWCurr());
1820 
1821  WorkVec.Add(1, TmpDir);
1822  WorkVec.Add(4, TmpArm1.Cross(TmpDir));
1823 
1824  /* In 7 non c'e' nulla */
1825  WorkVec.Add(10, (Omega1.Cross(TmpArm1)).Cross(TmpDir));
1826 
1827  WorkVec.Sub(7, TmpDir);
1828  WorkVec.Sub(10, TmpArm2.Cross(TmpDir));
1829 
1830  /* In 7 non c'e' nulla */
1831  WorkVec.Sub(16, (Omega2.Cross(TmpArm2)).Cross(TmpDir));
1832 
1833  return WorkVec;
1834 }
1835 
1836 /* AbsoluteInternalForce - end */
1837 
1838 
1839 /* FollowerInternalForce - begin */
1840 
1841 /* Costruttore non banale */
1842 
1843 FollowerInternalForce::FollowerInternalForce(unsigned int uL,
1844  const StructNode* pN1, const StructNode* pN2,
1845  const TplDriveCaller<Vec3>* pDC,
1846  const Vec3& TmpArm1, const Vec3& TmpArm2,
1847  flag fOut)
1848 : Elem(uL, fOut),
1849 StructuralInternalForce(uL, pN1, pN2, pDC, fOut),
1850 Arm1(TmpArm1), Arm2(TmpArm2)
1851 #ifdef USE_NETCDF
1852 ,
1853 Var_A1(0),
1854 Var_A2(0)
1855 #endif // USE_NETCDF
1856 {
1857  NO_OP;
1858 }
1859 
1860 
1861 FollowerInternalForce::~FollowerInternalForce(void)
1862 {
1863  NO_OP;
1864 }
1865 
1866 
1867 void
1868 FollowerInternalForce::WorkSpaceDim(
1869  integer* piNumRows,
1870  integer* piNumCols) const
1871 {
1872  *piNumRows = 12;
1873  *piNumCols = 6;
1874 }
1875 
1876 
1877 void
1878 FollowerInternalForce::InitialWorkSpaceDim(
1879  integer* piNumRows,
1880  integer* piNumCols) const
1881 {
1882  *piNumRows = 24;
1883  *piNumCols = 12;
1884 }
1885 
1886 
1887 /* Contributo al file di restart */
1888 std::ostream&
1889 FollowerInternalForce::Restart(std::ostream& out) const
1890 {
1891  Force::Restart(out) << ", follower internal, "
1892  << pNode1->GetLabel()
1893  << ", position, reference, node, ",
1894  Arm1.Write(out, ", ") << ", "
1895  << pNode2->GetLabel()
1896  << ", position, reference, node, ",
1897  Arm2.Write(out, ", ") << ", ";
1898  return f.pGetDriveCaller()->Restart(out) << ';' << std::endl;
1899 }
1900 
1901 
1902 VariableSubMatrixHandler&
1903 FollowerInternalForce::AssJac(VariableSubMatrixHandler& WorkMat,
1904  doublereal dCoef,
1905  const VectorHandler& /* XCurr */ ,
1906  const VectorHandler& /* XPrimeCurr */ )
1907 {
1908  DEBUGCOUT("Entering FollowerInternalForce::AssJac()" << std::endl);
1909 
1910  FullSubMatrixHandler& WM = WorkMat.SetFull();
1911 
1912  /* Dimensiona e resetta la matrice di lavoro */
1913  integer iNumRows = 0;
1914  integer iNumCols = 0;
1915  WorkSpaceDim(&iNumRows, &iNumCols);
1916  WM.ResizeReset(iNumRows, iNumCols);
1917 
1918  integer iFirstRotationIndex1 = pNode1->iGetFirstPositionIndex() + 3;
1919  integer iFirstMomentumIndex1 = pNode1->iGetFirstMomentumIndex();
1920  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
1921  WM.PutRowIndex(iCnt, iFirstMomentumIndex1 + iCnt); /* forza */
1922  WM.PutRowIndex(3 + iCnt, iFirstMomentumIndex1 + 3 + iCnt); /* coppia */
1923  WM.PutColIndex(iCnt, iFirstRotationIndex1 + iCnt); /* rotazione */
1924 
1925  WM.PutRowIndex(6 + iCnt, iFirstMomentumIndex1 + iCnt); /* forza */
1926  WM.PutRowIndex(9 + iCnt, iFirstMomentumIndex1 + 3 + iCnt); /* coppia */
1927  WM.PutColIndex(3 + iCnt, iFirstRotationIndex1 + iCnt); /* rotazione */
1928  }
1929 
1930  /* Dati */
1931  Vec3 TmpDir(pNode1->GetRRef()*(f.Get()*dCoef));
1932  Vec3 TmpArm1(pNode1->GetRRef()*Arm1);
1933  Vec3 TmpArm2(pNode2->GetRRef()*Arm2);
1934 
1935  /* | F/\ 0 | | F |
1936  * | | Delta_g_1 = | |
1937  * | (d1/\F)/\ 0 | | d1/\F |
1938  * | | Delta_g_2 = | |
1939  * | -F/\d2/\ d2/\F/\ | | d2/\F |
1940  */
1941 
1942  WM.Add(1, 1, Mat3x3(MatCross, TmpDir));
1943  WM.Add(4, 1, Mat3x3(MatCross, TmpArm1.Cross(TmpDir)));
1944  WM.Sub(7, 1, Mat3x3(MatCross, TmpDir));
1945  WM.Sub(7, 1, Mat3x3(MatCrossCross, TmpArm2, TmpDir));
1946  WM.Add(7, 4, Mat3x3(MatCrossCross, TmpDir, TmpArm2));
1947 
1948  return WorkMat;
1949 }
1950 
1951 
1952 /* Assembla il residuo */
1953 SubVectorHandler&
1954 FollowerInternalForce::AssRes(SubVectorHandler& WorkVec,
1955  doublereal /* dCoef */ ,
1956  const VectorHandler& /* XCurr */ ,
1957  const VectorHandler& /* XPrimeCurr */ )
1958 {
1959  DEBUGCOUT("Entering FollowerInternalForce::AssRes()" << std::endl);
1960 
1961  integer iNumRows;
1962  integer iNumCols;
1963  WorkSpaceDim(&iNumRows, &iNumCols);
1964  WorkVec.ResizeReset(iNumRows);
1965 
1966  /* Indici delle incognite del nodo */
1967  integer iFirstMomentumIndex1 = pNode1->iGetFirstMomentumIndex();
1968  integer iFirstMomentumIndex2 = pNode2->iGetFirstMomentumIndex();
1969  for (integer iCnt = 1; iCnt <= 6; iCnt++) {
1970  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex1 + iCnt);
1971  WorkVec.PutRowIndex(6 + iCnt, iFirstMomentumIndex2 + iCnt);
1972  }
1973 
1974  /* Dati */
1975  Vec3 TmpDir(f.Get());
1976  Vec3 F(pNode1->GetRCurr()*TmpDir);
1977  Vec3 M1(pNode1->GetRCurr()*Arm1.Cross(TmpDir));
1978  Vec3 M2(F.Cross(pNode2->GetRCurr()*Arm2));
1979 
1980  WorkVec.Add(1, F);
1981  WorkVec.Add(4, M1);
1982  WorkVec.Sub(7, F);
1983  WorkVec.Add(10, M2);
1984 
1985  return WorkVec;
1986 }
1987 
1988 /* Inverse Dynamics*/
1989 SubVectorHandler&
1990 FollowerInternalForce::AssRes(SubVectorHandler& WorkVec,
1991  const VectorHandler& /* XCurr */ ,
1992  const VectorHandler& /* XPrimeCurr */ ,
1993  const VectorHandler& /* XPrimePrimeCurr */ ,
1994  InverseDynamics::Order iOrder)
1995 {
1996  DEBUGCOUT("Entering FollowerInternalForce::AssRes()" << std::endl);
1997 
1998  integer iNumRows;
1999  integer iNumCols;
2000  WorkSpaceDim(&iNumRows, &iNumCols);
2001  WorkVec.ResizeReset(iNumRows);
2002 
2003  /* Indici delle incognite del nodo */
2004  integer iFirstMomentumIndex1 = pNode1->iGetFirstPositionIndex();
2005  integer iFirstMomentumIndex2 = pNode2->iGetFirstPositionIndex();
2006  for (integer iCnt = 1; iCnt <= 6; iCnt++) {
2007  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex1 + iCnt);
2008  WorkVec.PutRowIndex(6+iCnt, iFirstMomentumIndex2 + iCnt);
2009  }
2010 
2011  /* Dati */
2012  Vec3 TmpDir(f.Get());
2013  Vec3 F(pNode1->GetRCurr()*TmpDir);
2014  Vec3 M1(pNode1->GetRCurr()*Arm1.Cross(TmpDir));
2015  Vec3 M2(F.Cross(pNode2->GetRCurr()*Arm2));
2016 
2017  WorkVec.Add(1, F);
2018  WorkVec.Add(4, M1);
2019  WorkVec.Sub(7, F);
2020  WorkVec.Add(10, M2);
2021 
2022  return WorkVec;
2023 }
2024 
2025 void
2026 FollowerInternalForce::OutputPrepare(OutputHandler& OH)
2027 {
2028  if (bToBeOutput()) {
2029 #ifdef USE_NETCDF
2030  if (OH.UseNetCDF(OutputHandler::FORCES)) {
2031  ASSERT(OH.IsOpen(OutputHandler::NETCDF));
2032 
2033  std::ostringstream os;
2034  os << "elem.force." << GetLabel();
2035  (void)OH.CreateVar(os.str(), "internal follower");
2036 
2037  // joint sub-data
2038  os << '.';
2039  if (fToBeOutput() & StructuralForce::OUTPUT_REL) {
2040  Var_F = OH.CreateVar<Vec3>(os.str() + "f", "N",
2041  "local force components (x, y, z)");
2042  } else {
2043  Var_F = OH.CreateVar<Vec3>(os.str() + "F", "N",
2044  "global force components (x, y, z)");
2045  }
2046 
2047  Var_A1 = OH.CreateVar<Vec3>(os.str() + "Arm1", "m",
2048  "node 1 arm in global frame (x, y, z)");
2049 
2050  Var_A2 = OH.CreateVar<Vec3>(os.str() + "Arm2", "m",
2051  "node 2 arm in global frame (x, y, z)");
2052  }
2053 #endif // USE_NETCDF
2054  }
2055 }
2056 
2057 void
2058 FollowerInternalForce::Output(OutputHandler& OH) const
2059 {
2060  if (bToBeOutput()) {
2061 #ifdef USE_NETCDF
2062  if (OH.UseNetCDF(OutputHandler::FORCES)) {
2063  if (fToBeOutput() & StructuralForce::OUTPUT_REL) {
2064  Var_F->put_rec(f.Get().pGetVec(), OH.GetCurrentStep());
2065  } else {
2066  Var_F->put_rec((pNode1->GetRCurr()*f.Get()).pGetVec(), OH.GetCurrentStep());
2067  }
2068  Var_A1->put_rec((pNode1->GetXCurr() + pNode1->GetRCurr()*Arm1).pGetVec(), OH.GetCurrentStep());
2069  Var_A2->put_rec((pNode2->GetXCurr() + pNode2->GetRCurr()*Arm2).pGetVec(), OH.GetCurrentStep());
2070  }
2071 #endif // USE_NETCDF
2072 
2073  if (OH.UseText(OutputHandler::FORCES)) {
2074  Vec3 F;
2075  if (fToBeOutput() & StructuralForce::OUTPUT_REL) {
2076  F = f.Get();
2077  } else {
2078  F = pNode1->GetRCurr()*f.Get();
2079  }
2080 
2081  OH.Forces()
2082  << GetLabel()
2083  << " " << pNode1->GetLabel()
2084  << " " << F
2085  << " " << pNode1->GetXCurr() + pNode1->GetRCurr()*Arm1
2086  << " " << pNode2->GetLabel()
2087  << " " << -F
2088  << " " << pNode2->GetXCurr() + pNode2->GetRCurr()*Arm2
2089  << std::endl;
2090  }
2091 
2092  /* TODO: NetCDF */
2093  }
2094 }
2095 
2096 
2097 /* Contributo allo jacobiano durante l'assemblaggio iniziale */
2098 VariableSubMatrixHandler&
2099 FollowerInternalForce::InitialAssJac(VariableSubMatrixHandler& WorkMat,
2100  const VectorHandler& /* XCurr */ )
2101 {
2102  DEBUGCOUT("Entering FollowerInternalForce::InitialAssJac()" << std::endl);
2103 
2104  FullSubMatrixHandler& WM = WorkMat.SetFull();
2105 
2106  integer iNumRows;
2107  integer iNumCols;
2108  WorkSpaceDim(&iNumRows, &iNumCols);
2109 
2110  /* Dimensiona e resetta la matrice di lavoro */
2111  WM.ResizeReset(iNumRows, iNumCols);
2112 
2113  integer iFirstPositionIndex1 = pNode1->iGetFirstPositionIndex();
2114  integer iFirstVelocityIndex1 = iFirstPositionIndex1 + 6;
2115 
2116  integer iFirstPositionIndex2 = pNode2->iGetFirstPositionIndex();
2117  integer iFirstVelocityIndex2 = iFirstPositionIndex2 + 6;
2118 
2119  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
2120  WM.PutRowIndex(iCnt, iFirstPositionIndex1 + iCnt);
2121  WM.PutRowIndex(3 + iCnt, iFirstPositionIndex1 + 3 + iCnt);
2122  WM.PutRowIndex(6 + iCnt, iFirstVelocityIndex1 + iCnt);
2123  WM.PutRowIndex(9 + iCnt, iFirstVelocityIndex1 + 3 + iCnt);
2124 
2125  WM.PutColIndex(iCnt, iFirstPositionIndex1 +3+iCnt);
2126  WM.PutColIndex(3 + iCnt, iFirstVelocityIndex1 + 3 + iCnt);
2127 
2128  WM.PutRowIndex(12 + iCnt, iFirstPositionIndex2 + iCnt);
2129  WM.PutRowIndex(15 + iCnt, iFirstPositionIndex2 + 3 + iCnt);
2130  WM.PutRowIndex(18 + iCnt, iFirstVelocityIndex2 + iCnt);
2131  WM.PutRowIndex(21 + iCnt, iFirstVelocityIndex2 + 3 + iCnt);
2132 
2133  WM.PutColIndex(12 + iCnt, iFirstPositionIndex2 + 3 + iCnt);
2134  WM.PutColIndex(15 + iCnt, iFirstVelocityIndex2 + 3 + iCnt);
2135  }
2136 
2137  /* Dati */
2138  Vec3 TmpArm1(pNode1->GetRRef()*Arm1);
2139  Vec3 TmpArm2(pNode2->GetRRef()*Arm2);
2140  Vec3 TmpDir = pNode1->GetRRef()*f.Get();
2141  const Vec3& Omega1(pNode1->GetWRef());
2142  const Vec3& Omega2(pNode2->GetWRef());
2143 
2144  WM.Add(1, 1, Mat3x3(MatCross, TmpDir));
2145  WM.Add(4, 1, Mat3x3(MatCross, TmpArm1.Cross(TmpDir)));
2146 
2147  WM.Add(7, 1, Mat3x3(MatCrossCross, Omega1, TmpDir));
2148  WM.Add(7, 4, Mat3x3(MatCross, TmpDir));
2149  WM.Add(10, 1, Mat3x3(MatCrossCross, Omega1, TmpArm1.Cross(TmpDir)));
2150  WM.Add(10, 4, Mat3x3(MatCross, TmpArm1.Cross(TmpDir)));
2151 
2152  WM.Sub(13, 1, Mat3x3(MatCross, TmpDir));
2153  WM.Add(16, 1, Mat3x3(MatCrossCross, TmpArm2, TmpDir));
2154  WM.Sub(16, 7, Mat3x3(MatCrossCross, TmpDir, TmpArm2));
2155 
2156  WM.Sub(19, 1, Mat3x3(MatCrossCross, Omega1, TmpDir));
2157  WM.Sub(19, 4, Mat3x3(MatCross, TmpDir));
2158 
2159  WM.Add(22, 1, Mat3x3(MatCrossCross, TmpArm2, Omega1)*Mat3x3(MatCross, TmpDir)
2160  - Mat3x3(MatCrossCross, Omega2.Cross(TmpArm2), TmpDir));
2161  WM.Add(22, 4, Mat3x3(MatCrossCross, TmpArm2, TmpDir));
2162  WM.Add(22, 7, Mat3x3(MatCrossCross, TmpDir, Omega2)*Mat3x3(MatCross, TmpArm2)
2163  - Mat3x3(MatCrossCross, Omega1.Cross(TmpDir), TmpArm2));
2164  WM.Add(22, 10, Mat3x3(MatCrossCross, TmpDir, TmpArm2));
2165 
2166  return WorkMat;
2167 }
2168 
2169 
2170 /* Contributo al residuo durante l'assemblaggio iniziale */
2171 SubVectorHandler&
2172 FollowerInternalForce::InitialAssRes(SubVectorHandler& WorkVec,
2173  const VectorHandler& /* XCurr */ )
2174 {
2175  DEBUGCOUT("Entering FollowerInternalForce::InitialAssRes()" << std::endl);
2176 
2177  integer iNumRows;
2178  integer iNumCols;
2179  InitialWorkSpaceDim(&iNumRows, &iNumCols);
2180  WorkVec.ResizeReset(iNumRows);
2181 
2182  /* Indici delle incognite del nodo */
2183  integer iFirstPositionIndex1 = pNode1->iGetFirstPositionIndex();
2184  integer iFirstVelocityIndex1 = iFirstPositionIndex1 + 6;
2185 
2186  integer iFirstPositionIndex2 = pNode2->iGetFirstPositionIndex();
2187  integer iFirstVelocityIndex2 = iFirstPositionIndex2 + 6;
2188 
2189  for (integer iCnt = 1; iCnt <= 6; iCnt++) {
2190  WorkVec.PutRowIndex(iCnt, iFirstPositionIndex1 + iCnt);
2191  WorkVec.PutRowIndex(6 + iCnt, iFirstVelocityIndex1 + iCnt);
2192 
2193  WorkVec.PutRowIndex(12 + iCnt, iFirstPositionIndex2 + iCnt);
2194  WorkVec.PutRowIndex(18 + iCnt, iFirstVelocityIndex2 + iCnt);
2195  }
2196 
2197  /* Dati */
2198  Vec3 TmpDir(pNode1->GetRCurr()*f.Get());
2199  Vec3 TmpArm1(pNode1->GetRCurr()*Arm1);
2200  Vec3 TmpArm2(pNode2->GetRCurr()*Arm2);
2201  const Vec3& Omega1(pNode1->GetWCurr());
2202  const Vec3& Omega2(pNode2->GetWCurr());
2203 
2204  WorkVec.Add(1, TmpDir);
2205  WorkVec.Add(4, TmpArm1.Cross(TmpDir));
2206  WorkVec.Add(7, Omega1.Cross(TmpDir));
2207  WorkVec.Add(10, (Omega1.Cross(TmpArm1)).Cross(TmpDir)
2208  + TmpArm1.Cross(Omega1.Cross(TmpDir)));
2209 
2210  WorkVec.Sub(13, TmpDir);
2211  WorkVec.Sub(16, TmpArm2.Cross(TmpDir));
2212  WorkVec.Sub(19, Omega1.Cross(TmpDir));
2213  WorkVec.Sub(22, (Omega2.Cross(TmpArm2)).Cross(TmpDir)
2214  + TmpArm2.Cross(Omega1.Cross(TmpDir)));
2215 
2216  return WorkVec;
2217 }
2218 
2219 /* FollowerInternalForce - end */
2220 
2221 
2222 /* AbsoluteInternalCouple - begin */
2223 
2224 /* Costruttore non banale */
2225 
2226 AbsoluteInternalCouple::AbsoluteInternalCouple(unsigned int uL,
2227  const StructNode* pN1, const StructNode* pN2,
2228  const TplDriveCaller<Vec3>* pDC,
2229  flag fOut)
2230 : Elem(uL, fOut),
2231 StructuralInternalForce(uL, pN1, pN2, pDC, fOut)
2232 {
2233  NO_OP;
2234 }
2235 
2236 
2237 AbsoluteInternalCouple::~AbsoluteInternalCouple(void)
2238 {
2239  NO_OP;
2240 }
2241 
2242 void
2243 AbsoluteInternalCouple::WorkSpaceDim(
2244  integer* piNumRows,
2245  integer* piNumCols) const
2246 {
2247  *piNumRows = 6;
2248  *piNumCols = 1;
2249 }
2250 
2251 
2252 void
2253 AbsoluteInternalCouple::InitialWorkSpaceDim(
2254  integer* piNumRows,
2255  integer* piNumCols) const
2256 {
2257  *piNumRows = 6;
2258  *piNumCols = 1;
2259 }
2260 
2261 
2262 /* Contributo al file di restart */
2263 std::ostream&
2264 AbsoluteInternalCouple::Restart(std::ostream& out) const
2265 {
2266  out << " couple: " << GetLabel() << ", absolute internal, "
2267  << pNode1->GetLabel() << ", "
2268  << pNode2->GetLabel() << ", ";
2269  return f.pGetDriveCaller()->Restart(out) << ';' << std::endl;
2270 }
2271 
2272 
2273 /* Assembla il residuo */
2274 SubVectorHandler&
2275 AbsoluteInternalCouple::AssRes(SubVectorHandler& WorkVec,
2276  doublereal /* dCoef */ ,
2277  const VectorHandler& /* XCurr */ ,
2278  const VectorHandler& /* XPrimeCurr */ )
2279 {
2280  DEBUGCOUT("Entering AbsoluteInternalCouple::AssRes()" << std::endl);
2281 
2282  integer iNumRows;
2283  integer iNumCols;
2284  WorkSpaceDim(&iNumRows, &iNumCols);
2285  WorkVec.ResizeReset(iNumRows);
2286 
2287  /* Indici delle incognite del nodo */
2288  integer iFirstMomentumIndex1 = pNode1->iGetFirstMomentumIndex() + 3;
2289  integer iFirstMomentumIndex2 = pNode2->iGetFirstMomentumIndex() + 3;
2290  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
2291  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex1 + iCnt);
2292  WorkVec.PutRowIndex(3 + iCnt, iFirstMomentumIndex2 + iCnt);
2293  }
2294 
2295  /* Dati */
2296  Vec3 F(f.Get());
2297 
2298  WorkVec.Add(1, F);
2299  WorkVec.Sub(4, F);
2300 
2301  return WorkVec;
2302 }
2303 
2304 /* Inverse Dynamics*/
2305 SubVectorHandler&
2306 AbsoluteInternalCouple::AssRes(SubVectorHandler& WorkVec,
2307  const VectorHandler& /* XCurr */ ,
2308  const VectorHandler& /* XPrimeCurr */ ,
2309  const VectorHandler& /* XPrimePrimeCurr */ ,
2310  InverseDynamics::Order iOrder)
2311 {
2312  DEBUGCOUT("Entering AbsoluteInternalCouple::AssRes()" << std::endl);
2313 
2314  integer iNumRows;
2315  integer iNumCols;
2316  WorkSpaceDim(&iNumRows, &iNumCols);
2317  WorkVec.ResizeReset(iNumRows);
2318 
2319  /* Indici delle incognite del nodo */
2320  integer iFirstMomentumIndex1 = pNode1->iGetFirstPositionIndex() + 3;
2321  integer iFirstMomentumIndex2 = pNode2->iGetFirstPositionIndex() + 3;
2322  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
2323  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex1 + iCnt);
2324  WorkVec.PutRowIndex(3 + iCnt, iFirstMomentumIndex2 + iCnt);
2325  }
2326 
2327  /* Dati */
2328  Vec3 F(f.Get());
2329 
2330  WorkVec.Add(1, F);
2331  WorkVec.Sub(4, F);
2332 
2333  return WorkVec;
2334 }
2335 
2336 void
2337 AbsoluteInternalCouple::OutputPrepare(OutputHandler& OH)
2338 {
2339  if (bToBeOutput()) {
2340 #ifdef USE_NETCDF
2341  if (OH.UseNetCDF(OutputHandler::FORCES)) {
2342  ASSERT(OH.IsOpen(OutputHandler::NETCDF));
2343 
2344  std::ostringstream os;
2345  os << "elem.couple." << GetLabel();
2346  (void)OH.CreateVar(os.str(), "internal absolute");
2347 
2348  // joint sub-data
2349  os << '.';
2350  Var_F = OH.CreateVar<Vec3>(os.str() + "M", "Nm",
2351  "global couple components (x, y, z)");
2352  }
2353 #endif // USE_NETCDF
2354  }
2355 }
2356 
2357 void
2358 AbsoluteInternalCouple::Output(OutputHandler& OH) const
2359 {
2360  if (bToBeOutput()) {
2361 #ifdef USE_NETCDF
2362  if (OH.UseNetCDF(OutputHandler::FORCES)) {
2363  Var_F->put_rec(f.Get().pGetVec(), OH.GetCurrentStep());
2364  }
2365 #endif // USE_NETCDF
2366 
2367  if (OH.UseText(OutputHandler::FORCES)) {
2368  Vec3 F(f.Get());
2369  OH.Forces()
2370  << GetLabel()
2371  << " " << pNode1->GetLabel()
2372  << " " << F
2373  << " " << pNode2->GetLabel()
2374  << " " << -F
2375  << std::endl;
2376  }
2377 
2378  /* TODO: NetCDF */
2379  }
2380 }
2381 
2382 
2383 /* Contributo al residuo durante l'assemblaggio iniziale */
2384 SubVectorHandler&
2385 AbsoluteInternalCouple::InitialAssRes(SubVectorHandler& WorkVec,
2386  const VectorHandler& /* XCurr */ )
2387 {
2388  DEBUGCOUT("Entering AbsoluteInternalCouple::InitialAssRes()" << std::endl);
2389 
2390  WorkVec.ResizeReset(6);
2391 
2392  /* Indici delle incognite del nodo */
2393  integer iFirstPositionIndex1 = pNode1->iGetFirstPositionIndex() + 3;
2394  integer iFirstPositionIndex2 = pNode2->iGetFirstPositionIndex() + 3;
2395  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
2396  WorkVec.PutRowIndex(iCnt, iFirstPositionIndex1 + iCnt);
2397  WorkVec.PutRowIndex(3 + iCnt, iFirstPositionIndex2 + iCnt);
2398  }
2399 
2400  /* Dati */
2401  Vec3 F(f.Get());
2402 
2403  WorkVec.Add(1, F);
2404  WorkVec.Sub(4, F);
2405 
2406  return WorkVec;
2407 }
2408 
2409 /* AbsoluteInternalCouple - end */
2410 
2411 
2412 /* FollowerInternalCouple - begin */
2413 
2414 FollowerInternalCouple::FollowerInternalCouple(unsigned int uL,
2415  const StructNode* pN1, const StructNode* pN2,
2416  const TplDriveCaller<Vec3>* pDC,
2417  flag fOut)
2418 : Elem(uL, fOut),
2419 StructuralInternalForce(uL, pN1, pN2, pDC, fOut)
2420 {
2421  NO_OP;
2422 }
2423 
2424 
2425 FollowerInternalCouple::~FollowerInternalCouple(void)
2426 {
2427  NO_OP;
2428 }
2429 
2430 void
2431 FollowerInternalCouple::WorkSpaceDim(
2432  integer* piNumRows,
2433  integer* piNumCols) const
2434 {
2435  *piNumRows = 6;
2436  *piNumCols = 3;
2437 }
2438 
2439 
2440 void
2441 FollowerInternalCouple::InitialWorkSpaceDim(
2442  integer* piNumRows,
2443  integer* piNumCols) const
2444 {
2445  *piNumRows = 12;
2446  *piNumCols = 12;
2447 }
2448 
2449 
2450 /* Contributo al file di restart */
2451 std::ostream&
2452 FollowerInternalCouple::Restart(std::ostream& out) const
2453 {
2454  out << " couple: " << GetLabel() << ", follower internal, "
2455  << pNode1->GetLabel() << ", "
2456  << pNode2->GetLabel() << ", ";
2457  return f.pGetDriveCaller()->Restart(out) << ';' << std::endl;
2458 }
2459 
2460 
2461 VariableSubMatrixHandler&
2462 FollowerInternalCouple::AssJac(VariableSubMatrixHandler& WorkMat,
2463  doublereal dCoef,
2464  const VectorHandler& /* XCurr */ ,
2465  const VectorHandler& /* XPrimeCurr */ )
2466 {
2467  DEBUGCOUT("Entering FollowerInternalCouple::AssJac()" << std::endl);
2468 
2469  FullSubMatrixHandler& WM = WorkMat.SetFull();
2470 
2471  /* Dimensiona e resetta la matrice di lavoro */
2472  integer iNumRows = 0;
2473  integer iNumCols = 0;
2474  WorkSpaceDim(&iNumRows, &iNumCols);
2475  WM.ResizeReset(iNumRows, iNumCols);
2476 
2477  integer iFirstRotationIndex1 = pNode1->iGetFirstPositionIndex() + 3;
2478  integer iFirstMomentumIndex1 = pNode1->iGetFirstMomentumIndex() + 3;
2479  integer iFirstMomentumIndex2 = pNode2->iGetFirstMomentumIndex() + 3;
2480  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
2481  WM.PutRowIndex(iCnt, iFirstMomentumIndex1 + iCnt); /* coppia */
2482  WM.PutColIndex(iCnt, iFirstRotationIndex1 + iCnt); /* rotazione */
2483 
2484  WM.PutRowIndex(3 + iCnt, iFirstMomentumIndex2 + iCnt); /* coppia */
2485  }
2486 
2487  /* Dati */
2488  Mat3x3 MWedge(MatCross, pNode1->GetRRef()*(f.Get()*dCoef));
2489 
2490  /* | M /\| Delta_g = | M | */
2491 
2492  WM.Add(1, 1, MWedge);
2493  WM.Sub(4, 1, MWedge);
2494 
2495  return WorkMat;
2496 }
2497 
2498 
2499 /* Assembla il residuo */
2500 SubVectorHandler&
2501 FollowerInternalCouple::AssRes(SubVectorHandler& WorkVec,
2502  doublereal /* dCoef */ ,
2503  const VectorHandler& /* XCurr */ ,
2504  const VectorHandler& /* XPrimeCurr */ )
2505 {
2506  DEBUGCOUT("Entering FollowerInternalCouple::AssRes()" << std::endl);
2507 
2508  integer iNumRows;
2509  integer iNumCols;
2510  WorkSpaceDim(&iNumRows, &iNumCols);
2511  WorkVec.ResizeReset(iNumRows);
2512 
2513  /* Indici delle incognite del nodo */
2514  integer iFirstMomentumIndex1 = pNode1->iGetFirstMomentumIndex() + 3;
2515  integer iFirstMomentumIndex2 = pNode2->iGetFirstMomentumIndex() + 3;
2516  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
2517  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex1 + iCnt);
2518  WorkVec.PutRowIndex(3 + iCnt, iFirstMomentumIndex2 + iCnt);
2519  }
2520 
2521  /* Dati */
2522  Vec3 M(pNode1->GetRCurr()*f.Get());
2523 
2524  WorkVec.Add(1, M);
2525  WorkVec.Sub(4, M);
2526 
2527  return WorkVec;
2528 }
2529 
2530 /* Inverse Dynamics*/
2531 SubVectorHandler&
2532 FollowerInternalCouple::AssRes(SubVectorHandler& WorkVec,
2533  const VectorHandler& /* XCurr */ ,
2534  const VectorHandler& /* XPrimeCurr */ ,
2535  const VectorHandler& /* XPrimePrimeCurr */ ,
2536  InverseDynamics::Order iOrder)
2537 {
2538  DEBUGCOUT("Entering FollowerInternalCouple::AssRes()" << std::endl);
2539 
2540  integer iNumRows;
2541  integer iNumCols;
2542  WorkSpaceDim(&iNumRows, &iNumCols);
2543  WorkVec.ResizeReset(iNumRows);
2544 
2545  /* Indici delle incognite del nodo */
2546  integer iFirstMomentumIndex1 = pNode1->iGetFirstPositionIndex() + 3;
2547  integer iFirstMomentumIndex2 = pNode2->iGetFirstPositionIndex() + 3;
2548  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
2549  WorkVec.PutRowIndex(iCnt, iFirstMomentumIndex1 + iCnt);
2550  WorkVec.PutRowIndex(3 + iCnt, iFirstMomentumIndex2 + iCnt);
2551  }
2552 
2553  /* Dati */
2554  Vec3 M(pNode1->GetRCurr()*f.Get());
2555 
2556  WorkVec.Add(1, M);
2557  WorkVec.Sub(4, M);
2558 
2559  return WorkVec;
2560 }
2561 
2562 void
2563 FollowerInternalCouple::OutputPrepare(OutputHandler& OH)
2564 {
2565  if (bToBeOutput()) {
2566 #ifdef USE_NETCDF
2567  if (OH.UseNetCDF(OutputHandler::FORCES)) {
2568  ASSERT(OH.IsOpen(OutputHandler::NETCDF));
2569 
2570  std::ostringstream os;
2571  os << "elem.couple." << GetLabel();
2572  (void)OH.CreateVar(os.str(), "internal follower");
2573 
2574  // joint sub-data
2575  os << '.';
2576  if (fToBeOutput() & StructuralForce::OUTPUT_REL) {
2577  Var_F = OH.CreateVar<Vec3>(os.str() + "m", "Nm",
2578  "local couple components (x, y, z)");
2579 
2580  } else {
2581  Var_F = OH.CreateVar<Vec3>(os.str() + "M", "Nm",
2582  "global couple components (x, y, z)");
2583  }
2584  }
2585 #endif // USE_NETCDF
2586  }
2587 }
2588 
2589 void
2590 FollowerInternalCouple::Output(OutputHandler& OH) const
2591 {
2592  if (bToBeOutput()) {
2593 #ifdef USE_NETCDF
2594  if (OH.UseNetCDF(OutputHandler::FORCES)) {
2595  if (fToBeOutput() & StructuralForce::OUTPUT_REL) {
2596  Var_F->put_rec(f.Get().pGetVec(), OH.GetCurrentStep());
2597  } else {
2598  Var_F->put_rec((pNode1->GetRCurr()*f.Get()).pGetVec(), OH.GetCurrentStep());
2599  }
2600  }
2601 #endif // USE_NETCDF
2602 
2603  if (OH.UseText(OutputHandler::FORCES)) {
2604  Vec3 F;
2605  if (fToBeOutput() & StructuralForce::OUTPUT_REL) {
2606  F = f.Get();
2607  } else {
2608  F = pNode1->GetRCurr()*f.Get();
2609  }
2610 
2611  OH.Forces()
2612  << GetLabel()
2613  << " " << pNode1->GetLabel()
2614  << " " << F
2615  << " " << pNode2->GetLabel()
2616  << " " << -F
2617  << std::endl;
2618  }
2619 
2620  /* TODO: NetCDF */
2621  }
2622 }
2623 
2624 
2625 /* Contributo allo jacobiano durante l'assemblaggio iniziale */
2626 VariableSubMatrixHandler&
2627 FollowerInternalCouple::InitialAssJac(VariableSubMatrixHandler& WorkMat,
2628  const VectorHandler& /* XCurr */ )
2629 {
2630  DEBUGCOUT("Entering FollowerInternalCouple::InitialAssJac()" << std::endl);
2631 
2632  FullSubMatrixHandler& WM = WorkMat.SetFull();
2633 
2634  /* Dimensiona e resetta la matrice di lavoro */
2635  WM.ResizeReset(12, 6);
2636 
2637  integer iFirstPositionIndex1 = pNode1->iGetFirstPositionIndex() + 3;
2638  integer iFirstVelocityIndex1 = iFirstPositionIndex1 + 6;
2639  integer iFirstPositionIndex2 = pNode2->iGetFirstPositionIndex() + 3;
2640  integer iFirstVelocityIndex2 = iFirstPositionIndex2 + 6;
2641  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
2642  WM.PutRowIndex(iCnt, iFirstPositionIndex1 + iCnt);
2643  WM.PutRowIndex(3 + iCnt, iFirstVelocityIndex1 + iCnt);
2644 
2645  WM.PutColIndex(iCnt, iFirstPositionIndex1 + iCnt);
2646  WM.PutColIndex(3 + iCnt, iFirstVelocityIndex1 + iCnt);
2647 
2648  WM.PutRowIndex(6 + iCnt, iFirstPositionIndex2 + iCnt);
2649  WM.PutRowIndex(9 + iCnt, iFirstVelocityIndex2 + iCnt);
2650  }
2651 
2652  /* Dati */
2653  Vec3 TmpDir(pNode1->GetRRef()*f.Get());
2654  const Vec3& Omega1(pNode1->GetWRef());
2655 
2656  /* | F/\ | | F |
2657  * | | Delta_g = | |
2658  * | (d/\F)/\ | | d/\F |
2659  */
2660 
2661  WM.Add(1, 1, Mat3x3(MatCross, TmpDir));
2662  WM.Add(4, 1, Mat3x3(MatCrossCross, Omega1, TmpDir));
2663  WM.Add(4, 4, Mat3x3(MatCross, TmpDir));
2664 
2665  WM.Sub(7, 1, Mat3x3(MatCross, TmpDir));
2666  WM.Sub(10, 1, Mat3x3(MatCrossCross, Omega1, TmpDir));
2667  WM.Sub(10, 4, Mat3x3(MatCross, TmpDir));
2668 
2669  return WorkMat;
2670 }
2671 
2672 
2673 /* Contributo al residuo durante l'assemblaggio iniziale */
2674 SubVectorHandler&
2675 FollowerInternalCouple::InitialAssRes(SubVectorHandler& WorkVec,
2676  const VectorHandler& /* XCurr */ )
2677 {
2678  DEBUGCOUT("Entering FollowerInternalCouple::InitialAssRes()" << std::endl);
2679 
2680  WorkVec.ResizeReset(12);
2681 
2682  /* Indici delle incognite del nodo */
2683  integer iFirstPositionIndex1 = pNode1->iGetFirstPositionIndex() + 3;
2684  integer iFirstVelocityIndex1 = iFirstPositionIndex1 + 6;
2685 
2686  integer iFirstPositionIndex2 = pNode2->iGetFirstPositionIndex() + 3;
2687  integer iFirstVelocityIndex2 = iFirstPositionIndex2 + 6;
2688 
2689  for (integer iCnt = 1; iCnt <= 3; iCnt++) {
2690  WorkVec.PutRowIndex(iCnt, iFirstPositionIndex1+iCnt);
2691  WorkVec.PutRowIndex(3 + iCnt, iFirstVelocityIndex1 + iCnt);
2692 
2693  WorkVec.PutRowIndex(6 + iCnt, iFirstPositionIndex2 + iCnt);
2694  WorkVec.PutRowIndex(9 + iCnt, iFirstVelocityIndex2 + iCnt);
2695  }
2696 
2697  /* Dati */
2698  Vec3 TmpDir(pNode1->GetRCurr()*f.Get());
2699  const Vec3& Omega1(pNode1->GetWCurr());
2700 
2701  WorkVec.Add(1, TmpDir);
2702  WorkVec.Add(4, Omega1.Cross(TmpDir));
2703 
2704  WorkVec.Sub(7, TmpDir);
2705  WorkVec.Sub(10, Omega1.Cross(TmpDir));
2706 
2707  return WorkVec;
2708 }
2709 
2710 /* FollowerInternalCouple - end */
2711 
2712 Elem *
2713 ReadStructuralForce(DataManager* pDM,
2714  MBDynParser& HP,
2715  unsigned int uLabel,
2716  bool bDisp,
2717  bool bCouple,
2718  bool bFollower,
2719  bool bInternal)
2720 {
2721  Elem *pEl = 0;
2722  const char *sType = bCouple ? "Couple" : "Force";
2723 
2724  /* nodo collegato */
2725  const StructDispNode* pDispNode = pDM->ReadNode<const StructDispNode, Node::STRUCTURAL>(HP);
2726  const StructNode* pNode = dynamic_cast<const StructNode *>(pDispNode);
2727  ReferenceFrame rf;
2728  if (pNode) {
2729  rf = ReferenceFrame(pNode);
2730  }
2731  Vec3 Arm(Zero3);
2732 
2733  // FIXME: legacy...
2734  Vec3 Dir(Zero3);
2735  bool bLegacy(false);
2736  bool bGotPosition(false);
2737 
2738  /* distanza dal nodo (vettore di 3 elementi) (solo se e' una forza) */
2739  if (!bDisp) {
2740  if (pNode == 0) {
2741  silent_cerr(sType << "(" << uLabel << "): "
2742  "invalid node type at line " << HP.GetLineData() << std::endl);
2743  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
2744  }
2745 
2746  if (HP.IsKeyWord("position")) {
2747  Arm = HP.GetPosRel(rf);
2748  bGotPosition = true;
2749  DEBUGCOUT("Arm is supplied" << std::endl);
2750 
2751  } else {
2752  if (!bCouple) {
2753  silent_cerr(sType << "(" << uLabel << "): "
2754  "\"position\" keyword expected "
2755  "at line " << HP.GetLineData() << "; "
2756  "still using deprecated syntax?"
2757  << std::endl);
2758 
2759  if (bFollower) {
2760  try {
2761  Dir = HP.GetUnitVecRel(rf);
2762  } catch (ErrNullNorm) {
2763  silent_cerr(sType << "(" << uLabel << ") has null direction" << std::endl);
2764  throw ErrNullNorm(MBDYN_EXCEPT_ARGS);
2765  }
2766 
2767  } else {
2768  try {
2769  Dir = HP.GetUnitVecAbs(rf);
2770  } catch (ErrNullNorm) {
2771  silent_cerr(sType << "(" << uLabel << ") has null direction" << std::endl);
2772  throw ErrNullNorm(MBDYN_EXCEPT_ARGS);
2773  }
2774  }
2775 
2776  Arm = HP.GetPosRel(rf);
2777 
2778  bLegacy = true;
2779  }
2780  }
2781  }
2782 
2783  const StructDispNode *pDispNode2 = 0;
2784  const StructNode *pNode2 = 0;
2785  Vec3 Arm2(Zero3);
2786  if (bInternal) {
2787  /* nodo collegato */
2788  pDispNode2 = pDM->ReadNode<const StructDispNode, Node::STRUCTURAL>(HP);
2789  pNode2 = dynamic_cast<const StructNode *>(pDispNode2);
2790 
2791  if (!bDisp) {
2792  if (pNode2 == 0) {
2793  silent_cerr(sType << "(" << uLabel << "): "
2794  "invalid node type at line " << HP.GetLineData() << std::endl);
2795  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
2796  }
2797 
2798  ReferenceFrame rf2(pNode2);
2799 
2800  /* distanza dal nodo (vettore di 3 elementi) ( solo se e' una forza) */
2801  if (HP.IsKeyWord("position")) {
2802  Arm2 = HP.GetPosRel(rf2, rf, Arm);
2803  bGotPosition = true;
2804  DEBUGCOUT("Node 2 arm is supplied" << std::endl);
2805 
2806  } else if (bLegacy) {
2807  Arm2 = HP.GetPosRel(rf2, rf, Arm);
2808  }
2809  }
2810  }
2811 
2812  if (bCouple && bInternal && !bGotPosition) {
2813  silent_cerr(sType << "(" << uLabel << ") "
2814  "line " << HP.GetLineData() << ": "
2815  "warning, the syntax changed; "
2816  "you may safely ignore this warning if you used "
2817  "the syntax documented for MBDyn >= 1.3.7"
2818  << std::endl);
2819  }
2820 
2821  TplDriveCaller<Vec3>* pDC = 0;
2822  if (bLegacy) {
2823  pDC = DC2TDC(HP.GetDriveCaller(), Dir);
2824 
2825  } else {
2826  if (bFollower) {
2827  pDC = ReadDCVecRel(pDM, HP, rf);
2828 
2829  } else {
2830  pDC = ReadDCVecAbs(pDM, HP, rf);
2831  }
2832  }
2833 
2834  bool bOutRel(false);
2835  if (bFollower && HP.IsKeyWord("output" "relative")) {
2836  bOutRel = true;
2837  }
2838 
2839  flag fOut = pDM->fReadOutput(HP, Elem::FORCE);
2840 
2841  if (bOutRel) {
2842  fOut |= StructuralForce::OUTPUT_REL;
2843  }
2844 
2845  /* Alloca la forza */
2846  if (!bCouple) {
2847  if (!bFollower) {
2848  if (!bInternal) {
2849  if (bDisp) {
2850  SAFENEWWITHCONSTRUCTOR(pEl,
2851  AbsoluteDispForce,
2852  AbsoluteDispForce(uLabel, pDispNode, pDC, fOut));
2853 
2854  } else {
2855  SAFENEWWITHCONSTRUCTOR(pEl,
2856  AbsoluteForce,
2857  AbsoluteForce(uLabel, pNode, pDC, Arm, fOut));
2858  }
2859 
2860  } else {
2861  if (bDisp) {
2862  SAFENEWWITHCONSTRUCTOR(pEl,
2863  AbsoluteInternalDispForce,
2864  AbsoluteInternalDispForce(uLabel, pDispNode, pDispNode2, pDC, fOut));
2865 
2866  } else {
2867  SAFENEWWITHCONSTRUCTOR(pEl,
2868  AbsoluteInternalForce,
2869  AbsoluteInternalForce(uLabel, pNode, pNode2, pDC, Arm, Arm2, fOut));
2870  }
2871  }
2872 
2873  } else {
2874  if (!bInternal) {
2875  SAFENEWWITHCONSTRUCTOR(pEl,
2876  FollowerForce,
2877  FollowerForce(uLabel, pNode, pDC, Arm, fOut));
2878 
2879  } else {
2880  SAFENEWWITHCONSTRUCTOR(pEl,
2881  FollowerInternalForce,
2882  FollowerInternalForce(uLabel, pNode, pNode2, pDC, Arm, Arm2, fOut));
2883  }
2884  }
2885 
2886  } else {
2887  if (!bFollower) {
2888  if (!bInternal) {
2889  SAFENEWWITHCONSTRUCTOR(pEl,
2890  AbsoluteCouple,
2891  AbsoluteCouple(uLabel, pNode, pDC, fOut));
2892 
2893  } else {
2894  SAFENEWWITHCONSTRUCTOR(pEl,
2895  AbsoluteInternalCouple,
2896  AbsoluteInternalCouple(uLabel, pNode, pNode2, pDC, fOut));
2897  }
2898 
2899  } else {
2900  if (!bInternal) {
2901  SAFENEWWITHCONSTRUCTOR(pEl,
2902  FollowerCouple,
2903  FollowerCouple(uLabel, pNode, pDC, fOut));
2904 
2905  } else {
2906  SAFENEWWITHCONSTRUCTOR(pEl,
2907  FollowerInternalCouple,
2908  FollowerInternalCouple(uLabel, pNode, pNode2, pDC, fOut));
2909  }
2910  }
2911  }
2912 
2913  std::ostream& os = pDM->GetLogFile();
2914 
2915  os << "structural";
2916  if (bInternal) {
2917  os << " internal";
2918  }
2919  if (bFollower) {
2920  os << " follower";
2921  } else {
2922  os << " absolute";
2923  if (bDisp) {
2924  os << " displacement";
2925  }
2926  }
2927  if (bCouple) {
2928  os << " couple";
2929  } else {
2930  os << " force";
2931  }
2932  os << ": " << uLabel << ' ' << pDispNode->GetLabel()
2933  << ' ' << Arm;
2934 
2935  if (pDispNode2 != 0) {
2936  os << ' ' << pDispNode2->GetLabel() << ' ' << Arm2;
2937  }
2938 
2939  if (bFollower && (pEl->fToBeOutput() & StructuralForce::OUTPUT_REL)) {
2940  os << " output_relative";
2941  }
2942 
2943  os << std::endl;
2944 
2945  return pEl;
2946 }
2947 
AbsoluteInternalDispForce::Restart
virtual std::ostream & Restart(std::ostream &out) const
Definition: strforce.cc:257
DataManager::fReadOutput
flag fReadOutput(MBDynParser &HP, const T &t) const
Definition: dataman.h:1064
FollowerForce::AssJac
virtual VariableSubMatrixHandler & AssJac(VariableSubMatrixHandler &WorkMat, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:762
FollowerInternalCouple::OutputPrepare
virtual void OutputPrepare(OutputHandler &OH)
Definition: strforce.cc:2563
AbsoluteCouple::OutputPrepare
virtual void OutputPrepare(OutputHandler &OH)
Definition: strforce.cc:1119
StructuralInternalForce::f
TplDriveOwner< Vec3 > f
Definition: strforce_impl.h:461
FollowerForce::OutputPrepare
virtual void OutputPrepare(OutputHandler &OH)
Definition: strforce.cc:869
FollowerForce::Output
virtual void Output(OutputHandler &OH) const
Definition: strforce.cc:900
FullSubMatrixHandler::PutColIndex
void PutColIndex(integer iSubCol, integer iCol)
Definition: submat.h:325
FollowerCouple::FollowerCouple
FollowerCouple(unsigned int uL, const StructNode *pN, const TplDriveCaller< Vec3 > *pDC, flag fOut)
Definition: strforce.cc:1188
FollowerCouple::OutputPrepare
virtual void OutputPrepare(OutputHandler &OH)
Definition: strforce.cc:1324
Zero3
const Vec3 Zero3(0., 0., 0.)
Vec3::Cross
Vec3 Cross(const Vec3 &v) const
Definition: matvec3.h:218
AbsoluteInternalForce::Restart
virtual std::ostream & Restart(std::ostream &out) const
Definition: strforce.cc:1544
Vec3::Write
std::ostream & Write(std::ostream &out, const char *sFill=" ") const
Definition: matvec3.cc:738
FollowerInternalCouple::InitialAssJac
virtual VariableSubMatrixHandler & InitialAssJac(VariableSubMatrixHandler &WorkMat, const VectorHandler &XCurr)
Definition: strforce.cc:2627
flag
long int flag
Definition: mbdyn.h:43
StructNode::GetRRef
virtual const Mat3x3 & GetRRef(void) const
Definition: strnode.h:1006
ToBeOutput::bToBeOutput
virtual bool bToBeOutput(void) const
Definition: output.cc:890
MBDYN_EXCEPT_ARGS
#define MBDYN_EXCEPT_ARGS
Definition: except.h:63
AbsoluteInternalDispForce::AbsoluteInternalDispForce
AbsoluteInternalDispForce(unsigned int uL, const StructDispNode *pN1, const StructDispNode *pN2, const TplDriveCaller< Vec3 > *pDC, flag fOut)
Definition: strforce.cc:214
Vec3
Definition: matvec3.h:98
AbsoluteCouple::InitialWorkSpaceDim
virtual void InitialWorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:1048
FollowerInternalCouple::WorkSpaceDim
void WorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:2431
VectorHandler::ResizeReset
virtual void ResizeReset(integer)
Definition: vh.cc:55
tpldrive_impl.h
MatCross
const MatCross_Manip MatCross
Definition: matvec3.cc:639
OutputHandler::UseNetCDF
bool UseNetCDF(int out) const
Definition: output.cc:491
FollowerCouple::AssRes
virtual SubVectorHandler & AssRes(SubVectorHandler &WorkVec, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:1269
VariableSubMatrixHandler::SetFull
FullSubMatrixHandler & SetFull(void)
Definition: submat.h:1168
FollowerInternalForce::InitialAssRes
virtual SubVectorHandler & InitialAssRes(SubVectorHandler &WorkVec, const VectorHandler &XCurr)
Definition: strforce.cc:2172
StructNode::GetRCurr
virtual const Mat3x3 & GetRCurr(void) const
Definition: strnode.h:1012
FollowerForce::WorkSpaceDim
void WorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:732
StructDispNode::GetNodeType
virtual Node::Type GetNodeType(void) const
Definition: strnode.cc:145
AbsoluteInternalDispForce::WorkSpaceDim
void WorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:240
StructuralInternalForce::pNode1
const StructNode * pNode1
Definition: strforce_impl.h:462
FullSubMatrixHandler::Add
void Add(integer iRow, integer iCol, const Vec3 &v)
Definition: submat.cc:209
AbsoluteInternalDispForce::Output
virtual void Output(OutputHandler &OH) const
Definition: strforce.cc:330
AbsoluteForce::OutputPrepare
virtual void OutputPrepare(OutputHandler &OH)
Definition: strforce.cc:579
AbsoluteForce::AssRes
virtual SubVectorHandler & AssRes(SubVectorHandler &WorkVec, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:519
VectorHandler::Sub
virtual void Sub(integer iRow, const Vec3 &v)
Definition: vh.cc:78
FollowerCouple::AssJac
virtual VariableSubMatrixHandler & AssJac(VariableSubMatrixHandler &WorkMat, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:1233
StructuralInternalForce::~StructuralInternalForce
virtual ~StructuralInternalForce(void)
Definition: strforce.cc:1478
AbsoluteDispForce::AssRes
virtual SubVectorHandler & AssRes(SubVectorHandler &WorkVec, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:97
FollowerForce::FollowerForce
FollowerForce(unsigned int uL, const StructNode *pN, const TplDriveCaller< Vec3 > *pDC, const Vec3 &TmpArm, flag fOut)
Definition: strforce.cc:709
Force
Definition: force.h:46
FollowerInternalForce::OutputPrepare
virtual void OutputPrepare(OutputHandler &OH)
Definition: strforce.cc:2026
StructuralForce::GetConnectedNodes
virtual void GetConnectedNodes(std::vector< const Node * > &connectedNodes) const
Definition: strforce.cc:423
AbsoluteInternalDispForce::AssRes
virtual SubVectorHandler & AssRes(SubVectorHandler &WorkVec, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:268
FollowerCouple::WorkSpaceDim
void WorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:1205
Force::Restart
virtual std::ostream & Restart(std::ostream &out) const
Definition: force.cc:52
FollowerInternalForce::AssJac
virtual VariableSubMatrixHandler & AssJac(VariableSubMatrixHandler &WorkMat, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:1903
StructNode::GetWRef
virtual const Vec3 & GetWRef(void) const
Definition: strnode.h:1024
NO_OP
#define NO_OP
Definition: myassert.h:74
AbsoluteDispForce::pNode
const StructDispNode * pNode
Definition: strforce_impl.h:46
FollowerInternalForce::Output
virtual void Output(OutputHandler &OH) const
Definition: strforce.cc:2058
FollowerInternalCouple::FollowerInternalCouple
FollowerInternalCouple(unsigned int uL, const StructNode *pN1, const StructNode *pN2, const TplDriveCaller< Vec3 > *pDC, flag fOut)
Definition: strforce.cc:2414
FollowerCouple::InitialAssRes
virtual SubVectorHandler & InitialAssRes(SubVectorHandler &WorkVec, const VectorHandler &XCurr)
Definition: strforce.cc:1425
AbsoluteInternalForce::OutputPrepare
virtual void OutputPrepare(OutputHandler &OH)
Definition: strforce.cc:1671
FollowerInternalForce::FollowerInternalForce
FollowerInternalForce(unsigned int uL, const StructNode *pN1, const StructNode *pN2, const TplDriveCaller< Vec3 > *pDC, const Vec3 &TmpArm1, const Vec3 &TmpArm2, flag fOut)
Definition: strforce.cc:1843
strforce_impl.h
FollowerInternalCouple::InitialWorkSpaceDim
virtual void InitialWorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:2441
FollowerForce::AssRes
virtual SubVectorHandler & AssRes(SubVectorHandler &WorkVec, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:804
AbsoluteForce::Output
virtual void Output(OutputHandler &OH) const
Definition: strforce.cc:603
FollowerInternalCouple::Restart
virtual std::ostream & Restart(std::ostream &out) const
Definition: strforce.cc:2452
AbsoluteInternalDispForce::pNode2
const StructDispNode * pNode2
Definition: strforce_impl.h:111
StructuralInternalForce::StructuralInternalForce
StructuralInternalForce(unsigned int uL, const StructNode *pN1, const StructNode *pN2, const TplDriveCaller< Vec3 > *pDC, flag fOut)
Definition: strforce.cc:1457
FollowerForce::InitialWorkSpaceDim
virtual void InitialWorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:740
AbsoluteInternalCouple::AssRes
virtual SubVectorHandler & AssRes(SubVectorHandler &WorkVec, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:2275
AbsoluteDispForce::Restart
virtual std::ostream & Restart(std::ostream &out) const
Definition: strforce.cc:87
AbsoluteDispForce::f
TplDriveOwner< Vec3 > f
Definition: strforce_impl.h:45
ReadDCVecAbs
TplDriveCaller< Vec3 > * ReadDCVecAbs(const DataManager *pDM, MBDynParser &HP, const ReferenceFrame &rf)
Definition: tpldrive_impl.cc:923
SubVectorHandler::PutRowIndex
virtual void PutRowIndex(integer iSubRow, integer iRow)=0
FollowerForce::Arm
const Vec3 Arm
Definition: strforce_impl.h:280
AbsoluteCouple::AssRes
virtual SubVectorHandler & AssRes(SubVectorHandler &WorkVec, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:1069
AbsoluteForce::AssJac
virtual VariableSubMatrixHandler & AssJac(VariableSubMatrixHandler &WorkMat, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:484
AbsoluteCouple::~AbsoluteCouple
~AbsoluteCouple(void)
Definition: strforce.cc:1033
FollowerCouple::~FollowerCouple
~FollowerCouple(void)
Definition: strforce.cc:1198
StructuralForce::~StructuralForce
virtual ~StructuralForce(void)
Definition: strforce.cc:417
FollowerInternalForce::InitialWorkSpaceDim
virtual void InitialWorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:1878
StructuralForce::StructuralForce
StructuralForce(unsigned int uL, const StructNode *pN, const TplDriveCaller< Vec3 > *pDC, flag fOut)
Definition: strforce.cc:398
FollowerForce::InitialAssRes
virtual SubVectorHandler & InitialAssRes(SubVectorHandler &WorkVec, const VectorHandler &XCurr)
Definition: strforce.cc:983
TplDriveOwner::pGetDriveCaller
TplDriveCaller< T > * pGetDriveCaller(void) const
Definition: tpldrive.h:105
AbsoluteForce::Restart
virtual std::ostream & Restart(std::ostream &out) const
Definition: strforce.cc:475
AbsoluteDispForce::WorkSpaceDim
void WorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:70
AbsoluteInternalDispForce::GetConnectedNodes
virtual void GetConnectedNodes(std::vector< const Node * > &connectedNodes) const
Definition: strforce.cc:385
AbsoluteInternalForce::WorkSpaceDim
void WorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:1523
MBDynParser::GetPosRel
Vec3 GetPosRel(const ReferenceFrame &rf)
Definition: mbpar.cc:1331
FollowerForce::Restart
virtual std::ostream & Restart(std::ostream &out) const
Definition: strforce.cc:751
MBDynParser::pDM
DataManager * pDM
Definition: mbpar.h:252
HighParser::IsKeyWord
virtual bool IsKeyWord(const char *sKeyWord)
Definition: parser.cc:910
AbsoluteDispForce::InitialWorkSpaceDim
virtual void InitialWorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:77
FollowerForce::~FollowerForce
~FollowerForce(void)
Definition: strforce.cc:725
OutputHandler::GetCurrentStep
long GetCurrentStep(void) const
Definition: output.h:116
AbsoluteForce::InitialWorkSpaceDim
virtual void InitialWorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:465
DEBUGCOUT
#define DEBUGCOUT(msg)
Definition: myassert.h:232
TplDriveOwner::Get
T Get(const doublereal &dVar) const
Definition: tpldrive.h:109
AbsoluteInternalForce::AssRes
virtual SubVectorHandler & AssRes(SubVectorHandler &WorkVec, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:1601
StructuralForce::f
TplDriveOwner< Vec3 > f
Definition: strforce_impl.h:180
AbsoluteDispForce::Output
virtual void Output(OutputHandler &OH) const
Definition: strforce.cc:156
StructuralForce::pNode
const StructNode * pNode
Definition: strforce_impl.h:181
AbsoluteForce::Arm
const Vec3 Arm
Definition: strforce_impl.h:211
AbsoluteDispForce::GetConnectedNodes
virtual void GetConnectedNodes(std::vector< const Node * > &connectedNodes) const
Definition: strforce.cc:201
StructDispNode::iGetFirstMomentumIndex
virtual integer iGetFirstMomentumIndex(void) const =0
AbsoluteCouple::WorkSpaceDim
void WorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:1040
StructDispNode::iGetFirstPositionIndex
virtual integer iGetFirstPositionIndex(void) const
Definition: strnode.h:452
FollowerInternalCouple::AssRes
virtual SubVectorHandler & AssRes(SubVectorHandler &WorkVec, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:2501
StructNode::GetWCurr
virtual const Vec3 & GetWCurr(void) const
Definition: strnode.h:1030
AbsoluteInternalForce::InitialWorkSpaceDim
virtual void InitialWorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:1533
AbsoluteForce::InitialAssJac
virtual VariableSubMatrixHandler & InitialAssJac(VariableSubMatrixHandler &WorkMat, const VectorHandler &XCurr)
Definition: strforce.cc:629
AbsoluteCouple::AbsoluteCouple
AbsoluteCouple(unsigned int uL, const StructNode *pN, const TplDriveCaller< Vec3 > *pDC, flag fOut)
Definition: strforce.cc:1023
AbsoluteInternalCouple::Output
virtual void Output(OutputHandler &OH) const
Definition: strforce.cc:2358
AbsoluteInternalForce::AssJac
virtual VariableSubMatrixHandler & AssJac(VariableSubMatrixHandler &WorkMat, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:1556
OutputHandler::IsOpen
bool IsOpen(int out) const
Definition: output.cc:395
TplDriveCaller::Restart
virtual std::ostream & Restart(std::ostream &out) const =0
FollowerInternalForce::Restart
virtual std::ostream & Restart(std::ostream &out) const
Definition: strforce.cc:1889
FollowerCouple::InitialAssJac
virtual VariableSubMatrixHandler & InitialAssJac(VariableSubMatrixHandler &WorkMat, const VectorHandler &XCurr)
Definition: strforce.cc:1387
ASSERT
#define ASSERT(expression)
Definition: colamd.c:977
AbsoluteForce::InitialAssRes
virtual SubVectorHandler & InitialAssRes(SubVectorHandler &WorkVec, const VectorHandler &XCurr)
Definition: strforce.cc:670
FollowerInternalCouple::AssJac
virtual VariableSubMatrixHandler & AssJac(VariableSubMatrixHandler &WorkMat, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:2462
DC2TDC
TplDriveCaller< T > * DC2TDC(DriveCaller *pDC, T &t)
Definition: tpldrive_impl.h:215
grad::Cross
VectorExpression< VectorCrossExpr< VectorLhsExpr, VectorRhsExpr >, 3 > Cross(const VectorExpression< VectorLhsExpr, 3 > &u, const VectorExpression< VectorRhsExpr, 3 > &v)
Definition: matvec.h:3248
SAFENEWWITHCONSTRUCTOR
#define SAFENEWWITHCONSTRUCTOR(pnt, item, constructor)
Definition: mynewmem.h:698
MBDynParser::GetUnitVecAbs
Vec3 GetUnitVecAbs(const ReferenceFrame &rf)
Definition: mbpar.cc:1716
StructDispNode::GetXCurr
virtual const Vec3 & GetXCurr(void) const
Definition: strnode.h:310
VectorHandler::Add
virtual void Add(integer iRow, const Vec3 &v)
Definition: vh.cc:63
FollowerCouple::Restart
virtual std::ostream & Restart(std::ostream &out) const
Definition: strforce.cc:1224
AbsoluteDispForce::InitialAssRes
virtual SubVectorHandler & InitialAssRes(SubVectorHandler &WorkVec, const VectorHandler &XCurr)
Definition: strforce.cc:179
FullSubMatrixHandler::ResizeReset
virtual void ResizeReset(integer, integer)
Definition: submat.cc:182
AbsoluteInternalCouple::WorkSpaceDim
void WorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:2243
AbsoluteInternalForce::AbsoluteInternalForce
AbsoluteInternalForce(unsigned int uL, const StructNode *pN1, const StructNode *pN2, const TplDriveCaller< Vec3 > *pDC, const Vec3 &TmpArm1, const Vec3 &TmpArm2, flag fOut)
Definition: strforce.cc:1498
AbsoluteInternalForce::InitialAssJac
virtual VariableSubMatrixHandler & InitialAssJac(VariableSubMatrixHandler &WorkMat, const VectorHandler &XCurr)
Definition: strforce.cc:1729
AbsoluteInternalCouple::AbsoluteInternalCouple
AbsoluteInternalCouple(unsigned int uL, const StructNode *pN1, const StructNode *pN2, const TplDriveCaller< Vec3 > *pDC, flag fOut)
Definition: strforce.cc:2226
MatCrossCross
const MatCrossCross_Manip MatCrossCross
Definition: matvec3.cc:640
AbsoluteInternalCouple::OutputPrepare
virtual void OutputPrepare(OutputHandler &OH)
Definition: strforce.cc:2337
Elem
Definition: elem.h:75
AbsoluteInternalCouple::InitialWorkSpaceDim
virtual void InitialWorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:2253
AbsoluteDispForce::~AbsoluteDispForce
~AbsoluteDispForce(void)
Definition: strforce.cc:64
FullSubMatrixHandler::PutRowIndex
void PutRowIndex(integer iSubRow, integer iRow)
Definition: submat.h:311
FollowerInternalForce::AssRes
virtual SubVectorHandler & AssRes(SubVectorHandler &WorkVec, doublereal dCoef, const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
Definition: strforce.cc:1954
FollowerInternalForce::InitialAssJac
virtual VariableSubMatrixHandler & InitialAssJac(VariableSubMatrixHandler &WorkMat, const VectorHandler &XCurr)
Definition: strforce.cc:2099
Vec3::pGetVec
const doublereal * pGetVec(void) const
Definition: matvec3.h:192
DataManager::GetLogFile
std::ostream & GetLogFile(void) const
Definition: dataman.h:326
ToBeOutput::fToBeOutput
virtual flag fToBeOutput(void) const
Definition: output.cc:884
AbsoluteInternalForce::InitialAssRes
virtual SubVectorHandler & InitialAssRes(SubVectorHandler &WorkVec, const VectorHandler &XCurr)
Definition: strforce.cc:1788
AbsoluteForce::~AbsoluteForce
~AbsoluteForce(void)
Definition: strforce.cc:452
AbsoluteInternalDispForce::OutputPrepare
virtual void OutputPrepare(OutputHandler &OH)
Definition: strforce.cc:309
StructuralInternalForce::pNode2
const StructNode * pNode2
Definition: strforce_impl.h:463
AbsoluteInternalDispForce::InitialWorkSpaceDim
virtual void InitialWorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:247
AbsoluteInternalDispForce::f
TplDriveOwner< Vec3 > f
Definition: strforce_impl.h:109
AbsoluteCouple::InitialAssRes
virtual SubVectorHandler & InitialAssRes(SubVectorHandler &WorkVec, const VectorHandler &XCurr)
Definition: strforce.cc:1165
FollowerInternalForce::WorkSpaceDim
void WorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:1868
AbsoluteInternalCouple::Restart
virtual std::ostream & Restart(std::ostream &out) const
Definition: strforce.cc:2264
FullSubMatrixHandler::Sub
void Sub(integer iRow, integer iCol, const Vec3 &v)
Definition: submat.cc:215
AbsoluteDispForce::AbsoluteDispForce
AbsoluteDispForce(unsigned int uL, const StructDispNode *pN, const TplDriveCaller< Vec3 > *pDC, flag fOut)
Definition: strforce.cc:47
MBDynParser::GetDriveCaller
DriveCaller * GetDriveCaller(bool bDeferred=false)
Definition: mbpar.cc:2033
AbsoluteCouple::Output
virtual void Output(OutputHandler &OH) const
Definition: strforce.cc:1141
MBDynParser::GetUnitVecRel
Vec3 GetUnitVecRel(const ReferenceFrame &rf)
Definition: mbpar.cc:1695
AbsoluteInternalDispForce::InitialAssRes
virtual SubVectorHandler & InitialAssRes(SubVectorHandler &WorkVec, const VectorHandler &XCurr)
Definition: strforce.cc:359
FollowerInternalCouple::InitialAssRes
virtual SubVectorHandler & InitialAssRes(SubVectorHandler &WorkVec, const VectorHandler &XCurr)
Definition: strforce.cc:2675
ReadStructuralForce
Elem * ReadStructuralForce(DataManager *pDM, MBDynParser &HP, unsigned int uLabel, bool bDisp, bool bCouple, bool bFollower, bool bInternal)
Definition: strforce.cc:2713
AbsoluteDispForce::OutputPrepare
virtual void OutputPrepare(OutputHandler &OH)
Definition: strforce.cc:135
AbsoluteCouple::Restart
virtual std::ostream & Restart(std::ostream &out) const
Definition: strforce.cc:1059
AbsoluteInternalForce::Output
virtual void Output(OutputHandler &OH) const
Definition: strforce.cc:1698
doublereal
double doublereal
Definition: colamd.c:52
FollowerCouple::InitialWorkSpaceDim
virtual void InitialWorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:1213
AbsoluteForce::AbsoluteForce
AbsoluteForce(unsigned int uL, const StructNode *pN, const TplDriveCaller< Vec3 > *pDC, const Vec3 &TmpArm, flag fOut)
Definition: strforce.cc:435
integer
long int integer
Definition: colamd.c:51
IncludeParser::GetLineData
virtual HighParser::ErrOut GetLineData(void) const
Definition: parsinc.cc:697
WithLabel::GetLabel
unsigned int GetLabel(void) const
Definition: withlab.cc:62
AbsoluteInternalDispForce::pNode1
const StructDispNode * pNode1
Definition: strforce_impl.h:110
DataManager::ReadNode
Node * ReadNode(MBDynParser &HP, Node::Type type) const
Definition: dataman3.cc:2309
AbsoluteInternalCouple::InitialAssRes
virtual SubVectorHandler & InitialAssRes(SubVectorHandler &WorkVec, const VectorHandler &XCurr)
Definition: strforce.cc:2385
StructuralInternalForce::GetConnectedNodes
virtual void GetConnectedNodes(std::vector< const Node * > &connectedNodes) const
Definition: strforce.cc:1484
AbsoluteForce::WorkSpaceDim
void WorkSpaceDim(integer *piNumRows, integer *piNumCols) const
Definition: strforce.cc:458
FollowerForce::InitialAssJac
virtual VariableSubMatrixHandler & InitialAssJac(VariableSubMatrixHandler &WorkMat, const VectorHandler &XCurr)
Definition: strforce.cc:938
VectorHandler::Resize
virtual void Resize(integer iNewSize)=0
FollowerInternalCouple::Output
virtual void Output(OutputHandler &OH) const
Definition: strforce.cc:2590
OutputHandler::Forces
std::ostream & Forces(void) const
Definition: output.h:450
FollowerCouple::Output
virtual void Output(OutputHandler &OH) const
Definition: strforce.cc:1353
R
Mat3x3 R
Definition: mbdyn_util_oct.cc:74
ReadDCVecRel
TplDriveCaller< Vec3 > * ReadDCVecRel(const DataManager *pDM, MBDynParser &HP, const ReferenceFrame &rf)
Definition: tpldrive_impl.cc:917
OutputHandler::UseText
bool UseText(int out) const
Definition: output.cc:446