1 /* $Header: /var/cvs/mbdyn/mbdyn/mbdyn-1.0/mbdyn/struct/strnode.h,v 1.82 2017/01/12 14:46:44 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 #ifndef STRNODE_H
33 #define STRNODE_H
34 
35 #include "node.h"
36 #include "matvec3.h"
37 #include "rbk.h"
38 #include "invdyn.h"
39 
40 #ifdef USE_AUTODIFF
41 #ifdef USE_MULTITHREAD
42 #include "veciter.h"
43 #endif
44 #include "gradient.h"
45 #include "matvec.h"
46 #endif
47 
48 extern const char* psStructNodeNames[];
49 
50 
51 class StructDispNode;
52 
53 class StructNode;
54 
55 /* StructDispNode - begin */
56 
57 /* Nodo strutturale: possiede i gradi di liberta' di:
58  *  - spostamento assoluto,
59  * inoltre, se dinamico, i gdl di:
60  *  - quantita' di moto,
61  * Il nodo di per se' non ha caratteristiche inerziali, che gli vengono date
62  * dagli elementi ad esso collegati. In particolare, gli elementi Body, corpo
63  * rigido, sono responsabili dell'attribuzione di inerzia ai nodi. Altri
64  * contributi possono giungere da travi con matrice di inerzia consistente
65  * (non ancora implementate). */
66 
67 class StructDispNode : public Node, public RigidBodyKinematics {
68 public:
69 	class ErrGeneric : public MBDynErrBase {
70   	public:
ErrGeneric(MBDYN_EXCEPT_ARGS_DECL)71  		ErrGeneric(MBDYN_EXCEPT_ARGS_DECL) : MBDynErrBase(MBDYN_EXCEPT_ARGS_PASSTHRU) {};
72 	};
73 
74 	enum Type {
75 		UNKNOWN = -1,
76 
77 		DYNAMIC = 0,
78 		STATIC,
79 
80 		LASTSTRUCTDISPNODETYPE
81 	};
82 
83 	enum Output {
84 		OUTPUT_ACCELERATIONS = (ToBeOutput::OUTPUT_PRIVATE << 0),
85 		OUTPUT_INERTIA = (ToBeOutput::OUTPUT_PRIVATE << 1)
86 	};
87 
88 protected:
89 	mutable Vec3 XPrev;   /* Posizione al passo precedente */
90 	mutable Vec3 XCurr;   /* Posizione corrente */
91 
92 	mutable Vec3 VPrev;   /* Velocita' al passo precedente */
93 	mutable Vec3 VCurr;   /* Velocita' corrente */
94 
95 	mutable Vec3 XPPCurr;   /* Accelerazione lineare  corrente */
96 	mutable Vec3 XPPPrev;   /* Accelerazione lineare  al passo prec. */
97 
98 	const StructNode *pRefNode;	/* Reference node for relative prediction
99 					WARNING: used only if the relative macro is
100 					active (not default, see configuration options!) */
101 
102 #ifdef USE_NETCDF
103 	NcVar	*Var_X,
104 		*Var_Phi,
105 		*Var_XP,
106 		*Var_Omega,
107 		*Var_XPP,
108 		*Var_OmegaP;
109 #endif /* USE_NETCDF */
110 
111 	OrientationDescription od;
112 
113 	/* Rigidezze fittizie usate nell'assemblaggio dei vincoli */
114 	doublereal dPositionStiffness;
115 	doublereal dVelocityStiffness;
116 
117 	// FIXME: reference motion, for relative kinematics
118 	const RigidBodyKinematics *pRefRBK;
119 
120 	// makes sense also for dummy nodes, as they may inherit
121 	// accelerations from the parent node
122 	bool bOutputAccels;
123 
124 #ifdef USE_AUTODIFF
125 	/*
126 	 * Returns the dof index -1 of VCurr during initial assembly
127 	 */
128 	virtual integer iGetInitialFirstIndexPrime() const=0;
129 #endif
130 
131 public:
132 	/* Costruttore definitivo */
133 	StructDispNode(unsigned int uL,
134 		const DofOwner* pDO,
135 		const Vec3& X0,
136 		const Vec3& V0,
137 		const StructNode *pRN,
138 		const RigidBodyKinematics *pRBK,
139 		doublereal dPosStiff,
140 		doublereal dVelStiff,
141 		OrientationDescription od,
142 		flag fOut);
143 
144 	/* Distruttore (per ora e' banale) */
145 	virtual ~StructDispNode(void);
146 
147 	/* Tipo di nodo */
148 	virtual Node::Type GetNodeType(void) const;
149 
150 	/* FIXME: rigid-body kinematics */
151 	const RigidBodyKinematics *pGetRBK(void) const;
152 
153 	// RBK
154 	const Vec3& GetX(void) const;
155 	const Mat3x3& GetR(void) const;
156 	const Vec3& GetV(void) const;
157 	const Vec3& GetW(void) const;
158 	const Vec3& GetXPP(void) const;
159 	const Vec3& GetWP(void) const;
160 
161 	/* Ritorna il primo indice (-1) di posizione */
162 	virtual inline integer iGetFirstPositionIndex(void) const;
163 
164 	/* Ritorna il primo indice (-1) di Quantita' di moto */
165 	virtual integer iGetFirstMomentumIndex(void) const = 0;
166 
167 	/* Tipo di nodo strutturale */
168 	virtual StructDispNode::Type GetStructDispNodeType(void) const = 0;
169 
170 	/* Contributo del nodo strutturale al file di restart */
171 	virtual std::ostream& Restart(std::ostream& out) const;
172 
173 	virtual std::ostream& DescribeDof(std::ostream& out,
174 		const char *prefix = "",
175 		bool bInitial = false) const;
176 
177 	virtual void DescribeDof(std::vector<std::string>& desc,
178 		bool bInitial = false,
179 		int i = -1) const;
180 
181 	virtual std::ostream& DescribeEq(std::ostream& out,
182 		const char *prefix = "",
183 		bool bInitial = false) const;
184 
185 	virtual void DescribeEq(std::vector<std::string>& desc,
186 		bool bInitial = false,
187 		int i = -1) const;
188 
189 	/* Restituisce il valore del dof iDof;
190 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
191 	virtual const doublereal& dGetDofValue(int iDof, int iOrder = 0) const;
192 
193 	/* Restituisce il valore del dof iDof al passo precedente;
194 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
195 	virtual const doublereal& dGetDofValuePrev(int iDof, int iOrder = 0) const;
196 
197 	/* Setta il valore del dof iDof a dValue;
198 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
199 	virtual void SetDofValue(const doublereal& dValue,
200 		unsigned int iDof, unsigned int iOrder = 0);
201 
202 	virtual DofOrder::Order GetDofType(unsigned int) const;
203 
204 	virtual inline const Vec3& GetXPrev(void) const;
205 	virtual inline const Vec3& GetXCurr(void) const;
206 
207 	virtual inline const Vec3& GetVPrev(void) const;
208 	virtual inline const Vec3& GetVCurr(void) const;
209 
210 	virtual inline const Vec3& GetXPPPrev(void) const;
211 	virtual inline const Vec3& GetXPPCurr(void) const;
212 
213 #ifdef USE_AUTODIFF
214 	inline void GetXCurr(grad::Vector<doublereal, 3>& X, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
215 
216 	template <grad::index_type N_SIZE>
217 	inline void GetXCurr(grad::Vector<grad::Gradient<N_SIZE>, 3>& X, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
218 
219 	inline void GetVCurr(grad::Vector<doublereal, 3>& V, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
220 
221 	template <grad::index_type N_SIZE>
222 	inline void GetVCurr(grad::Vector<grad::Gradient<N_SIZE>, 3>& V, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
223 #endif
224 
225 	virtual inline const doublereal& dGetPositionStiffness(void) const;
226 	virtual inline const doublereal& dGetVelocityStiffness(void) const;
227 
228 	virtual bool ComputeAccelerations(bool b);
229 	virtual inline bool bComputeAccelerations(void) const;
230 	virtual inline bool bOutputAccelerations(void) const;
231 	virtual void OutputAccelerations(bool bOut);
232 
233 	virtual void OutputPrepare(OutputHandler &OH);
234 
235 	/* Output del nodo strutturale (da mettere a punto) */
236 	virtual void Output(OutputHandler& OH) const;
237 
238 #if 0
239 	/* Output della soluzione perturbata (modi ...) */
240 	virtual void Output(OutputHandler& OH,
241 		const VectorHandler& X, const VectorHandler& XP) const;
242 #endif
243 
244 	/* Aggiorna dati in base alla soluzione */
245 	virtual void Update(const VectorHandler& X,
246 		const VectorHandler& XP);
247 
248 	/* Aggiorna dati durante l'iterazione fittizia iniziale */
249 	virtual void DerivativesUpdate(const VectorHandler& X,
250 		const VectorHandler& XP);
251 
252 	/* Ritorna il numero di dofs usato nell'assemblaggio iniziale */
253 	virtual inline unsigned int iGetInitialNumDof(void) const;
254 
255 	/* Aggiorna dati in base alla soluzione durante l'assemblaggio iniziale */
256 	virtual void InitialUpdate(const VectorHandler& X);
257 
258 	/* Inverse Dynamics: */
259 	/* Do Update on node position, velocity or acceleration
260 	 * depending on iOrder */
261 	void Update(const VectorHandler& X, InverseDynamics::Order iOrder);
262 
263 	/* Funzioni di inizializzazione, ereditate da DofOwnerOwner */
264 	virtual void SetInitialValue(VectorHandler& X);
265 	virtual void SetValue(DataManager *pDM,
266 		VectorHandler& X, VectorHandler& XP,
267 		SimulationEntity::Hints *ph = 0);
268 
269 	/* Elaborazione vettori e dati prima e dopo la predizione
270 	 * per MultiStepIntegrator */
271 	virtual void BeforePredict(VectorHandler& X, VectorHandler& XP,
272 		VectorHandler& XPrev,
273 		VectorHandler& XPPrev) const;
274 	virtual void AfterPredict(VectorHandler& X, VectorHandler& XP);
275 
276 	/* Inverse Dynamics: reset orientation parameters */
277 	virtual void AfterConvergence(const VectorHandler& X,
278 			const VectorHandler& XP,
279 			const VectorHandler& XPP);
280 
281 	/* Metodi per l'estrazione di dati "privati".
282 	 * Si suppone che l'estrattore li sappia interpretare.
283 	 * Come default non ci sono dati privati estraibili */
284 	virtual unsigned int iGetNumPrivData(void) const;
285 
286 	/* Maps a string (possibly with substrings) to a private data;
287 	 * returns a valid index ( > 0 && <= iGetNumPrivData()) or 0
288 	 * in case of unrecognized data; error must be handled by caller */
289 	virtual unsigned int iGetPrivDataIdx(const char *s) const;
290 
291 	/* Returns the current value of a private data
292 	 * with 0 < i <= iGetNumPrivData() */
293 	virtual doublereal dGetPrivData(unsigned int i) const;
294 };
295 
296 /* Ritorna il numero di dofs usato nell'assemblaggio iniziale */
297 inline unsigned int
iGetInitialNumDof(void)298 StructDispNode::iGetInitialNumDof(void) const
299 {
300 	return 6;
301 }
302 
303 inline const Vec3&
GetXPrev(void)304 StructDispNode::GetXPrev(void) const
305 {
306 	return XPrev;
307 }
308 
309 inline const Vec3&
GetXCurr(void)310 StructDispNode::GetXCurr(void) const
311 {
312 	return XCurr;
313 }
314 
315 inline const Vec3&
GetVPrev(void)316 StructDispNode::GetVPrev(void) const
317 {
318 	return VPrev;
319 }
320 
321 inline const Vec3&
GetVCurr(void)322 StructDispNode::GetVCurr(void) const
323 {
324 	return VCurr;
325 }
326 
327 inline const Vec3&
GetXPPPrev(void)328 StructDispNode::GetXPPPrev(void) const
329 {
330 	return XPPPrev;
331 }
332 
333 inline const Vec3&
GetXPPCurr(void)334 StructDispNode::GetXPPCurr(void) const
335 {
336 	return XPPCurr;
337 }
338 
339 #ifdef USE_AUTODIFF
340 inline void
GetXCurr(grad::Vector<doublereal,3> & X,doublereal dCoef,enum grad::FunctionCall func,grad::LocalDofMap * pDofMap)341 StructDispNode::GetXCurr(grad::Vector<doublereal, 3>& X, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const
342 {
343 	X = XCurr;
344 }
345 
346 template <grad::index_type N_SIZE>
347 inline void
GetXCurr(grad::Vector<grad::Gradient<N_SIZE>,3> & X,doublereal dCoef,enum grad::FunctionCall func,grad::LocalDofMap * pDofMap)348 StructDispNode::GetXCurr(grad::Vector<grad::Gradient<N_SIZE>, 3>& X, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const
349 {
350 	using namespace grad;
351 
352 	index_type iFirstDofIndex;
353 
354 	switch (func) {
355 	case INITIAL_ASS_JAC:
356 		GRADIENT_ASSERT(dCoef == 1.);
357 	case INITIAL_DER_JAC:
358 	case REGULAR_JAC:
359 		iFirstDofIndex = iGetFirstIndex();
360 		break;
361 
362 	default:
363 		GRADIENT_ASSERT(false);
364 	}
365 
366 	for (index_type i = 1; i <= 3; ++i) {
367 		Gradient<N_SIZE>& g = X(i);
368 		g.SetValuePreserve(XCurr(i));
369 		g.DerivativeResizeReset(pDofMap,
370 								iFirstDofIndex + 1,
371 								iFirstDofIndex + 4,
372 								MapVectorBase::GLOBAL,
373 								0.);
374 		g.SetDerivativeGlobal(iFirstDofIndex + i, -dCoef);
375 	}
376 }
377 
378 inline void
GetVCurr(grad::Vector<doublereal,3> & V,doublereal dCoef,enum grad::FunctionCall func,grad::LocalDofMap * pDofMap)379 StructDispNode::GetVCurr(grad::Vector<doublereal, 3>& V, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const
380 {
381 	V = VCurr;
382 }
383 
384 template <grad::index_type N_SIZE>
385 inline void
GetVCurr(grad::Vector<grad::Gradient<N_SIZE>,3> & V,doublereal dCoef,enum grad::FunctionCall func,grad::LocalDofMap * pDofMap)386 StructDispNode::GetVCurr(grad::Vector<grad::Gradient<N_SIZE>, 3>& V, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const
387 {
388 	using namespace grad;
389 
390 	index_type iFirstDofIndex;
391 
392 	switch (func) {
393 	case INITIAL_ASS_JAC:
394 		GRADIENT_ASSERT(dCoef == 1.);
395 		iFirstDofIndex = iGetInitialFirstIndexPrime();
396 		break;
397 
398 	case INITIAL_DER_JAC:
399 	case REGULAR_JAC:
400 		iFirstDofIndex = iGetFirstIndex();
401 		break;
402 
403 	default:
404 		GRADIENT_ASSERT(false);
405 	}
406 
407 	for (index_type i = 1; i <= 3; ++i) {
408 		Gradient<N_SIZE>& g = V(i);
409 		g.SetValuePreserve(VCurr(i));
410 		g.DerivativeResizeReset(pDofMap,
411 								iFirstDofIndex + 1,
412 								iFirstDofIndex + 4,
413 								MapVectorBase::GLOBAL,
414 								0.);
415 		g.SetDerivativeGlobal(iFirstDofIndex + i, -1.);
416 	}
417 }
418 #endif
419 
420 inline const doublereal&
dGetPositionStiffness(void)421 StructDispNode::dGetPositionStiffness(void) const
422 {
423 	return dPositionStiffness;
424 }
425 
426 inline const doublereal&
dGetVelocityStiffness(void)427 StructDispNode::dGetVelocityStiffness(void) const
428 {
429 	return dVelocityStiffness;
430 }
431 
432 inline bool
bComputeAccelerations(void)433 StructDispNode::bComputeAccelerations(void) const
434 {
435 	return false;
436 }
437 
438 inline bool
bOutputAccelerations(void)439 StructDispNode::bOutputAccelerations(void) const
440 {
441 	return bOutputAccels;
442 }
443 
444 inline void
OutputAccelerations(bool bOut)445 StructDispNode::OutputAccelerations(bool bOut)
446 {
447 	bOutputAccels = bOut;
448 }
449 
450 /* Ritorna il primo indice (-1) di posizione */
451 inline integer
iGetFirstPositionIndex(void)452 StructDispNode::iGetFirstPositionIndex(void) const
453 {
454 	return DofOwnerOwner::iGetFirstIndex();
455 }
456 
457 /* StructDispNode - end */
458 
459 /* DynamicStructDispNode - begin */
460 
461 /* Nodo strutturale per problemi dinamici: possiede i gradi di liberta' di:
462  *  - spostamento assoluto,
463  *  - quantita' di moto,
464  * Il nodo di per se' non ha caratteristiche inerziali, che gli vengono date
465  * dagli elementi ad esso collegati. In particolare, gli elementi Body, corpo
466  * rigido, sono responsabili dell'attribuzione di inerzia ai nodi. Altri
467  * contributi possono giungere da travi con matrice di inerzia consistente
468  * (non ancora implementate). Fa eccezione il caso di un nodo incastrato.
469  * In questo caso non e' necessario attribuirgli inerzia perche' il vincolo
470  * di incastro, ClampJoint, si occupa di rendere non singolare la matrice
471  * jacobiana. */
472 
473 
474 /* Forward declaration */
475 class AutomaticStructDispElem;
476 
477 class DynamicStructDispNode : virtual public StructDispNode {
478 protected:
479 	/* Acceleration and angular acceleration; DynamicStructNode uses them
480 	 * only for output; ModalNode uses them to store actual unknowns */
481 
482 	// "mutable" because it can be set in iGetPrivDataIdx
483 	mutable bool bComputeAccels;
484 	mutable AutomaticStructDispElem *pAutoStr;
485 
486 #ifdef USE_AUTODIFF
487 	virtual inline integer iGetInitialFirstIndexPrime() const;
488 #endif
489 
490 public:
491 	/* Costruttore definitivo (da mettere a punto) */
492 	/* I dati sono passati a mezzo di reference, quindi i relativi oggetti
493 	 * devono essere creati da chi costruisce il nodo, ovvero la funzione
494 	 * DataManager::ReadStructNode(). Non e' il modo piu' efficiente ma e'
495 	 * comodo e sicuro */
496 	DynamicStructDispNode(unsigned int uL,
497 		const DofOwner* pDO,
498 		const Vec3& X0,
499 		const Vec3& V0,
500 	        const StructNode *pRN,
501 		const RigidBodyKinematics *pRBK,
502 		doublereal dPosStiff,
503 		doublereal dVelStiff,
504 		OrientationDescription od,
505 		flag fOut);
506 
507 	/* Distruttore (per ora e' banale) */
508 	virtual ~DynamicStructDispNode(void);
509 
510 	/* Tipo di nodo strutturale */
511 	virtual StructDispNode::Type GetStructDispNodeType(void) const;
512 
513 	virtual inline void SetAutoStr(const AutomaticStructDispElem *p);
514 
515 	/* rigid-body kinematics */
516 	const Vec3& GetXPP(void) const;
517 
518 	/* Ritorna il numero di dofs (comune a tutto cio' che possiede dof) */
519 	virtual inline unsigned int iGetNumDof(void) const;
520 
521 	virtual std::ostream& DescribeDof(std::ostream& out,
522 		const char *prefix = "",
523 		bool bInitial = false) const;
524 
525 	virtual void DescribeDof(std::vector<std::string>& desc,
526 		bool bInitial = false,
527 		int i = -1) const;
528 
529 	virtual std::ostream& DescribeEq(std::ostream& out,
530 		const char *prefix = "",
531 		bool bInitial = false) const;
532 
533 	virtual void DescribeEq(std::vector<std::string>& desc,
534 		bool bInitial = false,
535 		int i = -1) const;
536 
537 	/* Ritorna il primo indice (-1) di quantita' di moto */
538 	virtual inline integer iGetFirstMomentumIndex(void) const;
539 
540 	/* Usato dalle forze astratte, dai bulk ecc., per assemblare le forze
541 	 * al posto giusto */
542 	virtual integer iGetFirstRowIndex(void) const;
543 
544 	virtual void AddInertia(const doublereal& dm) const;
545 
546    	/* Accesso ai suoi dati */
547 	virtual const Vec3& GetBCurr(void) const;
548 	virtual const Vec3& GetBPCurr(void) const;
549 
550 	virtual void AfterConvergence(const VectorHandler& X,
551 		const VectorHandler& XP);
552 
553 	/* Elaborazione vettori e dati prima e dopo la predizione
554 	 * per MultiStepIntegrator */
555 	virtual void BeforePredict(VectorHandler& X, VectorHandler& XP,
556 		VectorHandler& XPrev,
557 		VectorHandler& XPPrev) const;
558 
559 	/* Restituisce il valore del dof iDof;
560 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
561 	virtual const doublereal& dGetDofValue(int iDof, int iOrder = 0) const;
562 
563 	/* Restituisce il valore del dof iDof al passo precedente;
564 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
565 	virtual const doublereal& dGetDofValuePrev(int iDof, int iOrder = 0) const;
566 
567 	/* Setta il valore del dof iDof a dValue;
568 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
569 	virtual void SetDofValue(const doublereal& dValue,
570 		unsigned int iDof, unsigned int iOrder = 0);
571 
572 	/* Aggiorna dati in base alla soluzione */
573 	virtual void Update(const VectorHandler& X,
574 		const VectorHandler& XP);
575 
576 	virtual inline bool bComputeAccelerations(void) const;
577 	virtual bool ComputeAccelerations(bool b);
578 	virtual void SetOutputFlag(flag f = flag(1));
579 };
580 
581 inline void
SetAutoStr(const AutomaticStructDispElem * p)582 DynamicStructDispNode::SetAutoStr(const AutomaticStructDispElem *p)
583 {
584 	pAutoStr = const_cast<AutomaticStructDispElem *>(p);
585 }
586 
587 inline bool
bComputeAccelerations(void)588 DynamicStructDispNode::bComputeAccelerations(void) const
589 {
590 	return bComputeAccels;
591 }
592 
593 /* Ritorna il numero di dofs (comune a tutto cio' che possiede dof) */
594 inline unsigned int
iGetNumDof(void)595 DynamicStructDispNode::iGetNumDof(void) const
596 {
597 	return 6;
598 }
599 
600 /* Ritorna il primo indice (-1) di quantita' di moto */
601 inline integer
iGetFirstMomentumIndex(void)602 DynamicStructDispNode::iGetFirstMomentumIndex(void) const
603 {
604 	return DofOwnerOwner::iGetFirstIndex() + 3;
605 }
606 
607 #ifdef USE_AUTODIFF
608 inline integer
iGetInitialFirstIndexPrime()609 DynamicStructDispNode::iGetInitialFirstIndexPrime() const
610 {
611 	// FIXME: Is it correct this way?
612 	return iGetFirstIndex() + 3;
613 }
614 #endif
615 
616 /* DynamicStructDispNode - end */
617 
618 
619 /* StaticStructDispNode - begin */
620 
621 /* Nodo strutturale per problemi statici: possiede i gradi di liberta' di:
622  *  - spostamento assoluto,
623  * Il nodo puo' essere usato:
624  * - in problemi statici e quasi statici
625  * - quando e' vincolato da un incastro
626  * - per punti geometrici statici, di cui si intende trascurare la dinamica,
627  *   la cui non-singolarita' sia garantita da elementi elastici
628  *   o da vincoli */
629 
630 /* Numero di dof del tipo di nodo - usato anche dal DofManager (?) */
631 class StaticStructDispNode : virtual public StructDispNode {
632 protected:
633 
634 #ifdef USE_AUTODIFF
635 	virtual inline integer iGetInitialFirstIndexPrime() const;
636 #endif
637 
638 public:
639 	/* Costruttore definitivo */
640 	StaticStructDispNode(unsigned int uL,
641 		const DofOwner* pDO,
642 		const Vec3& X0,
643 		const Vec3& V0,
644 	        const StructNode *pRN,
645 		const RigidBodyKinematics *pRBK,
646 		doublereal dPosStiff,
647 		doublereal dVelStiff,
648 		OrientationDescription od,
649 		flag fOut);
650 
651 	/* Distruttore (per ora e' banale) */
652 	virtual ~StaticStructDispNode(void);
653 
654 	/* Tipo di nodo strutturale */
655 	virtual StructDispNode::Type GetStructDispNodeType(void) const;
656 
657 	/* Ritorna il numero di dofs (comune a tutto cio' che possiede dof) */
658 	virtual inline unsigned int iGetNumDof(void) const;
659 
660 	/* Ritorna il primo indice (-1) di quantita' di moto */
661 	virtual inline integer iGetFirstMomentumIndex(void) const;
662 };
663 
664 /* Ritorna il numero di dofs (comune a tutto cio' che possiede dof) */
665 inline unsigned int
iGetNumDof(void)666 StaticStructDispNode::iGetNumDof(void) const
667 {
668 	return 3;
669 }
670 
671 
672 /* Ritorna il primo indice (-1) di quantita' di moto */
673 inline integer
iGetFirstMomentumIndex(void)674 StaticStructDispNode::iGetFirstMomentumIndex(void) const
675 {
676 	return DofOwnerOwner::iGetFirstIndex();
677 }
678 
679 #ifdef USE_AUTODIFF
680 inline integer
iGetInitialFirstIndexPrime()681 StaticStructDispNode::iGetInitialFirstIndexPrime() const
682 {
683 	// FIXME: Is it correct this way?
684 	return iGetFirstIndex() + 3;
685 }
686 #endif
687 
688 /* StaticStructDispNode - end */
689 
690 
691 /* StructNode - begin */
692 
693 /* Nodo strutturale: possiede i gradi di liberta' di:
694  *  - spostamento assoluto,
695  *  - parametri di rotazione incrementali,
696  * inoltre, se dinamico, i gdl di:
697  *  - quantita' di moto,
698  *  - momento della quantita' di moto rispetto al polo mobile.
699  * Il nodo di per se' non ha caratteristiche inerziali, che gli vengono date
700  * dagli elementi ad esso collegati. In particolare, gli elementi Body, corpo
701  * rigido, sono responsabili dell'attribuzione di inerzia ai nodi. Altri
702  * contributi possono giungere da travi con matrice di inerzia consistente
703  * (non ancora implementate). */
704 
705 class StructNode : virtual public StructDispNode
706 #ifdef USE_AUTODIFF
707 , public grad::AlignedAlloc
708 #endif
709 {
710 public:
711 	class ErrGeneric : public MBDynErrBase {
712   	public:
ErrGeneric(MBDYN_EXCEPT_ARGS_DECL)713  		ErrGeneric(MBDYN_EXCEPT_ARGS_DECL) : MBDynErrBase(MBDYN_EXCEPT_ARGS_PASSTHRU) {};
714 	};
715 
716 	enum Type {
717 		UNKNOWN = -1,
718 
719 		DYNAMIC = 0,
720 		STATIC,
721 		MODAL,
722 		DUMMY,
723 
724 		LASTSTRUCTNODETYPE
725 	};
726 
727 protected:
728 	mutable Mat3x3 RPrev;   /* Matrice di rotazione da zero al passo prec. */
729 	Mat3x3 RRef;            /* Matrice di rotazione predetta al passo corr. */
730 	mutable Mat3x3 RCurr;   /* Matrice di rotazione all'iterazione corrente */
731 
732 	mutable Vec3 gRef;
733 	mutable Vec3 gCurr;     /* parametri e derivate correnti */
734 	mutable Vec3 gPRef;
735 	mutable Vec3 gPCurr;
736 
737 	/* Valgono le relazioni:
738 	 *        RCurr = RDelta*RRef                (1)
739 	 *        RDelta = RCurr*RRef^T              (2)
740 	 * In base a questo, dal momento che la matrice RDelta e' richiesta
741 	 * solo in fase di aggiornamento ed e' usata dal nodo stesso, conviene
742 	 * non conservarla e calcolarla in base alla relazione (2).
743 	 */
744 
745 	mutable Vec3 WPrev;   /* Velocita' angolare al passo precedente */
746 	Vec3 WRef;            /* Velocita' angolare predetta al passo corrente */
747 	mutable Vec3 WCurr;   /* Velocita' angolare corrente */
748 
749 	mutable Vec3 WPCurr;    /* Accelerazione angolare corrente */
750 	mutable Vec3 WPPrev;    /* Accelerazione angolare al passo prec. */
751 
752 	// FIXME: qui o in StructDispNode
753 	// const StructNode *pRefNode;	/* Reference node for relative prediction */
754 
755 	/* Rigidezze fittizie usate nell'assemblaggio dei vincoli */
756 	bool bOmegaRot;       /* Flag di velocita' angolare solidale col nodo */
757 
758 	// reference motion, for relative kinematics
759 	// FIXME: qui o in StructDispNode
760 	// const RigidBodyKinematics *pRefRBK;
761 
762 	// makes sense also for dummy nodes, as they may inherit
763 	// accelerations from the parent node
764 	// FIXME: qui o in StructDispNode
765 	// bool bOutputAccels;
766 
767 #ifdef USE_AUTODIFF
768 #ifdef USE_MULTITHREAD
769 	mutable InUse gradInUse;
770 #endif
771 	mutable bool bUpdateRotation;
772 	mutable doublereal dCoefGrad;
773 
774 	static const grad::index_type iNumADVars = 3; // Account for the initial assembly phase
775 	mutable grad::Matrix<grad::Gradient<iNumADVars>, 3, 3> RCurr_grad;
776 	mutable grad::Vector<grad::Gradient<iNumADVars>, 3> WCurr_grad;
777 
778 	template <typename T>
779 	inline void UpdateRotation(const Mat3x3& RRef, const Vec3& WRef, const grad::Vector<T, 3>& g, const grad::Vector<T, 3>& gP, grad::Matrix<T, 3, 3>& RCurr, grad::Vector<T, 3>& WCurr, enum grad::FunctionCall func) const;
780 
781 	void UpdateRotation(doublereal dCoef, enum grad::FunctionCall func) const;
782 
783 	inline void GetgCurr(grad::Vector<grad::Gradient<iNumADVars>, 3>& g, doublereal dCoef, enum grad::FunctionCall func) const;
784 
785 	inline void GetgPCurr(grad::Vector<grad::Gradient<iNumADVars>, 3>& gP, doublereal dCoef, enum grad::FunctionCall func) const;
786 #endif
787 
788 public:
789 	/* Costruttore definitivo */
790 	StructNode(unsigned int uL,
791 		const DofOwner* pDO,
792 		const Vec3& X0,
793 		const Mat3x3& R0,
794 		const Vec3& V0,
795 		const Vec3& W0,
796 		const StructNode *pRN,
797 		const RigidBodyKinematics *pRBK,
798 		doublereal dPosStiff,
799 		doublereal dVelStiff,
800 		bool bOmRot,
801 		OrientationDescription ood,
802 		flag fOut);
803 
804 	/* Distruttore (per ora e' banale) */
805 	virtual ~StructNode(void);
806 
807 	// RBK
808 	const Mat3x3& GetR(void) const;
809 	const Vec3& GetW(void) const;
810 	const Vec3& GetWP(void) const;
811 
812 	/* Contributo del nodo strutturale al file di restart */
813 	virtual std::ostream& Restart(std::ostream& out) const;
814 
815 	virtual std::ostream& DescribeDof(std::ostream& out,
816 		const char *prefix = "",
817 		bool bInitial = false) const;
818 
819 	virtual void DescribeDof(std::vector<std::string>& desc,
820 		bool bInitial = false,
821 		int i = -1) const;
822 
823 	virtual std::ostream& DescribeEq(std::ostream& out,
824 		const char *prefix = "",
825 		bool bInitial = false) const;
826 
827 	virtual void DescribeEq(std::vector<std::string>& desc,
828 		bool bInitial = false,
829 		int i = -1) const;
830 
831 	/* Tipo di nodo strutturale */
832 	virtual StructNode::Type GetStructNodeType(void) const = 0;
833 
834 	/* Restituisce il valore del dof iDof;
835 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
836 	virtual const doublereal& dGetDofValue(int iDof, int iOrder = 0) const;
837 
838 	/* Restituisce il valore del dof iDof al passo precedente;
839 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
840 	virtual const doublereal& dGetDofValuePrev(int iDof, int iOrder = 0) const;
841 
842 	/* Setta il valore del dof iDof a dValue;
843 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
844 	virtual void SetDofValue(const doublereal& dValue,
845 		unsigned int iDof, unsigned int iOrder = 0);
846 
847 	/* Ritorna il numero di dofs usato nell'assemblaggio iniziale */
848 	virtual inline unsigned int iGetInitialNumDof(void) const;
849 
850 	/* Restituiscono i dati correnti */
851 	/* Attenzione: restituiscono un reference ai dati veri, per limitare
852 	 * l'overhead. Tuttavia, una loro modifica e' permanente. Valutare
853 	 * quindi la possibilita' di far passare un const Mat3x3& ecc., in modo
854 	 * da obbligare il chiamante a farsi una copia dei dati */
855 	virtual inline const Vec3& GetgRef(void) const;
856 	virtual inline const Vec3& GetgCurr(void) const;
857 
858 	virtual inline const Vec3& GetgPRef(void) const;
859 	virtual inline const Vec3& GetgPCurr(void) const;
860 
861 	virtual inline const Mat3x3& GetRPrev(void) const;
862 	virtual inline const Mat3x3& GetRRef(void) const;
863 	virtual inline const Mat3x3& GetRCurr(void) const;
864 
865 	virtual inline const Vec3& GetWPrev(void) const;
866 	virtual inline const Vec3& GetWRef(void) const;
867 	virtual inline const Vec3& GetWCurr(void) const;
868 
869 	virtual inline const Vec3& GetWPCurr(void) const;
870 	virtual inline const Vec3& GetWPPrev(void) const;
871 
872 #ifdef USE_AUTODIFF
873 	inline void GetgCurr(grad::Vector<doublereal, 3>& g, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
874 
875 	template <grad::index_type N_SIZE>
876 	inline void GetgCurr(grad::Vector<grad::Gradient<N_SIZE>, 3>& g, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
877 
878                 inline void GetgPCurr(grad::Vector<doublereal, 3>& gP, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
879 
880 	template <grad::index_type N_SIZE>
881 	inline void GetgPCurr(grad::Vector<grad::Gradient<N_SIZE>, 3>& gP, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
882 
883 	inline void GetRCurr(grad::Matrix<doublereal, 3, 3>& R, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
884 
885 	template <grad::index_type N_SIZE>
886 	inline void GetRCurr(grad::Matrix<grad::Gradient<N_SIZE>, 3, 3>& R, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
887 
888 	inline void GetWCurr(grad::Vector<doublereal, 3>& W, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
889 
890 	template <grad::index_type N_SIZE>
891 	inline void GetWCurr(grad::Vector<grad::Gradient<N_SIZE>, 3>& W, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
892 #endif
893 
894 	virtual inline bool bOmegaRotates(void) const;
895 
896 	virtual void OutputPrepare(OutputHandler &OH);
897 
898 	/* Output del nodo strutturale (da mettere a punto) */
899 	virtual void Output(OutputHandler& OH) const;
900 
901 #if 0
902 	/* Output della soluzione perturbata (modi ...) */
903 	virtual void Output(OutputHandler& OH,
904 		const VectorHandler& X, const VectorHandler& XP) const;
905 #endif
906 
907 	/* Aggiorna dati in base alla soluzione */
908 	virtual void Update(const VectorHandler& X,
909 		const VectorHandler& XP);
910 
911 	/* Aggiorna dati durante l'iterazione fittizia iniziale */
912 	virtual void DerivativesUpdate(const VectorHandler& X,
913 		const VectorHandler& XP);
914 
915 	/* Aggiorna dati in base alla soluzione durante l'assemblaggio iniziale */
916 	virtual void InitialUpdate(const VectorHandler& X);
917 
918 	/* Inverse Dynamics: */
919 	/* Do Update on node position, velocity or acceleration
920 	 * depending on iOrder */
921 	void Update(const VectorHandler& X, InverseDynamics::Order iOrder);
922 
923 	/* Funzioni di inizializzazione, ereditate da DofOwnerOwner */
924 	virtual void SetInitialValue(VectorHandler& X);
925 	virtual void SetValue(DataManager *pDM,
926 		VectorHandler& X, VectorHandler& XP,
927 		SimulationEntity::Hints *ph = 0);
928 
929 	/* Elaborazione vettori e dati prima e dopo la predizione
930 	 * per MultiStepIntegrator */
931 	virtual void BeforePredict(VectorHandler& X, VectorHandler& XP,
932 		VectorHandler& XPrev,
933 		VectorHandler& XPPrev) const;
934 	virtual void AfterPredict(VectorHandler& X, VectorHandler& XP);
935 
936 	/*Inverse Dynamics: reset orientation parameters*/
937 	virtual void AfterConvergence(const VectorHandler& X,
938 			const VectorHandler& XP,
939 			const VectorHandler& XPP);
940 
941 	/*
942 	 * Metodi per l'estrazione di dati "privati".
943 	 * Si suppone che l'estrattore li sappia interpretare.
944 	 * Come default non ci sono dati privati estraibili
945 	 */
946 	virtual unsigned int iGetNumPrivData(void) const;
947 
948 	/*
949 	 * Maps a string (possibly with substrings) to a private data;
950 	 * returns a valid index ( > 0 && <= iGetNumPrivData()) or 0
951 	 * in case of unrecognized data; error must be handled by caller
952 	 */
953 	virtual unsigned int iGetPrivDataIdx(const char *s) const;
954 
955 	/*
956 	 * Returns the current value of a private data
957 	 * with 0 < i <= iGetNumPrivData()
958 	 */
959 	virtual doublereal dGetPrivData(unsigned int i) const;
960 }; /* End class StructNode */
961 
962 /* Ritorna il numero di dofs usato nell'assemblaggio iniziale */
963 inline unsigned int
iGetInitialNumDof(void)964 StructNode::iGetInitialNumDof(void) const
965 {
966 	return 12;
967 }
968 
969 
970 /* Restituiscono i dati correnti */
971 /* Attenzione: restituiscono un reference ai dati veri, per limitare
972  * l'overhead. Tuttavia, una loro modifica e' permanente. Valutare
973  * quindi la possibilita' di far passare un const Mat3x3& ecc., in modo
974  * da obbligare il chiamante a farsi una copia dei dati */
975 inline const Vec3&
GetgRef(void)976 StructNode::GetgRef(void) const
977 {
978 	return gRef;
979 }
980 
981 inline const Vec3&
GetgCurr(void)982 StructNode::GetgCurr(void) const
983 {
984 	return gCurr;
985 }
986 
987 inline const Vec3&
GetgPRef(void)988 StructNode::GetgPRef(void) const
989 {
990 	return gPRef;
991 }
992 
993 inline const Vec3&
GetgPCurr(void)994 StructNode::GetgPCurr(void) const
995 {
996 	return gPCurr;
997 }
998 
999 inline const Mat3x3&
GetRPrev(void)1000 StructNode::GetRPrev(void) const
1001 {
1002 	return RPrev;
1003 }
1004 
1005 inline const Mat3x3&
GetRRef(void)1006 StructNode::GetRRef(void) const
1007 {
1008 	return RRef;
1009 }
1010 
1011 inline const Mat3x3&
GetRCurr(void)1012 StructNode::GetRCurr(void) const
1013 {
1014 	return RCurr;
1015 }
1016 
1017 inline const Vec3&
GetWPrev(void)1018 StructNode::GetWPrev(void) const
1019 {
1020 	return WPrev;
1021 }
1022 
1023 inline const Vec3&
GetWRef(void)1024 StructNode::GetWRef(void) const
1025 {
1026 	return WRef;
1027 }
1028 
1029 inline const Vec3&
GetWCurr(void)1030 StructNode::GetWCurr(void) const
1031 {
1032 	return WCurr;
1033 }
1034 
1035 inline const Vec3&
GetWPPrev(void)1036 StructNode::GetWPPrev(void) const
1037 {
1038 	return WPPrev;
1039 }
1040 
1041 inline const Vec3&
GetWPCurr(void)1042 StructNode::GetWPCurr(void) const
1043 {
1044 	return WPCurr;
1045 }
1046 
1047 #ifdef USE_AUTODIFF
GetgCurr(grad::Vector<doublereal,3> & g,doublereal dCoef,enum grad::FunctionCall func,grad::LocalDofMap * pDofMap)1048 inline void StructNode::GetgCurr(grad::Vector<doublereal, 3>& g, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const
1049 {
1050 	g = gCurr;
1051 }
1052 
1053 template <grad::index_type N_SIZE>
GetgCurr(grad::Vector<grad::Gradient<N_SIZE>,3> & g,doublereal dCoef,enum grad::FunctionCall func,grad::LocalDofMap * pDofMap)1054 inline void StructNode::GetgCurr(grad::Vector<grad::Gradient<N_SIZE>, 3>& g, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const
1055 {
1056 	using namespace grad;
1057 
1058 	index_type iFirstDofIndex;
1059 
1060 	switch (func) {
1061 	case INITIAL_ASS_JAC:
1062 		GRADIENT_ASSERT(dCoef == 1.);
1063 
1064 	case INITIAL_DER_JAC:
1065 	case REGULAR_JAC:
1066 		iFirstDofIndex = iGetFirstIndex();
1067 		break;
1068 
1069 	default:
1070 		GRADIENT_ASSERT(false);
1071 	}
1072 
1073 	for (index_type i = 1; i <= 3; ++i) {
1074 		Gradient<N_SIZE>& g_i = g(i);
1075 		g_i.SetValuePreserve(gCurr(i));
1076 		g_i.DerivativeResizeReset(pDofMap,
1077 								  iFirstDofIndex + 4,
1078 								  iFirstDofIndex + 7,
1079 								  MapVectorBase::GLOBAL,
1080 								  0.);
1081 		g_i.SetDerivativeGlobal(iFirstDofIndex + i + 3, -dCoef);
1082 	}
1083 }
1084 
GetgPCurr(grad::Vector<doublereal,3> & gP,doublereal,enum grad::FunctionCall func,grad::LocalDofMap *)1085 inline void StructNode::GetgPCurr(grad::Vector<doublereal, 3>& gP, doublereal, enum grad::FunctionCall func, grad::LocalDofMap*) const
1086 {
1087         GRADIENT_ASSERT(grad::INITIAL_DER_RES || func == grad::REGULAR_RES);
1088 
1089         gP = gPCurr;
1090 }
1091 
1092 template <grad::index_type N_SIZE>
GetgPCurr(grad::Vector<grad::Gradient<N_SIZE>,3> & gP,doublereal,enum grad::FunctionCall func,grad::LocalDofMap * pDofMap)1093 inline void StructNode::GetgPCurr(grad::Vector<grad::Gradient<N_SIZE>, 3>& gP, doublereal, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const
1094 {
1095 	using namespace grad;
1096 
1097 	index_type iFirstDofIndex;
1098 
1099 	switch (func) {
1100 	case INITIAL_ASS_JAC:
1101 	case INITIAL_DER_JAC:
1102 	case REGULAR_JAC:
1103 		iFirstDofIndex = iGetFirstIndex();
1104 		break;
1105 
1106 	default:
1107 		GRADIENT_ASSERT(false);
1108 	}
1109 
1110 	for (index_type i = 1; i <= 3; ++i) {
1111 		Gradient<N_SIZE>& g_i = gP(i);
1112 		g_i.SetValuePreserve(gPCurr(i));
1113 		g_i.DerivativeResizeReset(pDofMap,
1114                                           iFirstDofIndex + 4,
1115                                           iFirstDofIndex + 7,
1116                                           MapVectorBase::GLOBAL,
1117                                           0.);
1118 		g_i.SetDerivativeGlobal(iFirstDofIndex + i + 3, -1.);
1119 	}
1120 }
1121 
GetRCurr(grad::Matrix<doublereal,3,3> & R,doublereal dCoef,enum grad::FunctionCall func,grad::LocalDofMap * pDofMap)1122 inline void StructNode::GetRCurr(grad::Matrix<doublereal, 3, 3>& R, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const
1123 {
1124 	R = RCurr;
1125 }
1126 
1127 template <grad::index_type N_SIZE>
GetRCurr(grad::Matrix<grad::Gradient<N_SIZE>,3,3> & R,doublereal dCoef,enum grad::FunctionCall func,grad::LocalDofMap * pDofMap)1128 inline void StructNode::GetRCurr(grad::Matrix<grad::Gradient<N_SIZE>, 3, 3>& R, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const
1129 {
1130 	using namespace grad;
1131 
1132 	UpdateRotation(dCoef, func);
1133 
1134 	index_type iFirstDofIndex;
1135 
1136 	switch (func) {
1137 	case INITIAL_ASS_JAC:
1138 		GRADIENT_ASSERT(dCoef == 1.);
1139 	case INITIAL_DER_JAC:
1140 	case REGULAR_JAC:
1141 		iFirstDofIndex = iGetFirstIndex() + 4;
1142 		break;
1143 
1144 	default:
1145 		GRADIENT_ASSERT(false);
1146 	}
1147 
1148     for (index_type i = 1; i <= 3; ++i) {
1149         for (index_type j = 1; j <= 3; ++j) {
1150         	const Gradient<iNumADVars>& RCurr_ij = RCurr_grad(i, j);
1151         	Gradient<N_SIZE>& R_ij = R(i, j);
1152 
1153             R_ij.SetValuePreserve(RCurr_ij.dGetValue());
1154             R_ij.DerivativeResizeReset(pDofMap,
1155             						   iFirstDofIndex + RCurr_ij.iGetStartIndexLocal(),
1156             						   iFirstDofIndex + RCurr_ij.iGetEndIndexLocal(),
1157             						   MapVectorBase::GLOBAL,
1158             						   0.);
1159 
1160             for (index_type k = RCurr_ij.iGetStartIndexLocal(); k < RCurr_ij.iGetEndIndexLocal(); ++k) {
1161             	R_ij.SetDerivativeGlobal(iFirstDofIndex + k, RCurr_ij.dGetDerivativeLocal(k));
1162             }
1163         }
1164     }
1165 }
1166 
GetWCurr(grad::Vector<doublereal,3> & W,doublereal dCoef,enum grad::FunctionCall func,grad::LocalDofMap * pDofMap)1167 inline void StructNode::GetWCurr(grad::Vector<doublereal, 3>& W, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const
1168 {
1169 	W = WCurr;
1170 }
1171 
1172 template <grad::index_type N_SIZE>
GetWCurr(grad::Vector<grad::Gradient<N_SIZE>,3> & W,doublereal dCoef,enum grad::FunctionCall func,grad::LocalDofMap * pDofMap)1173 inline void StructNode::GetWCurr(grad::Vector<grad::Gradient<N_SIZE>, 3>& W, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const
1174 {
1175 	using namespace grad;
1176 
1177 	UpdateRotation(dCoef, func);
1178 
1179 	index_type iFirstDofIndex;
1180 
1181 	switch (func) {
1182 	case INITIAL_ASS_JAC:
1183 		GRADIENT_ASSERT(dCoef == 1.);
1184 		/*
1185 		 * XCurr = [X1,		#1
1186 		 * 			X2,		#2
1187 		 * 			X3,		#3
1188 		 * 			g1,		#4
1189 		 * 			g2,		#5
1190 		 * 			g3,		#6
1191 		 * 			V1,		#7
1192 		 * 			V2,		#8
1193 		 * 			V3,		#9
1194 		 * 			W1,		#10
1195 		 * 			W2,		#11
1196 		 * 			W3];	#12
1197 		 */
1198 		iFirstDofIndex = iGetFirstIndex() + 10;
1199 		break;
1200 
1201 	case REGULAR_JAC:
1202 		/*
1203 		 * XCurr = [X1,		#1
1204 		 * 			X2,		#2
1205 		 * 			X3,		#3
1206 		 * 			g1		#4
1207 		 * 			g2		#5
1208 		 * 			g3];	#6
1209 		 */
1210 		iFirstDofIndex = iGetFirstIndex() + 4;
1211 		break;
1212 
1213 	default:
1214 		GRADIENT_ASSERT(false);
1215 	}
1216 
1217     for (int i = 1; i <= 3; ++i) {
1218     	Gradient<N_SIZE>& W_i = W(i);
1219     	const Gradient<iNumADVars>& WCurr_i = WCurr_grad(i);
1220 
1221         W_i.SetValuePreserve(WCurr_i.dGetValue());
1222         W_i.DerivativeResizeReset(pDofMap,
1223         						  iFirstDofIndex + WCurr_i.iGetStartIndexLocal(),
1224         						  iFirstDofIndex + WCurr_i.iGetEndIndexLocal(),
1225         						  MapVectorBase::GLOBAL,
1226         						  0.);
1227 
1228         for (index_type j = WCurr_i.iGetStartIndexLocal(); j < WCurr_i.iGetEndIndexLocal(); ++j) {
1229         	W_i.SetDerivativeGlobal(iFirstDofIndex + j, WCurr_i.dGetDerivativeLocal(j));
1230         }
1231     }
1232 }
1233 #endif
1234 
1235 inline bool
bOmegaRotates(void)1236 StructNode::bOmegaRotates(void) const
1237 {
1238 	return bOmegaRot;
1239 }
1240 
1241 /* StructNode - end */
1242 
1243 
1244 /* DynamicStructNode - begin */
1245 
1246 /* Nodo strutturale per problemi dinamici: possiede i gradi di liberta' di:
1247  *  - spostamento assoluto,
1248  *  - parametri di rotazione incrementali,
1249  *  - quantita' di moto,
1250  *  - momento della quantita' di moto rispetto al polo mobile.
1251  * Il nodo di per se' non ha caratteristiche inerziali, che gli vengono date
1252  * dagli elementi ad esso collegati. In particolare, gli elementi Body, corpo
1253  * rigido, sono responsabili dell'attribuzione di inerzia ai nodi. Altri
1254  * contributi possono giungere da travi con matrice di inerzia consistente
1255  * (non ancora implementate). Fa eccezione il caso di un nodo incastrato.
1256  * In questo caso non e' necessario attribuirgli inerzia perche' il vincolo
1257  * di incastro, ClampJoint, si occupa di rendere non singolare la matrice
1258  * jacobiana. */
1259 
1260 
1261 /* Forward declaration */
1262 class AutomaticStructElem;
1263 
1264 class DynamicStructNode
1265 : virtual public StructDispNode,
1266 public DynamicStructDispNode,
1267 public StructNode
1268 {
1269 protected:
1270 
1271 #ifdef USE_AUTODIFF
1272 	virtual inline integer iGetInitialFirstIndexPrime() const;
1273 #endif
1274 
1275 public:
1276 	/* Costruttore definitivo (da mettere a punto) */
1277 	/* I dati sono passati a mezzo di reference, quindi i relativi oggetti
1278 	 * devono essere creati da chi costruisce il nodo, ovvero la funzione
1279 	 * DataManager::ReadStructNode(). Non e' il modo piu' efficiente ma e'
1280 	 * comodo e sicuro */
1281 	DynamicStructNode(unsigned int uL,
1282 		const DofOwner* pDO,
1283 		const Vec3& X0,
1284 		const Mat3x3& R0,
1285 		const Vec3& V0,
1286 		const Vec3& W0,
1287 	        const StructNode *pRN,
1288 		const RigidBodyKinematics *pRBK,
1289 		doublereal dPosStiff,
1290 		doublereal dVelStiff,
1291 		bool bOmRot,
1292 		OrientationDescription ood,
1293 		flag fOut);
1294 
1295 	/* Distruttore (per ora e' banale) */
1296 	virtual ~DynamicStructNode(void);
1297 
1298 	/* Tipo di nodo strutturale */
1299 	virtual StructNode::Type GetStructNodeType(void) const;
1300 
1301 	/* rigid-body kinematics */
1302 	const Vec3& GetWP(void) const;
1303 
1304 	/* Ritorna il numero di dofs (comune a tutto cio' che possiede dof) */
1305 	virtual inline unsigned int iGetNumDof(void) const;
1306 
1307 	virtual std::ostream& DescribeDof(std::ostream& out,
1308 		const char *prefix = "",
1309 		bool bInitial = false) const;
1310 
1311 	virtual void DescribeDof(std::vector<std::string>& desc,
1312 		bool bInitial = false,
1313 		int i = -1) const;
1314 
1315 	virtual std::ostream& DescribeEq(std::ostream& out,
1316 		const char *prefix = "",
1317 		bool bInitial = false) const;
1318 
1319 	virtual void DescribeEq(std::vector<std::string>& desc,
1320 		bool bInitial = false,
1321 		int i = -1) const;
1322 
1323 	/* Ritorna il primo indice (-1) di quantita' di moto */
1324 	virtual inline integer iGetFirstMomentumIndex(void) const;
1325 
1326 	/* Usato dalle forze astratte, dai bulk ecc., per assemblare le forze
1327 	 * al posto giusto */
1328 	virtual integer iGetFirstRowIndex(void) const;
1329 
1330 	virtual void AddInertia(const doublereal& dm, const Vec3& dS,
1331 		const Mat3x3& dJ) const;
1332 
1333    	/* Accesso ai suoi dati */
1334 	virtual const Vec3& GetGCurr(void) const;
1335 	virtual const Vec3& GetGPCurr(void) const;
1336 
1337 	virtual void AfterConvergence(const VectorHandler& X,
1338 		const VectorHandler& XP);
1339 
1340 	/* Elaborazione vettori e dati prima e dopo la predizione
1341 	 * per MultiStepIntegrator */
1342 	virtual void BeforePredict(VectorHandler& X, VectorHandler& XP,
1343 		VectorHandler& XPrev,
1344 		VectorHandler& XPPrev) const;
1345 
1346 	/* Restituisce il valore del dof iDof;
1347 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
1348 	virtual const doublereal& dGetDofValue(int iDof, int iOrder = 0) const;
1349 
1350 	/* Restituisce il valore del dof iDof al passo precedente;
1351 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
1352 	virtual const doublereal& dGetDofValuePrev(int iDof, int iOrder = 0) const;
1353 
1354 	/* Setta il valore del dof iDof a dValue;
1355 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
1356 	virtual void SetDofValue(const doublereal& dValue,
1357 		unsigned int iDof, unsigned int iOrder = 0);
1358 
1359 	/* Aggiorna dati in base alla soluzione */
1360 	virtual void Update(const VectorHandler& X,
1361 		const VectorHandler& XP);
1362 };
1363 
1364 /* Ritorna il numero di dofs (comune a tutto cio' che possiede dof) */
1365 inline unsigned int
iGetNumDof(void)1366 DynamicStructNode::iGetNumDof(void) const
1367 {
1368 	return 12;
1369 }
1370 
1371 /* Ritorna il primo indice (-1) di quantita' di moto */
1372 inline integer
iGetFirstMomentumIndex(void)1373 DynamicStructNode::iGetFirstMomentumIndex(void) const
1374 {
1375 	return DofOwnerOwner::iGetFirstIndex() + 6;
1376 }
1377 
1378 #ifdef USE_AUTODIFF
1379 inline integer
iGetInitialFirstIndexPrime()1380 DynamicStructNode::iGetInitialFirstIndexPrime() const
1381 {
1382 	return iGetFirstIndex() + 6;
1383 }
1384 #endif
1385 
1386 /* DynamicStructNode - end */
1387 
1388 
1389 /* StaticStructNode - begin */
1390 
1391 /* Nodo strutturale per problemi statici: possiede i gradi di liberta' di:
1392  *  - spostamento assoluto,
1393  *  - parametri di rotazione incrementali
1394  * Il nodo puo' essere usato:
1395  * - in problemi statici e quasi statici
1396  * - quando e' vincolato da un incastro
1397  * - per punti geometrici statici, di cui si intende trascurare la dinamica,
1398  *   la cui non-singolarita' sia garantita da elementi elastici
1399  *   o da vincoli */
1400 
1401 class StaticStructNode
1402 : virtual public StructDispNode,
1403 public StaticStructDispNode,
1404 public StructNode
1405 {
1406 protected:
1407 #ifdef USE_AUTODIFF
1408 	virtual inline integer iGetInitialFirstIndexPrime() const;
1409 #endif
1410 
1411 public:
1412 	/* Costruttore definitivo */
1413 	StaticStructNode(unsigned int uL,
1414 		const DofOwner* pDO,
1415 		const Vec3& X0,
1416 		const Mat3x3& R0,
1417 		const Vec3& V0,
1418 		const Vec3& W0,
1419 	        const StructNode *pRN,
1420 		const RigidBodyKinematics *pRBK,
1421 		doublereal dPosStiff,
1422 		doublereal dVelStiff,
1423 		bool bOmRot,
1424 		OrientationDescription ood,
1425 		flag fOut);
1426 
1427 	/* Distruttore (per ora e' banale) */
1428 	virtual ~StaticStructNode(void);
1429 
1430 	/* Tipo di nodo strutturale */
1431 	virtual StructNode::Type GetStructNodeType(void) const;
1432 
1433 	/* Ritorna il numero di dofs (comune a tutto cio' che possiede dof) */
1434 	virtual inline unsigned int iGetNumDof(void) const;
1435 
1436 	/* Ritorna il primo indice (-1) di quantita' di moto */
1437 	virtual inline integer iGetFirstMomentumIndex(void) const;
1438 };
1439 
1440 /* Ritorna il numero di dofs (comune a tutto cio' che possiede dof) */
1441 inline unsigned int
iGetNumDof(void)1442 StaticStructNode::iGetNumDof(void) const
1443 {
1444 	return 6;
1445 }
1446 
1447 
1448 /* Ritorna il primo indice (-1) di quantita' di moto */
1449 inline integer
iGetFirstMomentumIndex(void)1450 StaticStructNode::iGetFirstMomentumIndex(void) const
1451 {
1452 	return DofOwnerOwner::iGetFirstIndex();
1453 }
1454 
1455 #ifdef USE_AUTODIFF
1456 inline integer
iGetInitialFirstIndexPrime()1457 StaticStructNode::iGetInitialFirstIndexPrime() const
1458 {
1459 	return iGetFirstIndex() + 6;
1460 }
1461 #endif
1462 
1463 /* StaticStructNode - end */
1464 
1465 
1466 /* classe ModalNode derivato da Dynamic */
1467 
1468 class ModalNode : virtual public StructDispNode, public DynamicStructNode {
1469 protected:
1470 #ifdef USE_AUTODIFF
1471 	virtual inline integer iGetInitialFirstIndexPrime() const;
1472 #endif
1473 
1474 public:
1475 	/* Costruttore definitivo (da mettere a punto) */
1476 	/* I dati sono passati a mezzo di reference, quindi i relativi oggetti
1477 	 * devono essere creati da chi costruisce il nodo, ovvero la funzione
1478 	 * DataManager::ReadStructNode(). Non e' il modo piu' efficiente ma e'
1479 	 * comodo e sicuro */
1480 	ModalNode(unsigned int uL,
1481 		const DofOwner* pDO,
1482 		const Vec3& X0,
1483 		const Mat3x3& R0,
1484 		const Vec3& V0,
1485 		const Vec3& W0,
1486 		const RigidBodyKinematics *pRBK,
1487 		doublereal dPosStiff,
1488 		doublereal dVelStiff,
1489 		bool bOmRot,
1490 		OrientationDescription ood,
1491 		flag fOut);
1492 
1493 	/* Distruttore (per ora e' banale) */
1494 	virtual ~ModalNode(void);
1495 
1496 	/* Tipo di nodo strutturale */
1497 	virtual StructNode::Type GetStructNodeType(void) const;
1498 
1499 	virtual std::ostream& DescribeDof(std::ostream& out,
1500 		const char *prefix = "",
1501 		bool bInitial = false) const;
1502 
1503 	virtual void DescribeDof(std::vector<std::string>& desc,
1504 		bool bInitial = false,
1505 		int i = -1) const;
1506 
1507 	virtual std::ostream& DescribeEq(std::ostream& out,
1508 		const char *prefix = "",
1509 		bool bInitial = false) const;
1510 
1511 	virtual void DescribeEq(std::vector<std::string>& desc,
1512 		bool bInitial = false,
1513 		int i = -1) const;
1514 
1515 #if 0
1516 	/* Ritorna il primo indice (-1) di quantita' di moto */
1517 	virtual inline integer iGetFirstMomentumIndex(void) const;
1518 #endif
1519 
1520 	/* Usato dalle forze astratte, dai bulk ecc., per assemblare le forze
1521 	 * al posto giusto */
1522 	virtual integer iGetFirstRowIndex(void) const;
1523 
1524 	/* Aggiorna dati in base alla soluzione */
1525 	virtual void Update(const VectorHandler& X,
1526 		const VectorHandler& XP);
1527 
1528 	virtual void AfterConvergence(const VectorHandler& X,
1529 		const VectorHandler& XP);
1530 
1531 #ifdef USE_AUTODIFF
1532        using StructNode::GetXPPCurr;
1533        using StructNode::GetWPCurr;
1534 
1535        inline void
1536        GetXPPCurr(grad::Vector<doublereal, 3>& XPP, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
1537 
1538        template <grad::index_type N_SIZE>
1539        inline void
1540        GetXPPCurr(grad::Vector<grad::Gradient<N_SIZE>, 3>& XPP, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
1541 
1542        inline void
1543        GetWPCurr(grad::Vector<doublereal, 3>& XPP, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
1544 
1545        template <grad::index_type N_SIZE>
1546        inline void
1547        GetWPCurr(grad::Vector<grad::Gradient<N_SIZE>, 3>& WP, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const;
1548 #endif
1549 };
1550 
1551 
1552 #if 0
1553 /* Ritorna il primo indice (-1) di quantita' di moto */
1554 inline integer
1555 ModalNode::iGetFirstMomentumIndex(void) const
1556 {
1557 	return DofOwnerOwner::iGetFirstIndex() + 6;
1558 }
1559 #endif
1560 
1561 #ifdef USE_AUTODIFF
1562 inline integer
iGetInitialFirstIndexPrime()1563 ModalNode::iGetInitialFirstIndexPrime() const
1564 {
1565 	// FIXME: Don't know how it should be implemented!
1566 	silent_cerr("ModalNode::iGetInitialFirstIndexPrime() not supported yet!" << std::endl);
1567 	throw ErrGeneric(MBDYN_EXCEPT_ARGS);
1568 }
1569 
1570 inline void
GetXPPCurr(grad::Vector<doublereal,3> & XPP,doublereal,enum grad::FunctionCall func,grad::LocalDofMap *)1571 ModalNode::GetXPPCurr(grad::Vector<doublereal, 3>& XPP, doublereal, enum grad::FunctionCall func, grad::LocalDofMap*) const {
1572         GRADIENT_ASSERT(func == grad::INITIAL_DER_RES || func == grad::REGULAR_RES);
1573         XPP = XPPCurr;
1574 }
1575 
1576 template <grad::index_type N_SIZE>
1577 inline void
GetXPPCurr(grad::Vector<grad::Gradient<N_SIZE>,3> & XPP,doublereal dCoef,enum grad::FunctionCall func,grad::LocalDofMap * pDofMap)1578 ModalNode::GetXPPCurr(grad::Vector<grad::Gradient<N_SIZE>, 3>& XPP, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const {
1579 	using namespace grad;
1580 
1581 	index_type iFirstDofIndex;
1582 
1583 	switch (func) {
1584 	case INITIAL_DER_JAC:
1585 	case REGULAR_JAC:
1586 		iFirstDofIndex = iGetFirstIndex();
1587 		break;
1588 
1589 	default:
1590 		GRADIENT_ASSERT(false);
1591 	}
1592 
1593 	for (index_type i = 1; i <= 3; ++i) {
1594 		Gradient<N_SIZE>& g = XPP(i);
1595 		g.SetValuePreserve(XPPCurr(i));
1596 		g.DerivativeResizeReset(pDofMap,
1597                                         iFirstDofIndex + 7,
1598                                         iFirstDofIndex + 10,
1599                                         MapVectorBase::GLOBAL,
1600                                         0.);
1601 		g.SetDerivativeGlobal(iFirstDofIndex + i + 6, -1.);
1602 	}
1603 }
1604 
1605 inline void
GetWPCurr(grad::Vector<doublereal,3> & WP,doublereal,enum grad::FunctionCall func,grad::LocalDofMap *)1606 ModalNode::GetWPCurr(grad::Vector<doublereal, 3>& WP, doublereal, enum grad::FunctionCall func, grad::LocalDofMap*) const {
1607         GRADIENT_ASSERT(func == grad::INITIAL_DER_RES || func == grad::REGULAR_RES);
1608         WP = WPCurr;
1609 }
1610 
1611 template <grad::index_type N_SIZE>
1612 inline void
GetWPCurr(grad::Vector<grad::Gradient<N_SIZE>,3> & WP,doublereal dCoef,enum grad::FunctionCall func,grad::LocalDofMap * pDofMap)1613 ModalNode::GetWPCurr(grad::Vector<grad::Gradient<N_SIZE>, 3>& WP, doublereal dCoef, enum grad::FunctionCall func, grad::LocalDofMap* pDofMap) const {
1614 	using namespace grad;
1615 
1616 	index_type iFirstDofIndex;
1617 
1618 	switch (func) {
1619 	case INITIAL_DER_JAC:
1620 	case REGULAR_JAC:
1621 		iFirstDofIndex = iGetFirstIndex();
1622 		break;
1623 
1624 	default:
1625 		GRADIENT_ASSERT(false);
1626 	}
1627 
1628 	for (index_type i = 1; i <= 3; ++i) {
1629 		Gradient<N_SIZE>& g = WP(i);
1630 		g.SetValuePreserve(WPCurr(i));
1631 		g.DerivativeResizeReset(pDofMap,
1632                                         iFirstDofIndex + 10,
1633                                         iFirstDofIndex + 13,
1634                                         MapVectorBase::GLOBAL,
1635                                         0.);
1636 		g.SetDerivativeGlobal(iFirstDofIndex + i + 9, -1.);
1637 	}
1638 }
1639 
1640 #endif
1641 
1642 /* ModalNode - end */
1643 
1644 
1645 /* DummyStructNode - begin */
1646 
1647 class DummyStructNode : virtual public StructDispNode, public StructNode {
1648 public:
1649 	enum Type {
1650 		UNKNOWN = -1,
1651 
1652 		OFFSET = 0,
1653 		RELATIVEFRAME,
1654 		PIVOTRELATIVEFRAME,
1655 
1656 		LASTTYPE
1657 	};
1658 
1659 protected:
1660 	const StructNode* pNode;
1661 
1662 	virtual void Update_int(void) = 0;
1663 
1664 #ifdef USE_AUTODIFF
1665 	virtual inline integer iGetInitialFirstIndexPrime() const;
1666 #endif
1667 
1668 public:
1669 	/* Costruttore definitivo */
1670 	DummyStructNode(unsigned int uL,
1671 		const DofOwner* pDO,
1672 		const StructNode* pNode,
1673 		OrientationDescription ood,
1674 		flag fOut);
1675 
1676 	/* Distruttore (per ora e' banale) */
1677 	virtual ~DummyStructNode(void);
1678 
1679 	/* tipo */
1680 	virtual DummyStructNode::Type GetDummyType(void) const = 0;
1681 	virtual StructDispNode::Type GetStructDispNodeType(void) const;
1682 
1683 	/* Ritorna il numero di dofs (comune a tutto cio' che possiede dof) */
1684 	virtual inline unsigned int iGetNumDof(void) const;
1685 
1686 	/* Restituisce il valore del dof iDof;
1687 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
1688 	virtual const doublereal& dGetDofValue(int iDof, int iOrder = 0) const;
1689 
1690 	/* Restituisce il valore del dof iDof al passo precedente;
1691 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
1692 	virtual const doublereal& dGetDofValuePrev(int iDof, int iOrder = 0) const;
1693 
1694 	/* Setta il valore del dof iDof a dValue;
1695 	 * se differenziale, iOrder puo' essere = 1 per la derivata */
1696 	virtual void SetDofValue(const doublereal& dValue,
1697 		unsigned int iDof, unsigned int iOrder = 0);
1698 
1699 	/* Tipo di nodo strutturale */
1700 	virtual StructNode::Type GetStructNodeType(void) const;
1701 
1702 	/* Ritorna il numero di dofs usato nell'assemblaggio iniziale */
1703 	virtual inline unsigned int iGetInitialNumDof(void) const;
1704 
1705 	/* Ritorna il primo indice (-1) di posizione */
1706 	virtual inline integer iGetFirstPositionIndex(void) const;
1707 
1708 	/* Ritorna il primo indice (-1) di Quantita' di moto */
1709 	virtual inline integer iGetFirstMomentumIndex(void) const;
1710 
1711 	/* Aggiorna dati durante l'iterazione fittizia iniziale */
1712 	virtual void DerivativesUpdate(const VectorHandler& X,
1713 		const VectorHandler& XP);
1714 
1715 	/* Aggiorna dati in base alla soluzione durante l'assemblaggio iniziale */
1716 	virtual void InitialUpdate(const VectorHandler& X);
1717 
1718 	/* Funzioni di inizializzazione, ereditate da DofOwnerOwner */
1719 	virtual void SetInitialValue(VectorHandler& X);
1720 	virtual void SetValue(DataManager *pDM,
1721 		VectorHandler& X, VectorHandler& XP,
1722 		SimulationEntity::Hints *ph = 0);
1723 
1724 	/* Elaborazione vettori e dati prima e dopo la predizione
1725 	 * per MultiStepIntegrator */
1726 	virtual void BeforePredict(VectorHandler& X, VectorHandler& XP,
1727 		VectorHandler& XPrev,
1728 		VectorHandler& XPPrev) const;
1729 	virtual void AfterPredict(VectorHandler& X, VectorHandler& XP);
1730 
1731 	virtual inline bool bComputeAccelerations(void) const;
1732 	virtual bool ComputeAccelerations(bool b);
1733 };
1734 
1735 /* Ritorna il numero di dofs usato nell'assemblaggio iniziale */
1736 inline unsigned int
iGetInitialNumDof(void)1737 DummyStructNode::iGetInitialNumDof(void) const
1738 {
1739 	return 0;
1740 }
1741 
1742 /* Ritorna il primo indice (-1) di posizione */
1743 inline integer
iGetFirstPositionIndex(void)1744 DummyStructNode::iGetFirstPositionIndex(void) const
1745 {
1746 	silent_cerr("DummyStructNode(" << GetLabel() << ") has no dofs"
1747 		<< std::endl);
1748 	throw ErrGeneric(MBDYN_EXCEPT_ARGS);
1749 }
1750 
1751 /* Ritorna il primo indice (-1) di Quantita' di moto */
1752 inline integer
iGetFirstMomentumIndex(void)1753 DummyStructNode::iGetFirstMomentumIndex(void) const
1754 {
1755 	silent_cerr("DummyStructNode(" << GetLabel() << ") has no dofs"
1756 		<< std::endl);
1757 	throw ErrGeneric(MBDYN_EXCEPT_ARGS);
1758 }
1759 
1760 #ifdef USE_AUTODIFF
1761 inline integer
iGetInitialFirstIndexPrime()1762 DummyStructNode::iGetInitialFirstIndexPrime() const
1763 {
1764 	silent_cerr("DummyStructNode(" << GetLabel() << ") has no dofs"
1765 		<< std::endl);
1766 	throw ErrGeneric(MBDYN_EXCEPT_ARGS);
1767 }
1768 #endif
1769 
1770 /* Ritorna il numero di dofs (comune a tutto cio' che possiede dof) */
1771 inline unsigned int
iGetNumDof(void)1772 DummyStructNode::iGetNumDof(void) const
1773 {
1774 	return 0;
1775 }
1776 
1777 inline bool
bComputeAccelerations(void)1778 DummyStructNode::bComputeAccelerations(void) const
1779 {
1780 	return pNode->bComputeAccelerations();
1781 }
1782 
1783 /* DummyStructNode - end */
1784 
1785 
1786 /* OffsetDummyStructNode - begin */
1787 
1788 class OffsetDummyStructNode : virtual public StructDispNode, public DummyStructNode {
1789 protected:
1790 	Vec3 f;
1791 	Mat3x3 R;
1792 
1793 	void Update_int(void);
1794 
1795 public:
1796 	/* Costruttore definitivo */
1797 	OffsetDummyStructNode(unsigned int uL,
1798 		const DofOwner* pDO,
1799 		const StructNode* pNode,
1800 		const Vec3& f,
1801 		const Mat3x3& R,
1802 		OrientationDescription ood,
1803 		flag fOut);
1804 
1805 	/* Distruttore (per ora e' banale) */
1806 	virtual ~OffsetDummyStructNode(void);
1807 
1808 	/* tipo */
1809 	virtual DummyStructNode::Type GetDummyType(void) const;
1810 
1811 	/* Aggiorna dati in base alla soluzione */
1812 	virtual void Update(const VectorHandler& X,
1813 		const VectorHandler& XP);
1814 };
1815 
1816 /* OffsetDummyStrNode - end */
1817 
1818 
1819 /* RelFrameDummyStructNode - begin */
1820 
1821 class RelFrameDummyStructNode : virtual public StructDispNode, public DummyStructNode {
1822 protected:
1823 	const StructNode* pNodeRef;
1824 	const Mat3x3 RhT;
1825 	const Vec3 fhT;
1826 
1827 	void Update_int(void);
1828 
1829 public:
1830 	/* Costruttore definitivo */
1831 	RelFrameDummyStructNode(unsigned int uL,
1832 		const DofOwner* pDO,
1833 		const StructNode* pNode,
1834 		const StructNode* pNodeRef,
1835 		const Vec3& fh,
1836 		const Mat3x3& Rh,
1837 		OrientationDescription ood,
1838 		flag fOut);
1839 
1840 	/* Distruttore (per ora e' banale) */
1841 	virtual ~RelFrameDummyStructNode(void);
1842 
1843 	/* tipo */
1844 	virtual DummyStructNode::Type GetDummyType(void) const;
1845 
1846 	/* Aggiorna dati in base alla soluzione */
1847 	virtual void Update(const VectorHandler& X,
1848 		const VectorHandler& XP);
1849 
1850 	virtual inline bool bComputeAccelerations(void) const;
1851 	virtual bool ComputeAccelerations(bool b);
1852 };
1853 
1854 inline bool
bComputeAccelerations(void)1855 RelFrameDummyStructNode::bComputeAccelerations(void) const
1856 {
1857 	return pNode->bComputeAccelerations() && pNodeRef->bComputeAccelerations();
1858 }
1859 
1860 /* RelFrameDummyStructNode - end */
1861 
1862 /* PivotRelFrameDummyStructNode - begin */
1863 
1864 class PivotRelFrameDummyStructNode : virtual public StructDispNode, public RelFrameDummyStructNode {
1865 protected:
1866 	const StructNode* pNodeRef2;
1867 	const Mat3x3 Rh2;
1868 	const Vec3 fh2;
1869 
1870 	void Update_int(void);
1871 
1872 public:
1873 	/* Costruttore definitivo */
1874 	PivotRelFrameDummyStructNode(unsigned int uL,
1875 		const DofOwner* pDO,
1876 		const StructNode* pNode,
1877 		const StructNode* pNodeRef,
1878 		const Vec3& fh,
1879 		const Mat3x3& Rh,
1880 		const StructNode* pNodeRef2,
1881 		const Vec3& fh2,
1882 		const Mat3x3& Rh2,
1883 		OrientationDescription ood,
1884 		flag fOut);
1885 
1886 	/* Distruttore (per ora e' banale) */
1887 	virtual ~PivotRelFrameDummyStructNode(void);
1888 
1889 	/* tipo */
1890 	virtual DummyStructNode::Type GetDummyType(void) const;
1891 
1892 	/* Aggiorna dati in base alla soluzione */
1893 	virtual void Update(const VectorHandler& X,
1894 		const VectorHandler& XP);
1895 
1896 	virtual inline bool bComputeAccelerations(void) const;
1897 	virtual bool ComputeAccelerations(bool b);
1898 };
1899 
1900 inline bool
bComputeAccelerations(void)1901 PivotRelFrameDummyStructNode::bComputeAccelerations(void) const
1902 {
1903 	return pNode->bComputeAccelerations()
1904 		&& pNodeRef->bComputeAccelerations()
1905 		&& pNodeRef2->bComputeAccelerations();
1906 }
1907 
1908 /* PivotRelFrameDummyStructNode - end */
1909 
1910 class DataManager;
1911 class MBDynParser;
1912 
1913 extern Node*
1914 ReadStructNode(DataManager* pDM,
1915 	MBDynParser& HP,
1916 	DofOwner* pDO,
1917 	unsigned int uLabel);
1918 
1919 #endif /* STRNODE_H */
1920