1 //#*************************************************************************************** 2 //# 3 //# filename: element.h 4 //# 5 //# author: Gerstmayr Johannes 6 //# 7 //# generated: July 2004 8 //# description: 9 //# 10 //# remarks: 11 //# 12 //# Copyright (c) 2003-2013 Johannes Gerstmayr, Linz Center of Mechatronics GmbH, Austrian 13 //# Center of Competence in Mechatronics GmbH, Institute of Technical Mechanics at the 14 //# Johannes Kepler Universitaet Linz, Austria. All rights reserved. 15 //# 16 //# This file is part of HotInt. 17 //# HotInt is free software: you can redistribute it and/or modify it under the terms of 18 //# the HOTINT license. See folder 'licenses' for more details. 19 //# 20 //# bug reports are welcome!!! 21 //# WWW: www.hotint.org 22 //# email: bug_reports@hotint.org or support@hotint.org 23 //#*************************************************************************************** 24 25 #ifndef ELEMENT__H 26 #define ELEMENT__H 27 28 #include "mbs_interface.h" 29 30 #include "FieldVariableDescriptor.h" 31 #include "mbsload.h" // elements need to deal with loads 32 33 class ReadWriteElementDataVariableType; 34 35 class Constraint; 36 37 38 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 39 // ELEMENT 40 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 41 42 typedef enum {TElement=1, TConstraint=2, TBody=4,/* TFFRF=8, TBody2D=16, TBody3D=32,*/ 43 TCMS = 64, TCMSflag = 128, TController = 256, TFiniteElement = 512, TGCMS = 1024, TParticle = 2048, TIOElementDataModifier = 4096} TMBSElement; 44 // TCMSflag ... set if the Element is a BaseCMSElement 45 // TCMS ....... set if the Element is a finite element belonging to a BaseCMSElement (i.e. CMS *OR* GCMS) 46 // TGCMS ...... set if the Element is a finite element belonging to a GCMSElement 47 // TParticle... set if the Element is a particle (e.g., of a fluid, see the classes SPHParticle2D, SPHParticle3D, SPHMass2D, SPHMass3D) 48 49 //definition of element equations structure: 50 //this information is for element description and for decision of solvers which method to use 51 typedef enum {TET_algebraic=1, TET_algebraic_linear=2, 52 TET_first_order_ODE=4, TET_first_order_ODE_linear=8, 53 TET_second_order_ODE=16, TET_second_order_ODE_linear=32, 54 TET_discontinuous=64, TET_fixedpoint_iteration=128, 55 TET_delay_equation=256, TET_stochastic_process=512, 56 TET_constant_mass_matrix=1024, TET_constant_stiffness_matrix=2048, 57 TET_Lagrange_multipliers=4096 } TEquationType; 58 59 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 60 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 61 //new position, displacement, velocity ... access functions 62 typedef enum{TKAF_none=0, //dummy for initialization 63 TKAF_position=1, //GetPos 64 TKAF_displacement=2, //GetDisplacement 65 TKAF_ref_conf_position=131072, //GetRefConfPos 66 TKAF_ref_position=4194304, //GetRefPos 67 TKAF_node_position = 4, //GetNodePos 68 TKAF_node_ref_conf_position = 262144, //GetNodeRefPos 69 TKAF_global_to_local_position=8, 70 TKAF_velocity=16, //GetVel 71 TKAF_node_velocity=32, //GetNodeVel 72 TKAF_acceleration=64, 73 TKAF_angular_velocity=128, //GetAngularVel 74 TKAF_rotation_matrix=256, //GetRotMatrix 75 TKAF_rotation_matrix_P=512, //GetRotMatrixP 76 TKAF_D_pos_D_q=1024, //GetdPosdqT 77 TKAF_D_pos_D_x=2048, //GetdPosdx 78 TKAF_D_node_pos_D_q = 4096, //GetNodedPosdqT 79 TKAF_D_ang_vel_D_q =8192, //GetdAngVeldqpT 80 TKAF_D_rot_D_q=16384, //GetdRotdqT 81 TKAF_D_rot_v_D_q=32768, //GetdRotvdqT 82 TKAF_int_D_u_D_q =65536, //GetIntDuDq 83 TKAF_position_2D =524288, //GetPos2D 84 TKAF_velocity_2D =1048576, //GetVel2D 85 TKAF_D_pos_D_q_2D =2097152, //GetdPosdqT2D 86 87 88 //GetIntDuDqFCentrifugal 89 TKAF_maximum=8388608} TKinematicsAccessFunctions; //always update maximum, such that it is still the maximum, used in Constraint:IsSuitableElement 90 91 //flag T_compute, T_draw, T_reference_configuration and T_initial_values decides, where to retrieve the system coordinate 92 //flag T_draw_magnified decides, whether the displacment is magnified or not in the drawing function (e.g. GetPos_dc) 93 //flag T_reference_configuration returns zero under the assumption that in reference configuration the coordinates are zero (especially for rigid bodies and classical displacement based finite elements) 94 //flag are not implemented yet 95 typedef enum {TCD_compute=1, TCD_draw=2, TCD_reference_configuration=4, TCD_initial_values=8, TCD_draw_magnified=16, TCD_cached=32} TComputeDrawInitFlag; 96 97 class Element //$EDC$[beginclass,classname=Element] 98 { 99 100 public: Element()101 Element():ltg(), ltgdata(), constraintindices(), loads(), constraints(), constraints_nodouble(), drawelements(), sensors(), elements(), x_init(), data_init() 102 { 103 Element::ElementDefaultConstructorInitialization(); 104 //damping_m=0; 105 //altshape = 0; 106 //type = TElement; 107 //elementname = GetElementSpec(); 108 //rho = 0; 109 //materialnum = 0; 110 //draw_element = 1; 111 }; Element(MBS * mbsi)112 Element(MBS* mbsi):ltg(), ltgdata(), constraintindices(), loads(), constraints(), constraints_nodouble(), drawelements(), sensors(), elements(), x_init(), data_init() 113 { 114 Element::ElementDefaultConstructorInitialization(); 115 mbs=mbsi; 116 //damping_m=0; 117 //altshape = 0; 118 //type = TElement; 119 //elementname = GetElementSpec(); 120 //rho = 0; 121 //materialnum = 0; 122 //draw_element = 1; 123 }; Element(const Element & e)124 Element(const Element& e):ltg(), ltgdata(), constraintindices(), loads(), constraints(), constraints_nodouble(), drawelements(), sensors(), elements(), x_init(), data_init() 125 { 126 CopyFrom(e); 127 }; 128 Element& operator=(const Element& e) 129 { 130 if (this == &e) {return *this;} 131 CopyFrom(e); 132 return *this; 133 } 134 //To be overwritten in derived class: 135 virtual Element* GetCopy(); 136 137 virtual ~Element(); 138 139 //To be overwritten in derived class: 140 virtual void CopyFrom(const Element& e); 141 142 143 //this function assigns default values to the element variables ElementDefaultConstructorInitialization()144 virtual void ElementDefaultConstructorInitialization() 145 { 146 damping_m=0; 147 altshape = 0; 148 type = TElement; 149 elementname = GetElementSpec(); 150 //rho = 0; 151 materialnum = 0; 152 draw_element = 1; 153 mass = 0; 154 elnum = 0; 155 col = Vector3D(0.1,0.1,0.8); //default body color 156 //col = Vector3D(-1,-1,-1); //default body color 157 158 x_init = Vector(2*SOS()+ES()+IS()); //initial conditions 159 data_init = Vector(DataS()); //initial conditions for data variables 160 161 loads.SetLen(0); //$ DR 2012-10 162 163 //not initilized: TArray<int> ltg; //local to global dof reference vector; 164 //not initilized: TArray<int> ltgdata; //local to global data reference vector; 165 //not initilized: TArray<int> constraintindices; //index of this element in the constraint (mostly 1 or 2) 166 //not initilized: TArray<char> dependencies; //contains dependencies to coordinates 167 //not initilized: TArray<Constraint*> constraints; //this constraint list might be double! 168 //not initilized: TArray<Constraint*> constraints_nodouble; //one constraint is only added once to an element 169 //not initilized: TArray<int> sensors; //EDC[varaccess,EDCvarname="sensors",EDCfolder="Info",tooltiptext="attached sensors",readonly]//if the element equations of motion are dependent on sensor values (e.g. as an input), then add the sensor numbers here 170 //not initilized: TArray<int> elements; //dependent elements (mainly for contraints, but also for follower loads!) 171 //not initilized: TArray<int> drawelements; 172 } 173 Initialize()174 virtual void Initialize() //initialize after initial conditions are set!!! 175 { 176 }; 177 GetMBS()178 virtual const MBS* GetMBS() const {return mbs;} GetMBS()179 virtual MBS* GetMBS() {return mbs;} 180 //virtual UserOutput& UO() {return mbs->uout;}; 181 //virtual UserOutput& UO(int message_level = UO_LVL_all) { mbs->uout.SetLocalMessageLevel(message_level); return mbs->uout;}; //(AD) 182 //$ YV 2012-11-28 183 //virtual UserOutput& UO(int message_level = UO_LVL_all, int output_prec = -1) { mbs->uout.SetLocalMessageLevel(message_level); mbs->uout.SetOutputPrec(output_prec); return mbs->uout;}; //$ AD 2011-02 output_prec 184 virtual UserOutputInterface& UO(int message_level = UO_LVL_all, int output_prec = -1) const { return mbs->UO(message_level,output_prec); } MaxIndex()185 virtual int MaxIndex() const {return mbs->MaxIndex();} 186 187 //$ JG 2011-02: compute index 3 constraint drift, only for postprocessing! (evaluate sensors) 188 virtual double GetConstraintDrift(double t) const; 189 GetOwnNum()190 virtual int GetOwnNum() const {return elnum;} SetOwnNum(int i)191 virtual void SetOwnNum(int i) {elnum = i;} 192 GetElementName()193 virtual const mystr& GetElementName() const {return elementname;} GetElementName()194 virtual mystr& GetElementName() {return elementname;} SetElementName(const char * name)195 virtual void SetElementName(const char* name) {elementname = name;} GetElementSpec()196 virtual const char* GetElementSpec() const {return "Element";} //$ DR this is implemented already in many elements and used at multiple places GetElementSpecification()197 virtual mystr GetElementSpecification() //$EDC$[funcaccess,EDCvarname="element_type",tooltiptext="specification of element type. Once the element is added to the mbs, you MUST NOT change this type anymore!"] 198 { 199 mystr type = GetElementSpec(); //$ DR this is necessary for the skript language 200 return type; 201 } 202 IsType(TMBSElement te)203 virtual int IsType(TMBSElement te) const 204 { 205 return (te&type) != 0; 206 } SetType(TMBSElement te)207 virtual void SetType(TMBSElement te) {type = te;} 208 GetType()209 virtual TMBSElement GetType() const {return type;} //return element type (e.g. for knowing if element is a constraint, etc.) 210 //the following string returns a list of strings which contains main specification of the element 211 //virtual MyStrList GetType_String() const 212 //{ 213 // MyStrList strlist; 214 // if (IsType(TElement)) {strlist.Add("connector element");} 215 // if (IsType(TConstraint)) {strlist.Add("connector element");} 216 // if (IsType(TBody) && !IsType(TFiniteElement)) {strlist.Add("body");} 217 // if (IsType(TFiniteElement)) {strlist.Add("finite element");} 218 // if (IsType(TCMS)) {strlist.Add("finite element belonging to a CMS reference frame");} 219 // if (IsType(TCMSflag)) {strlist.Add("reference frame for component mode synthesis (CMS)\n");} 220 // if (IsType(TGCMS)) {strlist.Add("GCMS modally reduced element\n");} 221 // if (IsType(TController)) {strlist.Add("control element (with input and output)\n");} 222 // if (IsType(TParticle)) {strlist.Add("particle for sph simulation\n");} 223 224 // return str; 225 //} 226 AddType(TMBSElement te)227 virtual void AddType(TMBSElement te) 228 { 229 type = (TMBSElement)(type|te); 230 } IsRigid()231 virtual int IsRigid() const {return 0;} //$EDC$[funcaccess,readonly,EDCvarname="is_rigid",EDCfolder="Info",tooltiptext="flag if body is a rigid body"]//default value IsFiniteElement()232 virtual int IsFiniteElement() const {return 0;} //$EDC$[funcaccess,readonly,EDCvarname="is_finite_element",EDCfolder="Info",tooltiptext="flag if body is a finite element"] //if is of type "FiniteElement" 2D or 3D 233 234 ////return the information about the equations type of the element 235 //virtual TElementEquationType GetEquationsType() const 236 //{ 237 // TElementEquationType eet = 0; 238 // if (IS() != 0) eet += TET_algebraic; 239 // if (ES() != 0) eet += TET_first_order_ODE; 240 // if (SOS() != 0) eet += TET_second_order_ODE; 241 // if (IS() != 0) eet += TET_Lagrange_multipliers; 242 243 // return eet; 244 //} 245 //return stringlist containing the possible types of the equation that apply 246 //virtual MyStrList GetEquationsType_String() const 247 //{ 248 // MyStrList strlist; 249 // if ((GetEquationsType()&TET_algebraic != 0) && (GetEquationsType()&TET_algebraic_linear != 0)) {strlist.Add("linear algebraic equations");} 250 // else if ((GetEquationsType()&TET_algebraic != 0)) {strlist.Add("nonlinear algebraic equations");} 251 252 // if ((GetEquationsType()&TET_first_order_ODE != 0) && (GetEquationsType()&TET_first_order_ODE_linear != 0)) {strlist.Add("linear first order differential equations");} 253 // else if ((GetEquationsType()&TET_first_order_ODE != 0)) {strlist.Add("nonlinear first order differential equations");} 254 255 // if ((GetEquationsType()&TET_second_order_ODE != 0) && (GetEquationsType()&TET_second_order_ODE_linear != 0)) {strlist.Add("linear second order differential equations");} 256 // else if ((GetEquationsType()&TET_second_order_ODE != 0)) {strlist.Add("nonlinear second order differential equations");} 257 258 // if (GetEquationsType()&TET_discontinuous != 0) {strlist.Add("discontinuous equations (switches, jumps, etc.)");} 259 // if (GetEquationsType()&TET_fixedpoint_iteration != 0) {strlist.Add("fixed point iteration for discontinuous equations");} 260 261 // if (GetEquationsType()&TET_delay_equation != 0) {strlist.Add("delay equations (storage and retrieving of data)");} 262 // if (GetEquationsType()&TET_stochastic_process != 0) {strlist.Add("stochastic process (solution might be non-deterministic)");} 263 // if (GetEquationsType()&TET_constant_mass_matrix != 0) {strlist.Add("second order equations contain a constant mass matrix");} 264 // if (GetEquationsType()&TET_constant_stiffness_matrix != 0) {strlist.Add("second order equations contain a constant stiffness matrix");} 265 // if (GetEquationsType()&TET_Lagrange_multipliers != 0) {strlist.Add("Lagrange multipliers are used and forces added to other elements");} 266 //} 267 268 #pragma endregion 269 PerformNodeCheck()270 virtual int PerformNodeCheck() const {return 1;} //node number are not checked for CMS-elements ... 271 virtual int CheckConsistency(mystr& errorstr); //rv==0 --> OK, rv==1 --> can not compute, rv==2 --> can not draw and not compute 272 273 //the following functions are for reading/writing the whole element data including the necessary initialization after setting the element data 274 virtual void GetElementData(ElementDataContainer& edc); //fill in all element data 275 virtual int SetElementData(ElementDataContainer& edc); //set element data acc. to ElementDataContainer, return 0 if failed or values invalid! 276 virtual void GetElementDataAuto(ElementDataContainer& edc); //fill in all element data 277 virtual int SetElementDataAuto(ElementDataContainer& edc); //set element data acc. to ElementDataContainer, return 0 if failed or values invalid! 278 279 280 //the following functions are for read/write access to single variables or functions (mapped to EDC-style-variables) 281 //these functions should be implemented in all derived classes 282 virtual int ReadSingleElementData(ReadWriteElementDataVariableType& RWdata); //retrieve value of element variable/function named RWdata.variable_name with specified components (vector/matrix); return 1, if variable found and read; 0 if not found, -1 if no readaccess, -2 if range-fault 283 virtual int WriteSingleElementData(const ReadWriteElementDataVariableType& RWdata); //write value of element variable/function named RWdata.variable_name with specified components (vector/matrix); return 1, if variable found and written; 0 if not found, -1 if no writeaccess, -2 if range-fault 284 virtual int GetAvailableSpecialValues(TArrayDynamic<ReadWriteElementDataVariableType>& available_variables); //$ DR 2012-10 285 286 virtual int ReadSingleElementDataAuto(ReadWriteElementDataVariableType& RWdata); //automatically generated function from EDC_converter 287 virtual int WriteSingleElementDataAuto(const ReadWriteElementDataVariableType& RWdata); //automatically generated function from EDC_converter 288 virtual int GetAvailableSpecialValuesAuto(TArrayDynamic<ReadWriteElementDataVariableType>& available_variables); //automatically generated function from EDC_converter 289 290 //virtual int IsConstraint() const {return 0;} 291 PrecomputeEvalFunctions()292 virtual void PrecomputeEvalFunctions() {}; //PG: precomputation for EvalF(), EvalF2(), EvalM(), and also for EvalG() and AddElementCqTLambda() in case of constraints EvalF(Vector & f,double t)293 virtual void EvalF(Vector& f, double t) {}; //first order equations: \dot q=f(q,t), return vector: len(q) EvalG(Vector & f,double t)294 virtual void EvalG(Vector& f, double t) {}; //evaluate constraints: len(z) 295 //second order equations: M \ddot u = F2(u,\dot u,t), len(u) EvalM(Matrix & m,double t)296 virtual void EvalM(Matrix& m, double t) {}; //evaluate mass matrix 297 //Element-elastic forces+external forces+Cq^T*lambda 298 virtual void EvalF2(Vector& f, double t); 299 virtual void EvalMinvF2(Vector& f, double t); 300 //Jacobian for F2, compute element Jacobian m, LTGVector ref 301 virtual void JacobianF2(double t, Matrix& m, IVector& colref); 302 virtual void JacobianG(double t, Matrix& m, IVector& colref); GetKineticEnergy()303 virtual double GetKineticEnergy() {return 0;}; //compute kinetic energy of system 304 virtual double GetPotentialEnergy(); //compute potential (strain) energy of system 305 FastStiffnessMatrix()306 virtual int FastStiffnessMatrix() const {return 0;} StiffnessMatrix(Matrix & m)307 virtual void StiffnessMatrix(Matrix& m) {assert(0 && "ERROR: StiffnessMatrix() not defined"); }; //fill in sos x sos components, m might be larger 308 GyroscopicMatrix(SparseMatrix & gy)309 virtual void GyroscopicMatrix(SparseMatrix& gy)const // $ MSax 2013-07-25 : added 310 { 311 // gyroscopic matrix for rotordynamic elements. Used e.g. for computation of eigenmodes (campbell) 312 gy = SparseMatrix(); 313 gy.SetSize(SOS(),SOS()); 314 gy.FillWithZeros(); 315 }; 316 317 //to be activated, MBS->SetTransformJacApply must be set: TransformJacApply()318 virtual int TransformJacApply() const {return 0;}; //transform jacobian in solver (x = x_0 + A*J^(-1)*B * f) ApplyTransform(const Vector & v,Vector & Av,int mode)319 virtual void ApplyTransform(const Vector& v, Vector& Av, int mode) {}; //compute Av=A^T*v in mode==0 and Av=A*v in mode==1 320 321 //sparse jacobian operations: AddMSparse(SparseMatrix & m,double t)322 virtual void AddMSparse(SparseMatrix& m, double t) {}; //add sparse mass matrix into full system matrix AddKSparse(SparseMatrix & m,double t)323 virtual void AddKSparse(SparseMatrix& m, double t) {}; //add sparse stiffness matrix into full system matrix AddDSparse(SparseMatrix & m,double t)324 virtual void AddDSparse(SparseMatrix& m, double t) {}; //add sparse damping matrix into full system matrix (at sosmbs) UseSparseM()325 virtual int UseSparseM() const {return 0;}; UseSparseK()326 virtual int UseSparseK() const {return 0;}; 327 virtual void SetUseSparseM(int i=1) {}; 328 virtual void SetUseSparseK(int i=1) {}; 329 330 //floating frame of reference formulation: for FFRF elements GetI1(Vector & I1)331 virtual void GetI1(Vector& I1) {assert(0);}; //Shabana p. 209-211 GetIkl(int k,int l)332 virtual double GetIkl(int k, int l) {assert(0);return 0;}; GetIbarkl(int k,int l,Vector & I1)333 virtual void GetIbarkl(int k, int l, Vector& I1) {assert(0);}; GetSbar(Matrix & Sbar)334 virtual void GetSbar(Matrix& Sbar) {assert(0);}; GetSbarkl(int k,int l,Matrix & Sbar)335 virtual void GetSbarkl(int k, int l, Matrix& Sbar) {assert(0);}; EvalMff(Matrix & m,double t)336 virtual void EvalMff(Matrix& m, double t) {assert(0);}; //matrix without rigid body motion, only flexible part GetH(Matrix & H)337 virtual void GetH(Matrix& H) {assert(0);}; NCMSNodes()338 virtual int NCMSNodes() const {mbs->UO() << "Element::NCMSNodes() called, should be called for FFRFElement only!!\n"; return 0;} 339 340 //dir: normal == 0=X, 1=-X, 2=Y, 3=-Y, 4=Z, 5=-Z AddSurfacePressure(Vector & f,double pressure,int dir)341 virtual void AddSurfacePressure(Vector& f, double pressure, int dir) {assert(0);} GetSurfaceNormalD(int dir)342 virtual Vector3D GetSurfaceNormalD(int dir) { assert(0); return Vector3D(); } // AH: get normal to a specfied surface in order to plot surface pressure 343 IS()344 virtual int IS() const {return 0;}; //$EDC$[funcaccess,readonly,EDCvarname="n_algebraic_equations",EDCfolder="Info",tooltiptext="number of algebraic equations"] //implicit (algebraic) size SOS()345 virtual int SOS() const {return 0;}; //$EDC$[funcaccess,readonly,EDCvarname="n_second_order_ODEs",EDCfolder="Info",tooltiptext="number of second order ODEs in which the element contributes with some terms to the mass matrix and RHS"] //size of second order components (size of M and EvalF2=^=K) SOSowned()346 virtual int SOSowned() const {return SOS();}; //$EDC$[funcaccess,readonly,EDCvarname="n_second_order_ODE_unknowns",EDCfolder="Info", tooltiptext="number of second order ODE unknowns, which are owned (caused) by one element (not borrowed e.g. from another element or node)] //size of second order unknowns, len(u) or len(v) ES()347 virtual int ES() const {return 0;}; //$EDC$[funcaccess,readonly,EDCvarname="n_first_order_ODEs",EDCfolder="Info",tooltiptext="number of first order ODEs"] //size of first order explicit equations SS()348 virtual int SS() const {return 2*SOS()+ES()+IS();}; //$EDC$[funcaccess,readonly,EDCvarname="element_equation_size",EDCfolder="Info",tooltiptext="element system size=2*SOS()+ES()+IS()"]//system size DataS()349 virtual int DataS() const {return 0;} //$EDC$[funcaccess,readonly,EDCvarname="data_size",EDCfolder="Info",tooltiptext="number of data variables (e.g. plastic strains), for which there are no separate equations"] //Data size for non-state variables (contact condition, plastic strains, discrete states, etc.) 350 IS_RS()351 virtual int IS_RS() const {return IS();}; //implicit (algebraic) size SOSowned_RS()352 virtual int SOSowned_RS() const {return SOSowned();}; //for resort, the implicit element variables belong to SOS2 FlexDOF()353 virtual int FlexDOF() const {return SOS();} 354 355 //# Timeint specific derived functions: for discontinuities StartTimeStep()356 virtual void StartTimeStep() {}; //function is called when computation of time step is started EndTimeStep()357 virtual void EndTimeStep() {}; //function is called when computation of time step is started ComputationFinished()358 virtual void ComputationFinished() {}; //function is called when computation is finished (e.g. in order to free memory, write results, close connections, etc.) 359 360 //$!LA 2011-02-16: NonlinStep was replaced by PostNewtonStep 361 // virtual double NonlinStep(double t) {return 0;}; PostNewtonStep(double t)362 virtual double PostNewtonStep(double t) {return 0;}; 363 //$!LA 2011-02-16: FixNonlinStep was replaced by PostprocessingStep 364 // virtual void FixNonlinStep() {}; PostprocessingStep()365 virtual void PostprocessingStep() {}; GetError()366 virtual double GetError() const 367 { 368 double err = 0; 369 for (int i=1; i<=SS()-IS(); i++) 370 { 371 //err += fabs(XG(i)); 372 err += Sqr(XG(i)); 373 } 374 return err; 375 }; 376 SetGlobalInitConditions(Vector & x_glob)377 virtual void SetGlobalInitConditions(Vector& x_glob) 378 { 379 for (int i=1; i<=SS(); i++) 380 { 381 // (AD) changed () to .Get() 382 x_glob(ltg.Get(i)) = x_init(i); 383 // x_glob(ltg(i)) = x_init(i); 384 } 385 } 386 SetGlobalInitData(Vector & data_glob)387 virtual void SetGlobalInitData(Vector& data_glob) 388 { 389 for (int i=1; i<=DataS(); i++) 390 { 391 // (AD) changed () to .Get() 392 data_glob(ltgdata.Get(i)) = data_init(i); 393 // data_glob(ltgdata(i)) = data_init(i); 394 } 395 } 396 Dim()397 virtual int Dim() const {return 2;} //$EDC$[funcaccess,readonly,EDCvarname="element_dimension",EDCfolder="Info",tooltiptext="dimensionality of element (2D/3D)"]//default value GetXact(int i)398 virtual const double& GetXact(int i) const {return mbs->GetXact(i);} GetXact(int i)399 virtual double& GetXact(int i) {return mbs->GetXact(i);} 400 //virtual const Vector& GetXact() const {return mbs->GetXact();} 401 //virtual Vector& GetXact() {return mbs->GetXact();} 402 // SetInitConditions(const Vector & x0)403 virtual void SetInitConditions(const Vector& x0) {x_init = x0;} GetXInit()404 virtual const Vector& GetXInit() const {return x_init;} 405 GetXInit(int i)406 virtual const double& GetXInit(int i) const {return x_init.Get(i);} // $ MSax 2013-07-12 : added 407 SetDataInit(const Vector & data_initI)408 virtual void SetDataInit(const Vector& data_initI) {data_init = data_initI;} GetDataInit()409 virtual const Vector& GetDataInit() const {return data_init;} 410 411 //compute things which should be actualized after every change of xact Update()412 virtual void Update() {}; 413 //local to global transformation, for positions, velocities, first order ODE and algebraic variables 414 415 virtual double GetOutput(double t, int i=1) const {return 0;} 416 ElementBandwidth()417 virtual int ElementBandwidth() const {return SS();} 418 419 420 //local to global transformation variables 421 // (AD) changed () to .Get() & Elem() LTG(int iloc)422 virtual int LTG(int iloc) const {return ltg.Get(iloc);} LTG(int iloc)423 virtual int& LTG(int iloc) {return ltg.Elem(iloc);} 424 // virtual int LTG(int iloc) const {return ltg(iloc);} 425 // virtual int& LTG(int iloc) {return ltg(iloc);} GetLTGArray()426 virtual const TArray<int>& GetLTGArray() const {return ltg;} AddLTG(int gi)427 virtual void AddLTG(int gi) {ltg.Add(gi);} LTGlength()428 virtual int LTGlength() const {return ltg.Length();} LTGreset()429 virtual void LTGreset() {ltg.Flush();} 430 431 //local to global transformation variables for data elements: 432 // (AD) changed () to .Get() & Elem() LTGdata(int iloc)433 virtual int LTGdata(int iloc) const {return ltgdata.Get(iloc);} LTGdata(int iloc)434 virtual int& LTGdata(int iloc) {return ltgdata.Elem(iloc);} 435 // virtual int LTGdata(int iloc) const {return ltgdata(iloc);} 436 // virtual int& LTGdata(int iloc) {return ltgdata(iloc);} GetLTGdataArray()437 virtual const TArray<int>& GetLTGdataArray() const {return ltgdata;} AddLTGdata(int gi)438 virtual void AddLTGdata(int gi) {ltgdata.Add(gi);} LTGdataReset()439 virtual void LTGdataReset() {ltgdata.Flush();} 440 441 //tell loads the element reference pointer LinkLoads()442 virtual void LinkLoads() 443 { 444 for (int i=1; i <= loads.Length(); i++) 445 { 446 //loads(i)->SetElement(this); //$ DR 2012-10: loads moved from element to mbs, old code 447 //mbs->GetLoad(loads(i)).SetElement(this); //$ DR 2012-10: loads moved from element to mbs, old code 448 mbs->GetLoad(loads(i)).SetMBS(mbs); //$ DR 2012-10: loads moved from element to mbs 449 mbs->GetLoad(loads(i)).SetDim(this->Dim()); //$ DR 2012-10: loads moved from element to mbs 450 } 451 } 452 //Link element to other elements if necessary (e.g. constraints) LinkToElements()453 virtual void LinkToElements() {}; 454 455 virtual void BuildDependencies(); 456 457 //$ YV 2011-08-02: 458 // the following function can be implemented by elements, 459 // which wish to perform some actions just before the system is assembled PreAssemble()460 virtual void PreAssemble() {} 461 IsDependent(int i)462 virtual int IsDependent(int i) const 463 { 464 if (i == 0 || !GetMBS()->UseDependencies()) return 1; 465 #ifdef _DEBUG 466 if (i > dependencies.Length()) 467 { 468 mbs->UO() << "Error in IsDependent!!!!\n"; 469 mbs->UO() << "deplen=" << dependencies.Length() << ", i=" << i << "\n"; 470 return 1; 471 } 472 #endif 473 return dependencies(i); 474 }; 475 476 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 477 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 478 //get system coordinate of actual state (must be set in EvalF, EvalG, EvalF2, EvalM) XG_dc(int iloc,TComputeDrawInitFlag flag)479 virtual const double XG_dc(int iloc, TComputeDrawInitFlag flag) const 480 { 481 if(flag & TCD_compute) 482 { 483 return GetXact(ltg.Get(iloc)); 484 } 485 else if(flag & TCD_draw) 486 { 487 return GetDrawValue(ltg.Get(iloc)); 488 } 489 else if(flag & TCD_reference_configuration) 490 { 491 return 0; 492 } 493 else if(flag & TCD_initial_values) 494 { 495 return GetXInit(iloc); 496 } 497 else 498 { 499 assert("XG_dc called with illegal flag" && 0); 500 return 0; 501 } 502 } 503 504 505 //get system coordinate of actual state (must be set in EvalF, EvalG, EvalF2, EvalM) XG(int iloc)506 virtual const double& XG(int iloc) const {return GetXact(ltg.Get(iloc));} 507 //EDC Vector XG(); //$EDC[funcaccess,readonly,EDCvarname="SOS_ODE_GC_poslevel",EDCfolder="Debug",condition=SOS()>0,vecstart=1,vecend=SOS(),tooltiptext="second order size state variables at position level (XG(i))"]//default value XG(int iloc)508 virtual double& XG(int iloc) {return GetXact(ltg(iloc));} XGP(int iloc)509 virtual const double& XGP(int iloc) const {return GetXact(ltg.Get(iloc+SOS()));} 510 //EDC Vector XGP(); //$EDC[funcaccess,readonly,EDCvarname="SOS_ODE_GC_vellevel",EDCfolder="Debug",condition=SOS()>0,vecstart=1,vecend=SOS(),tooltiptext="second order size state variables at velocity level (XGP(i))"]//default value XGP(int iloc)511 virtual double& XGP(int iloc) {return GetXact(ltg(iloc+SOS()));} 512 //EDC Vector XG(); //$EDC[funcaccess,readonly,EDCvarname="FOS_ODE_GC",EDCfolder="Debug",condition=ES()>0,vecstart=2*SOS()+1,vecend=2*SOS()+ES(),tooltiptext="first order size state variables (XG(2*SOS()+i))"]//default value 513 //EDC Vector XG(); //$EDC[funcaccess,readonly,EDCvarname="AE_lagrange_multipliers",EDCfolder="Debug",condition=IS()>0,vecstart=2*SOS()+ES()+1,vecend=SS(),tooltiptext="first order size state variables (XG(2*SOS()+ES()+i))"]//default value 514 //virtual const double& XGPP(int iloc) const {return mbs->GetVelocityAndAcceleration(ltg.Get(iloc+SOS()));} //this retrieves the acceleration values; only with implicit time integration! // $ MSax 2013-07-16 : removed XGPP(int iloc)515 virtual double XGPP(int iloc) const // $ MSax 2013-07-16 : added, returns the acceleration in dynamic case, returns zero in static case 516 { 517 if(mbs->DoStaticComputation()) return 0; 518 double xgpp = mbs->GetVelocityAndAcceleration(ltg.Get(iloc+SOS())); //this retrieves the acceleration values; only with implicit time integration! 519 return xgpp; 520 } 521 522 523 //global functions XGG(int iglob)524 virtual const double& XGG(int iglob) const {return GetXact(iglob);} XGG(int iglob)525 virtual double& XGG(int iglob) {return GetXact(iglob);} 526 527 //functions with values for drawing: 528 //get coordinate of actual state (must be set in EvalF, EvalG, EvalF2, EvalM) 529 // (AD) changed () to .Get() XGD(int iloc)530 virtual const double& XGD(int iloc) const {return GetDrawValue(ltg.Get(iloc));} XGD(int iloc)531 virtual double& XGD(int iloc) {return GetDrawValue(ltg.Get(iloc));} 532 // virtual const double& XGD(int iloc) const {return GetDrawValue(ltg(iloc));} 533 //virtual double& XGD(int iloc) {return GetDrawValue(ltg(iloc));} 534 // (AD) changed () to .Get() XGPD(int iloc)535 virtual const double& XGPD(int iloc) const {return GetDrawValue(ltg.Get(iloc+SOS()));} 536 /// virtual const double& XGPD(int iloc) const {return GetDrawValue(ltg(iloc+SOS()));} 537 //virtual double& GetDrawValue(int iloc) {return mbs->GetDrawValue(iloc);} GetDrawValue(int iloc)538 virtual const double& GetDrawValue(int iloc) const {return mbs->GetDrawValue(iloc);} GetDrawValue(int iloc)539 virtual double& GetDrawValue(int iloc) {return mbs->GetDrawValue(iloc);} 540 541 // (AD) changed () to .Get() XData(int iloc)542 virtual const double& XData(int iloc) const {return GetMBS()->GetDataAct(ltgdata.Get(iloc));} 543 // virtual const double& XData(int iloc) const {return GetMBS()->GetDataAct(ltgdata(iloc));} XData(int iloc)544 virtual double& XData(int iloc) {return GetMBS()->GetDataAct(ltgdata(iloc));} 545 546 // (AD) changed () to .Get() XDataD(int iloc)547 virtual const double& XDataD(int iloc) const {return GetMBS()->GetDataDraw(ltgdata.Get(iloc));} 548 // virtual const double& XDataD(int iloc) const {return GetMBS()->GetDataDraw(ltgdata(iloc));} XDataD(int iloc)549 virtual double& XDataD(int iloc) {return GetMBS()->GetDataDraw(ltgdata(iloc));} 550 551 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 552 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 553 554 //Get rot matrix in 2D and 3D 555 //rot matrix A describes the transformation of the local(body) to the 556 //global coordinate system, such that r_g=R_g+A_bg*u_b GetRotMatrix()557 virtual Matrix3D GetRotMatrix() const {return Matrix3D(1);}; GetRotMatrixP()558 virtual Matrix3D GetRotMatrixP() const {return Matrix3D();}; GetRotMatrix(const Vector3D & ploc)559 virtual Matrix3D GetRotMatrix(const Vector3D& ploc) const {return Matrix3D(1);}; 560 //virtual const Matrix& GetRotMatrix2D() {return rot;}; 561 562 //Get roation matrix in 2D or 3D GetRotMatrixD()563 virtual Matrix3D GetRotMatrixD() const {return Matrix3D(1);}; 564 GetRotMatrixInit()565 virtual Matrix3D GetRotMatrixInit() const {return Matrix3D(1);}; 566 //virtual const Matrix& GetRotMatrix2DD() {return rot;}; 567 //get a reference position of the body in 3d 568 569 // in specialized classes: keep this function up to date for consistency check with constraints, see Constraint::IsSuitableElement(..) 570 // add all kinematic access functions available to those of the base class 571 virtual TKinematicsAccessFunctions GetKinematicsAccessFunctions(int mode = 1) const 572 { 573 return TKinematicsAccessFunctions(TKAF_position+TKAF_displacement+TKAF_velocity); 574 } 575 576 577 //virtual Vector3D GetPos_dc(const Vector3D& ploc, TComputeDrawInitFlag flag) const 578 //{ 579 // return GetRefPos_dc(flag)+ploc; 580 //} 581 582 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 583 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 584 GetRefPos()585 virtual Vector3D GetRefPos() const 586 { 587 if (Dim()==3) 588 { 589 return Vector3D(XG(1),XG(2),XG(3)); 590 } 591 else 592 { 593 return Vector3D(XG(1),XG(2),0.); 594 } 595 } GetRefPosInit()596 virtual Vector3D GetRefPosInit() const 597 { 598 if (Dim()==3) 599 { 600 return Vector3D(x_init(1),x_init(2),x_init(3)); 601 } 602 else 603 { 604 return Vector3D(x_init(1),x_init(2),0.); 605 } 606 } GetRefVel()607 virtual Vector3D GetRefVel() const 608 { 609 if (Dim()==3) 610 { 611 return Vector3D(XGP(1),XGP(2),XGP(3)); 612 } 613 else 614 { 615 return Vector3D(XGP(1),XGP(2),0.); 616 } 617 } 618 //$ PG 2013-3-28: added GetRefConfPos to Element (for calling Constraint::GetBody3D(..).GetRefConfPos(...)) GetRefConfPos(const Vector3D & p_loc)619 virtual Vector3D GetRefConfPos(const Vector3D& p_loc) const 620 { 621 assert(0); 622 return Vector3D(0); //maby no assert and return GetRefPosInit()+ploc; ??? 623 } GetPos(const Vector3D & ploc)624 virtual Vector3D GetPos(const Vector3D& ploc) const 625 { 626 return GetRefPos()+ploc; 627 } GetVel(const Vector3D & ploc)628 virtual Vector3D GetVel(const Vector3D& ploc) const 629 { 630 return GetRefVel(); 631 } 632 //$ PG 2013-6-25: added for availablility in Element::GetFieldVariableValue(..), since displacement is added via Element::GetAvailableFieldVariables(..) GetDisplacement(const Vector3D & ploc)633 virtual Vector3D GetDisplacement(const Vector3D& ploc) const 634 { 635 return ploc; 636 } 637 //$ YV: sensors should use field variables, and the function below should be eliminated GetAcceleration(const Vector3D & ploc)638 virtual Vector3D GetAcceleration(const Vector3D& ploc) const // global acceleration vector used in MBSSensor 639 { 640 assert(0); 641 return Vector3D(0.); 642 } 643 644 //draw-functions: GetPosD(const Vector3D & ploc)645 virtual Vector3D GetPosD(const Vector3D& ploc) const 646 { 647 return GetRefPosD()+ploc; 648 } GetVelD(const Vector3D & ploc)649 virtual Vector3D GetVelD(const Vector3D& ploc) const 650 { 651 return GetRefVelD(); 652 } 653 //$ PG 2013-6-25: added for availablility in Element::GetFieldVariableValue(..), since displacement is added via Element::GetAvailableFieldVariables(..) GetDisplacementD(const Vector3D & ploc)654 virtual Vector3D GetDisplacementD(const Vector3D& ploc) const 655 { 656 return ploc; 657 } 658 GetRefPosD()659 virtual Vector3D GetRefPosD() const 660 { 661 if (Dim()==3) 662 { 663 return Vector3D(XGD(1),XGD(2),XGD(3)); 664 } 665 else 666 { 667 return Vector3D(XGD(1),XGD(2),0.); 668 } 669 } GetRefVelD()670 virtual Vector3D GetRefVelD() const 671 { 672 if (Dim()==3) 673 { 674 return Vector3D(XGPD(1),XGPD(2),XGPD(3)); 675 } 676 else 677 { 678 return Vector3D(XGPD(1),XGPD(2),0.); 679 } 680 } GetDOFPosD(int idof)681 virtual Vector3D GetDOFPosD(int idof) const //returns postion of i-th DOF 682 { 683 return GetRefPosD(); 684 } GetDOFDirD(int idof)685 virtual Vector3D GetDOFDirD(int idof) const //returns direction of action of i-th DOF 686 { 687 return Vector3D(0.,0.,0.); 688 } 689 690 NNodes()691 virtual int NNodes() const {return 0; }; //$EDC$[funcaccess,readonly,EDCvarname="number_of_nodes",EDCfolder="Info",tooltiptext="number of nodes in finite elements"] NodeNum(int i)692 virtual const int& NodeNum(int i) const {assert(0); return altshape;} NodeNum(int i)693 virtual int& NodeNum(int i) {assert(0); return altshape;} 694 GetNode(int i)695 virtual const Node& GetNode(int i) const {return GetMBS()->GetNode(NodeNum(i));} //get local node number i GetNode(int i)696 virtual Node& GetNode(int i) { 697 int j = NodeNum(i); 698 return GetMBS()->GetNode(j); 699 } //get local node number i GetNodeLocPos(int i)700 virtual Vector3D GetNodeLocPos(int i) const {assert(0); return Vector3D(0.,0.,0.);} //local position of node 701 GetNodePos(int i)702 virtual Vector3D GetNodePos(int i) const {return GetPos(GetNodeLocPos(i));} 703 GetNodePosD(int i)704 virtual Vector3D GetNodePosD(int i) const {return GetPosD(GetNodeLocPos(i));} 705 GetNodeVel(int i)706 virtual Vector3D GetNodeVel(int i) const {return GetVel(GetNodeLocPos(i));} 707 GetNodeVelD(int i)708 virtual Vector3D GetNodeVelD(int i) const {return GetVelD(GetNodeLocPos(i));} 709 GetNodeLocPos2D(int i)710 virtual Vector2D GetNodeLocPos2D(int i) const {assert(0); return Vector2D(0.,0.);} //local position of node 711 GetNodePos2D(int i)712 virtual Vector2D GetNodePos2D(int i) const {return GetPos2D(GetNodeLocPos2D(i));} 713 GetNodePos2DD(int i)714 virtual Vector2D GetNodePos2DD(int i) const {return GetPos2DD(GetNodeLocPos2D(i));} 715 GetNodeVel2D(int i)716 virtual Vector2D GetNodeVel2D(int i) const {return GetVel2D(GetNodeLocPos2D(i));} 717 718 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 719 //2D: GetRefPosInit2D()720 virtual Vector2D GetRefPosInit2D() const 721 { 722 return Vector2D(x_init(1),x_init(2)); 723 } GetRefPos2D()724 virtual Vector2D GetRefPos2D() const 725 { 726 return Vector2D(XG(1),XG(2)); 727 } GetRefVel2D()728 virtual Vector2D GetRefVel2D() const 729 { 730 return Vector2D(XGP(1),XGP(2)); 731 } 732 GetPos2D(const Vector2D & ploc)733 virtual Vector2D GetPos2D(const Vector2D& ploc) const 734 { 735 return GetRefPos2D()+ploc; 736 } GetVel2D(const Vector2D & ploc)737 virtual Vector2D GetVel2D(const Vector2D& ploc) const 738 { 739 return GetRefVel2D(); 740 } 741 //$ PG 2013-6-25: added for availablility in Element::GetFieldVariableValue(..), since displacement is added via Element::GetAvailableFieldVariables(..) GetDisplacement2D(const Vector2D & ploc)742 virtual Vector2D GetDisplacement2D(const Vector2D& ploc) const 743 { 744 return ploc; 745 } 746 747 //draw-functions 2D: GetRefPos2DD()748 virtual Vector2D GetRefPos2DD() const 749 { 750 return Vector2D(XGD(1),XGD(2)); 751 } GetRefVel2DD()752 virtual Vector2D GetRefVel2DD() const 753 { 754 return Vector2D(XGPD(1),XGPD(2)); 755 } 756 GetPos2DD(const Vector2D & ploc)757 virtual Vector2D GetPos2DD(const Vector2D& ploc) const 758 { 759 return GetRefPos2DD()+ploc; 760 } GetVel2DD(const Vector2D & ploc)761 virtual Vector2D GetVel2DD(const Vector2D& ploc) const 762 { 763 return GetRefVel2DD(); 764 } 765 //$ PG 2013-6-25: added for availablility in Element::GetFieldVariableValue(..), since displacement is added via Element::GetAvailableFieldVariables(..) GetDisplacement2DD(const Vector2D & ploc)766 virtual Vector2D GetDisplacement2DD(const Vector2D& ploc) const 767 { 768 return ploc; 769 } 770 771 772 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 773 774 775 //Add load to Element 776 //$ DR 2012-10:[ loads moved from element to mbs 777 // old code: 778 //virtual void AddLoad(const MBSLoad& li); 779 //virtual int NLoads() const {return loads.Length();} //EDC$[funcaccess,readonly,EDCvarname="number_of_loads",EDCfolder="Info"] 780 //virtual const MBSLoad& GetLoad(int i) const {return *loads(i);} 781 //virtual MBSLoad& GetLoad(int i) {return *loads(i);} 782 //virtual void DeleteLoad(int i) 783 //{ 784 // delete loads(i); 785 // for (int j=i; j < NLoads(); j++) 786 // { 787 // loads(j) = loads(j+1); 788 // } 789 // loads.SetLen(loads.Length()-1); 790 //} 791 792 // new code: 793 virtual void AddLoad(const MBSLoad& li); NLoads()794 virtual int NLoads() const {return loads.Length();} //$EDC$[funcaccess,readonly,EDCvarname="number_of_loads",EDCfolder="Info"] GetLoad(int i)795 virtual const MBSLoad& GetLoad(int i) const {return mbs->GetLoad(loads(i));} GetLoad(int i)796 virtual MBSLoad& GetLoad(int i) {return mbs->GetLoad(loads(i));} 797 virtual void DeleteLoad(int i); GetLoadNr(int i)798 virtual int GetLoadNr(int i) {return loads(i);} GetLoadNrs()799 virtual TArray<int> GetLoadNrs() {return loads;} SetLoadNrs(TArray<int> load_nrs)800 virtual void SetLoadNrs(TArray<int> load_nrs) {loads = load_nrs;} 801 802 //$ DR 2012-10:] loads moved from element to mbs 803 804 805 806 //Add constraint, do not copy!!! AddConstraint(Constraint * c,int lind)807 virtual void AddConstraint(Constraint* c, int lind) 808 { 809 constraints.Add(c); 810 constraintindices.Add(lind); 811 812 int found = 0; 813 for (int i=1; i <= constraints_nodouble.Length(); i++) 814 { 815 if (constraints_nodouble(i) == c) {found = 1; break;} 816 } 817 if (!found) constraints_nodouble.Add(c); 818 } RemoveConstraints()819 virtual void RemoveConstraints() 820 { 821 constraints.SetLen(0); 822 constraints_nodouble.SetLen(0); 823 constraintindices.SetLen(0); 824 } AddSensor(int sensornum)825 virtual int AddSensor(int sensornum) {return sensors.Add(sensornum);} NSensors()826 virtual int NSensors() const {return sensors.Length();} 827 GetSensorNum(int i)828 virtual int GetSensorNum(int i) const {return sensors(i);} SetSensorNum(int i,int sensnum)829 virtual void SetSensorNum(int i, int sensnum) {sensors(i) = sensnum;} GetSensor(int i)830 virtual const Sensor& GetSensor(int i) const {return mbs->GetSensor(sensors(i));} GetSensor(int i)831 virtual Sensor & GetSensor(int i) {return mbs->GetSensor(sensors(i));} GetSpecialSensorValue(int nr,double time)832 virtual double GetSpecialSensorValue(int nr, double time) const {return 0;} // DR 2012-01-12: has to be overwritten in derived class 833 834 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 835 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 836 //access functions for constraints 837 virtual void AddElement(int en); NE()838 virtual int NE() const {return elements.Length();} 839 GetElnum(int i)840 virtual int GetElnum(int i) const {return elements(i);} SetElnum(int i,int elnum)841 virtual void SetElnum(int i, int elnum) {elements(i) = elnum;} 842 GetElem(int i)843 virtual const Element& GetElem(int i) const {return mbs->GetElement(elements(i));} GetElem(int i)844 virtual Element& GetElem(int i) {return mbs->GetElement(elements(i));} 845 GetDirectFeedThroughElements(TArray<int> & elnums)846 virtual void GetDirectFeedThroughElements(TArray<int>& elnums) const {} //add all elements which are depending on this element by direct feed through to list 847 848 // functions to catch the input - implemented in derived class RespondToKey(int key)849 virtual int RespondToKey(int key) {return false;} 850 851 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 852 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 853 NC()854 virtual int NC() const {return constraints.Length();} GetConstraint(int i)855 virtual Constraint* GetConstraint(int i) const {return constraints(i);} ExternalSetRefCoord(double val)856 virtual void ExternalSetRefCoord(double val) {} //assigns an external input value to an element data --> used in rotactor and in LinearTransformation element 857 GetSize()858 virtual const Vector3D& GetSize() const {assert(0); return col;} 859 GetMaterialNum()860 virtual const int& GetMaterialNum() const {return materialnum;} GetMaterialNum()861 virtual int& GetMaterialNum() {return materialnum;} 862 virtual void SetMaterialNum(int mnum); 863 virtual void AddMaterial(const Material& m); GetMaterial()864 virtual const Material& GetMaterial() const {return GetMBS()->GetMaterial(materialnum);} GetMaterial()865 virtual Material& GetMaterial() {return GetMBS()->GetMaterial(materialnum);} 866 867 //virtual double GetRho() const {return rho;} //old version, should vanish in new or revised elements! 868 virtual const double& Rho() const; //link to material data 869 virtual const double& Em() const; //link to material data 870 virtual const double& Nu() const; //link to material data 871 virtual double& Rho(); //link to material data 872 virtual double& Em(); //link to material data 873 virtual double& Nu(); //link to material data 874 GetMass()875 virtual double GetMass() const {return mass;} GetVolume()876 virtual double GetVolume() const {return mass/Rho();} GetMassDamping()877 virtual double GetMassDamping() const {return damping_m;} SetMassDamping(double d)878 virtual void SetMassDamping(double d) {damping_m = d;} 879 880 //for loads, rigid bodies: GetDuxDq(Vector & dudq)881 virtual void GetDuxDq(Vector& dudq) {assert(0);}; GetDuyDq(Vector & dudq)882 virtual void GetDuyDq(Vector& dudq) {assert(0);}; GetDuzDq(Vector & dudq)883 virtual void GetDuzDq(Vector& dudq) {assert(0);}; GetDrotxDq(Vector & dudq)884 virtual void GetDrotxDq(Vector& dudq) {assert(0);}; GetDrotyDq(Vector & dudq)885 virtual void GetDrotyDq(Vector& dudq) {assert(0);}; GetDrotzDq(Vector & dudq)886 virtual void GetDrotzDq(Vector& dudq) {assert(0);}; 887 //volume loads for flexible bodies: GetIntDuDq(Matrix & dudq)888 virtual void GetIntDuDq(Matrix& dudq) {assert(0);}; //in fact it is DuDq Transposed 889 // gravity load for flexible bodies GetIntRhoDuDq(Matrix & rhodudq)890 virtual void GetIntRhoDuDq(Matrix& rhodudq) //in fact it is Int_V Rho DuDq Transposed 891 { 892 GetIntDuDq(rhodudq); 893 rhodudq *= Rho(); 894 } GetIntDkappaDq2D(Vector & dudq)895 virtual void GetIntDkappaDq2D(Vector& dudq) {assert(0);}; 896 // centrifugal load for flexible bodies GetIntDuDqFCentrifugal(Matrix & dudq,const Vector3D & omega,const Vector3D & r0)897 virtual void GetIntDuDqFCentrifugal(Matrix& dudq, const Vector3D& omega, const Vector3D& r0) {assert(0);}; 898 SetAltShape(int truefalse)899 virtual void SetAltShape(int truefalse) {altshape = truefalse;} GetAltShape()900 virtual int GetAltShape() const {return altshape;} 901 NGeomElements()902 virtual int NGeomElements() {return drawelements.Length();} AddGeomElement(int i)903 virtual int AddGeomElement(int i) {return drawelements.Add(i);} GetGeomElementNum(int i)904 virtual const int& GetGeomElementNum(int i) const {return drawelements(i);} GetGeomElementNum(int i)905 virtual int& GetGeomElementNum(int i) {return drawelements(i);} DeleteGeomElement(int i)906 virtual void DeleteGeomElement(int i) {drawelements.Erase(i);} GetGeomElement(int i)907 virtual GeomElement* GetGeomElement(int i) {return GetMBS()->GetDrawElement(drawelements(i));} Add(const GeomElement & de)908 virtual int Add(const GeomElement& de) //adds GeomElement to MBS and adds index to drawelments list 909 { 910 int i = GetMBS()->Add(de, GetOwnNum()); 911 return drawelements.Add(i); 912 913 //GeomElement* den = de.GetCopy(); //old version! 914 //drawelements.Add(den); 915 } 916 917 virtual Box3D GetBoundingBox() const; 918 virtual Box3D GetBoundingBoxD() const; 919 virtual Box2D GetBoundingBox2D() const; 920 virtual Box2D GetBoundingBox2DD() const; 921 GetElementBox()922 virtual Box3D GetElementBox() const 923 { 924 return Box3D(GetRefPos(),GetRefPos()); 925 } 926 GetElementBoxD()927 virtual Box3D GetElementBoxD() const 928 { 929 return Box3D(GetRefPosD(),GetRefPosD()); 930 } 931 GetElementBox2D()932 virtual Box2D GetElementBox2D() const 933 { 934 return Box2D((const Vector2D&) GetRefPos(), (const Vector2D&) GetRefPos()); 935 } 936 GetElementBox2DD()937 virtual Box2D GetElementBox2DD() const 938 { 939 return Box2D((const Vector2D&) GetRefPosD(), (const Vector2D&) GetRefPosD()); 940 } 941 942 virtual void DrawElementAdd(); 943 944 // the derived classes may fill in this list to inform the system 945 // which field variables they can produce for contour plotting and for sensors 946 virtual void GetAvailableFieldVariables(TArray<FieldVariableDescriptor> & variables); 947 // computes a value of a given field variable at a given point; 948 // flagD is true for drawing (contour plotting) and false for sensing (computation), 949 // and the actual computation should be based either on XG() or XGD() 950 // three versions are available with different definitions of the point, in which the field variable is computed: 951 // 3D local coordinates, 2D local coordinates, and local node number 952 953 //$ PG 2013-6-25:[ return position by GetPos(), displacement by GetDisplacemenet(), and velocity by GetVel(), those functions have to be implemented in the inherited element classes, though. 954 //virtual double GetFieldVariableValue(const FieldVariableDescriptor & fvd, const Vector3D & local_position, bool flagD) { assert(0); return 0; } 955 //virtual double GetFieldVariableValue(const FieldVariableDescriptor & fvd, const Vector2D & local_position, bool flagD) { assert(0); return 0; } 956 virtual double GetFieldVariableValue(const FieldVariableDescriptor & fvd, const Vector3D & local_position, bool flagD); 957 virtual double GetFieldVariableValue(const FieldVariableDescriptor & fvd, const Vector2D & local_position, bool flagD); GetFieldVariableValue(const FieldVariableDescriptor & fvd,const double & local_position,bool flagD)958 virtual double GetFieldVariableValue(const FieldVariableDescriptor & fvd, const double & local_position, bool flagD) { assert(0); return 0; }; //$ DR 2013-06: added 959 //$ PG 2013-6-25:] GetFieldVariableValue(const FieldVariableDescriptor & fvd,int local_node_number,bool flagD)960 virtual double GetFieldVariableValue(const FieldVariableDescriptor & fvd, int local_node_number, bool flagD) { assert(0); return 0; } //$ YV 2012-06: added 961 962 DrawElementPreProc()963 virtual void DrawElementPreProc() {}; 964 965 virtual void DrawElement(); 966 virtual void DrawElement2D(); 967 virtual void DrawElementLocalFrame(); 968 virtual void DrawElementVelocityVector(); DrawElementFlag()969 virtual const int& DrawElementFlag() const {return draw_element;} DrawElementFlag()970 virtual int& DrawElementFlag() {return draw_element;} 971 GetCol()972 virtual const Vector3D& GetCol() const {return col;} GetCol()973 virtual Vector3D& GetCol() {return col;} 974 GetAngularMomentum()975 virtual Vector3D GetAngularMomentum() const { return GetAngularMomentum(Vector3D(0.,0.,0.)); } GetAngularMomentum(const Vector3D & p_ref)976 virtual Vector3D GetAngularMomentum(const Vector3D& p_ref) const { assert(0); return Vector3D(0.,0.,0.); } 977 978 //$ YV 2013-01-03: for constraint elements this function adds the global numbers of constrained 979 // dofs to the list (i.e. the old contents of the list remains in there!); 980 // the function is used by the eigenmodes solver GetGlobalConstrainedDOFs(TArray<int> & dofs)981 virtual void GetGlobalConstrainedDOFs(TArray<int>& dofs) const {} 982 983 //$ AD 2013-01-20: added functions to write name of ouput (load) to list available in InputOutputElement GetNOutputs()984 virtual int GetNOutputs() const {return 0;} SetOutputName(int output_nr,mystr & str)985 virtual void SetOutputName(int output_nr, mystr& str) {} 986 987 //$ AD 2013-07-08: added functions to maniputlate the connection lines between the IOElements InsertConNode(int list_idx,int input_nr,Vector2D & node)988 virtual void InsertConNode(int list_idx, int input_nr, Vector2D& node) {} DeleteConNode(int list_idx)989 virtual void DeleteConNode(int list_idx) {} 990 //$AD 2013-07-10: Change Position of a single element (MBS element, conNode, ...) MoveConNode2D(int list_idx,double delta_x,double delta_y)991 virtual void MoveConNode2D(int list_idx, double delta_x, double delta_y) {} MoveElement(double delta_x,double delta_y,double delta_z)992 virtual void MoveElement(double delta_x, double delta_y, double delta_z) {} 993 994 //$ AD 2013-02-22: added functions to manipulate an internal variable of InputOutputElement - increase or decrease the current value of IOGUIResponse Increase()995 virtual void Increase() {}; Decrease()996 virtual void Decrease() {}; 997 998 protected: 999 mystr elementname; //$EDC$[varaccess,EDCvarname="name",EDCfolder="",tooltiptext="name of the element"] 1000 mutable MBS* mbs; //$!JG 2011-02: switched to mutable in order to modify e.g. in GetConstraintDrift 1001 Vector x_init; //initial conditions 1002 Vector data_init; //initial conditions for data variables 1003 //Matrix rot; //temporal matrix, 2D 1004 1005 TArray<int> ltg; //local to global dof reference vector; 1006 TArray<int> ltgdata; //local to global data reference vector; 1007 1008 TArray<int> constraintindices; //index of this element in the constraint (mostly 1 or 2) 1009 TArray<char> dependencies; //contains dependencies to coordinates 1010 int elnum; //$EDC$[varaccess,EDCvarname="element_number",EDCfolder="",readonly,tooltiptext="number of the element in the mbs"] //number in MBS-system 1011 TMBSElement type; 1012 1013 1014 //TArray<MBSLoad*> loads; //$ DR 2012-10: loads moved from element to mbs 1015 TArray<int> loads; //$EDC$[varaccess,EDCvarname="loads",variable_length_vector,condition=!IsType(TConstraint),tooltiptext="Set loads attached to this element: 'nr_load1, nr_load2, ...' or empty"] 1016 TArray<Constraint*> constraints; //this constraint list might be double! 1017 TArray<Constraint*> constraints_nodouble; //one constraint is only added once to an element 1018 1019 TArray<int> sensors; //$EDC$[varaccess,EDCvarname="sensors",EDCfolder="Info",variable_length_vector,tooltiptext="attached sensors",readonly]//if the element equations of motion are dependent on sensor values (e.g. as an input), then add the sensor numbers here 1020 //this improves the Jacobian; otherwise, convergence might fail 1021 TArray<int> elements; //dependent elements (mainly for contraints, but also for follower loads!) 1022 1023 //density for gravitational force ==> erase!!!: 1024 //double rho; //JG+DR: should not be accessible! EDC$[varaccess,EDCvarname="density",EDCfolder="Physics",condition=GetMaterialNum()==0,tooltiptext="density of simple bodies, which do not obtain density from material"] 1025 1026 int materialnum; //JG+DR: should not be accessible! EDC$[varaccess,EDCvarname="material_number",EDCfolder="Physics",tooltiptext="material number which contains the main material properties of the body or element"] //material number in MBS which contains material information 1027 1028 double mass; //JG+DR: should not be accessible! EDC$[varaccess,EDCvarname="mass",EDCfolder="Physics",tooltiptext="for rigid bodies and FFRF bodies, not available for constraints"] //mass, for rigid bodies and FFRF bodies // DR 2012-07 1029 double damping_m; //JG+DR: should not be accessible! EDC[varaccess,EDCvarname="mass_prop_damping",EDCfolder="Physics",tooltiptext="only for bodies, not available for constraints"] //damping constant for mass proportional damping 1030 1031 Vector3D col; //$EDC$[varaccess,EDCvarname="RGB_color",EDCfolder="Graphics",tooltiptext="[red, green, blue] color of element, range = 0..1, use default color:[-1,-1,-1]"] //color for body 1032 1033 //geomelements are managed in MBS main class 1034 TArray<int> drawelements; //$EDC$[varaccess,EDCvarname="geom_elements",EDCfolder="Graphics",variable_length_vector,condition=!IsType(TConstraint),tooltiptext="Set Geometric elements to represent body 'geomelem1, geomelem2, ...' or empty"] 1035 int altshape; //$EDC$[varaccess,EDCvarname="use_alternative_shape",EDCfolder="Graphics",int_bool,condition=!IsType(TConstraint),tooltiptext="Graphical representation of element with geom-objects that are attached to the element"] 1036 1037 int draw_element; //$EDC$[varaccess,EDCvarname="show_element",EDCfolder="Graphics",int_bool,tooltiptext="Flag to draw element"] 1038 }; //$EDC$[endclass,Element] 1039 1040 1041 typedef enum {TRWElementDataNoAccess=0, TRWElementDataRead=1, TRWElementDataWrite=2, TRWElementDataReadWrite=3} TReadWriteElementDataVariableAccess; 1042 1043 class ReadWriteElementDataVariableType 1044 { 1045 public: ReadWriteElementDataVariableType()1046 ReadWriteElementDataVariableType():variable_name(),comp1(0),comp2(0),value(0.),tooltiptext()/*, RWaccess(1)*/ {RWaccess = TRWElementDataRead;}; // default constructor ReadWriteElementDataVariableType(const mystr & full_string)1047 ReadWriteElementDataVariableType(const mystr& full_string) { GetVariableNameAndComponents(full_string); }; 1048 ReadWriteElementDataVariableType(const mystr& name, int c1, int c2, double val, const mystr& tooltip, TReadWriteElementDataVariableAccess RWaccessI = TRWElementDataRead) // used in function GetAvailableSpecialValuesAuto 1049 { 1050 RWaccess = RWaccessI; 1051 variable_name = name; 1052 comp1 = c1; 1053 comp2 = c2; 1054 value = 0.; 1055 tooltiptext = tooltip; 1056 } 1057 //possibly needed in future: copy-constructor, destructor 1058 // convert an full string that CAN contain [] into desired parts ( variable_name / component1 / component2 ) 1059 // example "MassMatrix[2,2]" is converted to: { variable_name = "MassMatrix", comp1 = 2, comp2 = 2 } 1060 public: GetVariableNameAndComponents(const mystr & full_string)1061 int GetVariableNameAndComponents(const mystr& full_string) //initialization from string (for access by sensor) 1062 { 1063 //initialize with defaults: 1064 value = 0.; 1065 RWaccess = TRWElementDataRead; 1066 //tooltiptext = ""; //not needed! 1067 1068 // --> MYSTR function 1069 variable_name = full_string; 1070 int n_comp = variable_name.Trim_And_Get_Indices(comp1,comp2); 1071 1072 if(n_comp>=0) {return 1;} // parsed successfully 1073 else {return 0;} // error 1074 } 1075 1076 // YV 2012-11-5: comparison, e.g. for checking consistensy of a sensor 1077 // JG 2013-01-11: changed to function which checks range 1078 //int operator==(const ReadWriteElementDataVariableType & RWdata) const 1079 //{ 1080 // return variable_name == RWdata.variable_name && comp1 == RWdata.comp1 && comp2 == RWdata.comp2; 1081 //} 1082 //this function checks, if the RWdata provided fits into the name and range (comp1 for maximum range of comp1 and comp2 for maximum range of comp2) 1083 //and if RWaccess fits IsAvailable(const ReadWriteElementDataVariableType & RWdata)1084 int IsAvailable(const ReadWriteElementDataVariableType & RWdata) const 1085 { 1086 if (variable_name == RWdata.variable_name) 1087 { 1088 if (comp1 != 0) 1089 { 1090 if (RWdata.comp1 </*=*/ 0 || RWdata.comp1 > comp1) return 0; 1091 } 1092 if (comp2 != 0) 1093 { 1094 if (RWdata.comp2 </*=*/ 0 || RWdata.comp2 > comp2) return 0; 1095 } 1096 1097 if (((int)RWdata.RWaccess & (int)RWaccess) == 0) return 0; //does not work for case that (int)RWdata.RWaccess == 3, but should not occur 1098 1099 return 1; 1100 } 1101 1102 return 0; 1103 } 1104 1105 public: 1106 mystr variable_name; //name of variable for read/write access 1107 int comp1; //component for vector variables or matrix row 1108 int comp2; //second component for matrix variable (column) 1109 TReadWriteElementDataVariableAccess RWaccess; //binary flags: 0==no access, 1==read access, 2==write access, 3==read and write access 1110 double value; //read/write value 1111 mystr tooltiptext; //short description matching the variable name - to be filled ONLY by function GetAvailableSpecialValues 1112 double time; // time for evaluation 1113 }; 1114 1115 #endif