1 /* $Header: /var/cvs/mbdyn/mbdyn/mbdyn-1.0/mbdyn/aero/aeroext.cc,v 1.44 2017/01/12 14:45:58 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 #include "mbconfig.h" /* This goes first in every *.c,*.cc file */
33
34 #ifdef USE_EXTERNAL
35 #include "aeroext.h"
36 #include "myassert.h"
37 #include "except.h"
38 #include "dataman.h"
39 #include "external.h"
40
41 /* le label delle comunicazioni vengono costruite in questo modo:
42 *
43 * blockID*10 + 1 label e numero nodi aggiuntivi
44 * 2 posizioni
45 * 3 velocita'
46 * 4 forze
47 * 5 momenti
48 *
49 */
50
51 /* Aerodynamic External - begin */
AerodynamicExternal(unsigned int uLabel,const DofOwner * pDO,const DriveCaller * pDC,int NN,const StructNode ** ppN,const doublereal * pRefL,MPI::Intercomm * pIC,flag fOut,bool VF,bool MF)52 AerodynamicExternal::AerodynamicExternal(unsigned int uLabel,
53 const DofOwner *pDO,
54 const DriveCaller* pDC,
55 int NN,
56 const StructNode** ppN,
57 const doublereal* pRefL,
58 MPI::Intercomm* pIC,
59 flag fOut,
60 bool VF,
61 bool MF)
62 : Elem(uLabel, fOut),
63 AerodynamicElem(uLabel, pDO, fOut),
64 DriveOwner(pDC),
65 pdBuffer(NULL),
66 pdBufferVel(NULL),
67 NodeN(NN),
68 ppNode(ppN),
69 pRefLength(pRefL),
70 OffN(3),
71 pOffsetVectors(NULL),
72 pInterfComm(pIC),
73 VelFlag(VF),
74 MomFlag(MF),
75 SentFlag(false),
76 pForce(NULL),
77 pMoms(NULL),
78 pLabList(NULL)
79 {
80 /* costruisce la mat 3xN degli offset */
81 SAFENEWWITHCONSTRUCTOR(pOffsetVectors, Mat3xN, Mat3xN(OffN, 0.));
82 for (int i=1; i<=3; i++) pOffsetVectors->Put(i,i,1.);
83 ConstructAndInitialize();
84 }
85
AerodynamicExternal(unsigned int uLabel,const DofOwner * pDO,const DriveCaller * pDC,int NN,const StructNode ** ppN,const doublereal * pRefL,MPI::Intercomm * pIC,int ON,Mat3xN * OV,flag fOut,bool VF,bool MF)86 AerodynamicExternal::AerodynamicExternal(unsigned int uLabel,
87 const DofOwner* pDO,
88 const DriveCaller* pDC,
89 int NN,
90 const StructNode** ppN,
91 const doublereal* pRefL,
92 MPI::Intercomm* pIC,
93 int ON,
94 Mat3xN* OV,
95 flag fOut,
96 bool VF,
97 bool MF)
98 : Elem(uLabel, fOut),
99 AerodynamicElem(uLabel, pDO, fOut),
100 DriveOwner(pDC),
101 pdBuffer(NULL),
102 pdBufferVel(NULL),
103 NodeN(NN),
104 ppNode(ppN),
105 pRefLength(pRefL),
106 OffN(ON),
107 pOffsetVectors(OV),
108 pInterfComm(pIC),
109 VelFlag(VF),
110 MomFlag(MF),
111 SentFlag(false),
112 pForce(NULL),
113 pMoms(NULL),
114 pLabList(NULL)
115 {
116 ConstructAndInitialize();
117 }
118
~AerodynamicExternal(void)119 AerodynamicExternal::~AerodynamicExternal(void)
120 {
121 if (pdBuffer != NULL) SAFEDELETE(pdBuffer);
122 if (pdBufferVel != NULL) SAFEDELETE(pdBufferVel);
123 if (ppNode != NULL) SAFEDELETEARR(ppNode);
124 if (pRefLength != NULL) SAFEDELETEARR(pRefLength);
125 if (pForce != NULL) SAFEDELETEARR(pForce);
126 if (pMoms != NULL) SAFEDELETEARR(pMoms);
127 if (pLabList != NULL) SAFEDELETEARR(pLabList);
128 }
129
130 void
ConstructAndInitialize(void)131 AerodynamicExternal::ConstructAndInitialize(void)
132 {
133 /* dimensiona l'array per il buffer */
134 DEBUGCOUTFNAME("AerodynamicExternal::ConstructAndInitialize");
135 if (bToBeOutput()) {
136 SAFENEWWITHCONSTRUCTOR(pForce, Mat3xN, Mat3xN(NodeN+1,0.));
137 SAFENEWWITHCONSTRUCTOR(pMoms, Mat3xN, Mat3xN(NodeN+1,0.));
138 }
139 SAFENEWWITHCONSTRUCTOR(pdBuffer, MyVectorHandler,
140 MyVectorHandler(8+3*NodeN*(OffN+1)));
141 if (VelFlag || MomFlag) SAFENEWWITHCONSTRUCTOR(pdBufferVel, MyVectorHandler,
142 MyVectorHandler(3*NodeN*(OffN+1)));
143
144 /* invio delle posizioni iniziali dei nodi per il calcolo della matrice di interpolazione */
145 SAFENEWARR(pLabList, unsigned int, 2*NodeN);
146 Mat3xN MatTmp(OffN,0.);
147 Vec3 X;
148 for (int i=0; i < NodeN; i++) {
149 pLabList[i] = ppNode[i]->GetLabel();
150 pLabList[i+NodeN] = 1+OffN;
151 X = ppNode[i]->GetXCurr();
152 pdBuffer->Put(i*(OffN+1)*3+1, X);
153 for (int j=1; j <= OffN; j++) {
154 MatTmp.LeftMult(ppNode[i]->GetRCurr(), *pOffsetVectors);
155 MatTmp *= pRefLength[i];
156 pdBuffer->Put(i*(OffN+1)*3+j*3+1, X + MatTmp.GetVec(j));
157 }
158 }
159 /* il rank del processo con cui si comunica e' zero perche' il comuncatire del
160 codice di interfaccia e' sempre fatto da una macchina sola */
161 pInterfComm->Send(pLabList, 2*NodeN, MPI::UNSIGNED, 0,(this->GetLabel())*10+1);
162 pInterfComm->Send(pdBuffer->pdGetVec(), 3*NodeN*(OffN+1), MPI::DOUBLE, 0,(this->GetLabel())*10+2);
163 }
164
165
166 /* invia posizione e velocita' predetti */
167 void
AfterPredict(VectorHandler & X,VectorHandler & XP)168 AerodynamicExternal::AfterPredict(VectorHandler& X, VectorHandler& XP)
169 {
170 DEBUGCOUTFNAME("AerodynamicExternal::AfterPredict");
171 //std::cout << "AerodynamicExternal::AfterPredict" << std::endl;
172 Send(X, XP);
173 }
174
175 /* invia posizione e velocita' predetti */
176 void
Update(const VectorHandler & X,const VectorHandler & XP)177 AerodynamicExternal::Update(const VectorHandler& X, const VectorHandler& XP)
178 {
179 DEBUGCOUTFNAME("AerodynamicExternal::Update");
180 //std::cout << "AerodynamicExternal::Update" << std::endl;
181 Send(X, XP);
182 //std::cout << "done" << std::endl;
183 }
184
185 void
Send(const VectorHandler & X,const VectorHandler & XP)186 AerodynamicExternal::Send(const VectorHandler& X, const VectorHandler& XP)
187 {
188 if (!SentFlag) {
189 //std::cout << "AerodynamicExternal::Send " << SentFlag << std::endl;
190 Vec3 Vinf(0.);
191 Vec3 Xc(0.), V, W, Xt;
192 Mat3x3 R;
193 doublereal rho, c, p, T;
194 GetAirProps(Xc, rho, c, p, T);
195 fGetAirVelocity(Vinf,Xc);
196 /* current time */
197 pdBuffer->PutCoef(1, DriveOwner::dGet());
198 pdBuffer->Put(2, Vinf);
199 pdBuffer->PutCoef(5, rho);
200 pdBuffer->PutCoef(6,c);
201 pdBuffer->PutCoef(7,p);
202 pdBuffer->PutCoef(8,T);
203
204 Mat3xN MatTmp(OffN,0.);
205 for (int i=0; i < NodeN; i++) {
206 Xc = ppNode[i]->GetXCurr();
207 pdBuffer->Put(8+i*(OffN+1)*3+1, Xc);
208 //std::cout << "Block " << (this->GetLabel()) << " Node " << i << " Pos " << Xc << std::endl;
209 if (VelFlag) {
210 V = ppNode[i]->GetVCurr();
211 W = ppNode[i]->GetWCurr();
212 pdBufferVel->Put(i*(OffN+1)*3+1, V);
213 }
214 R = ppNode[i]->GetRCurr();
215 //std::cout << "Vel " << V << " R " << R << " W " << W << std::endl;
216 for (int j=1; j <= OffN; j++) {
217 MatTmp.LeftMult(R, *pOffsetVectors);
218 MatTmp *= pRefLength[i];
219 Xt = MatTmp.GetVec(j);
220 pdBuffer->Put(8+i*(OffN+1)*3+j*3+1, Xc + Xt);
221 if (VelFlag) pdBufferVel->Put(i*(OffN+1)*3+j*3+1, V + W.Cross(Xt));
222 //std::cout << " OffVec " << j << ": Pos " << Xt << " Vel " << V + W.Cross(Xt) << std::endl;
223 }
224 }
225 pInterfComm->Send(pdBuffer->pdGetVec(), 8+3*NodeN*(OffN+1), MPI::DOUBLE, 0,(this->GetLabel())*10+2);
226 if (VelFlag) pInterfComm->Send(pdBufferVel->pdGetVec(), 3*NodeN*(OffN+1), MPI::DOUBLE, 0,(this->GetLabel())*10+3);
227 SentFlag = true;
228 }
229 }
230
231 SubVectorHandler&
AssRes(SubVectorHandler & WorkVec,doublereal dCoef,const VectorHandler & XCurr,const VectorHandler & XPrimeCurr)232 AerodynamicExternal::AssRes(SubVectorHandler& WorkVec,
233 doublereal dCoef,
234 const VectorHandler& XCurr,
235 const VectorHandler& XPrimeCurr)
236 {
237 DEBUGCOUTFNAME("AerodynamicExternal::AssRes");
238 /* verifica il completamento dei send */
239 //std::cout << "AerodynamicExternal::AssRes " << SentFlag << std::endl;
240
241 Send(XCurr, XPrimeCurr);
242
243 WorkVec.ResizeReset(6*NodeN);
244
245 for (int i=0; i < NodeN; i++) {
246 integer iFirstIndex = ppNode[i]->iGetFirstMomentumIndex();
247 for (int iCnt = 1; iCnt <= 6; iCnt++) {
248 WorkVec.PutRowIndex(i*6 + iCnt, iFirstIndex + iCnt);
249 }
250 }
251
252 pInterfComm->Recv(pdBuffer->pdGetVec(), 3*NodeN*(OffN+1), MPI::DOUBLE, 0,(this->GetLabel())*10+4);
253 if(MomFlag) pInterfComm->Recv(pdBufferVel->pdGetVec(), 3*NodeN*(OffN+1), MPI::DOUBLE, 0,(this->GetLabel())*10+5);
254 /* calcola le forze su ciascun nodo */
255 /* le forze e i momenti sono gia' espressi nel sistema di riferimento global */
256 if (bToBeOutput()) {
257 pForce->Reset();
258 pMoms->Reset();
259 }
260 for (int i=0; i < NodeN; i++) {
261 Vec3 F(0.);
262 Vec3 Fo(0.);
263 Vec3 M(0.);
264 F[0] += pdBuffer->operator()(i*(OffN+1)*3+1);
265 F[1] += pdBuffer->operator()(i*(OffN+1)*3+2);
266 F[2] += pdBuffer->operator()(i*(OffN+1)*3+3);
267 if (MomFlag) {
268 M[0] += pdBufferVel->operator()(i*(OffN+1)*3+1);
269 M[1] += pdBufferVel->operator()(i*(OffN+1)*3+2);
270 M[2] += pdBufferVel->operator()(i*(OffN+1)*3+3);
271 }
272 for (int j=1; j <= OffN; j++) {
273 Fo[0] = pdBuffer->operator()(i*(OffN+1)*3+j*3+1);
274 Fo[1] = pdBuffer->operator()(i*(OffN+1)*3+j*3+2);
275 Fo[2] = pdBuffer->operator()(i*(OffN+1)*3+j*3+3);
276 F += Fo;
277 if (MomFlag) {
278 M[0] += pdBufferVel->operator()(i*(OffN+1)*3+j*3+1);
279 M[1] += pdBufferVel->operator()(i*(OffN+1)*3+j*3+2);
280 M[2] += pdBufferVel->operator()(i*(OffN+1)*3+j*3+3);
281 }
282 M -= Fo.Cross( ppNode[i]->GetRCurr() * pOffsetVectors->GetVec(j) * pRefLength[i]);
283 }
284 if (bToBeOutput()) {
285 pForce->AddVec(1,F);
286 pMoms->AddVec(1,M);
287 pForce->PutVec(i+2,F);
288 pMoms->PutVec(i+2,M);
289 }
290 WorkVec.Add(i*6+1, F);
291 WorkVec.Add(i*6+4, M);
292 }
293 //std::cout << "FT " << pForce->GetVec(1) << std::endl;
294 //std::cout << "MT " << pMoms->GetVec(1) << std::endl;
295
296 SentFlag = false;
297 return WorkVec;
298 }
299
300 void
Output(OutputHandler & OH) const301 AerodynamicExternal::Output(OutputHandler& OH) const
302 {
303 if (bToBeOutput()) {
304 integer lab = -1*GetLabel();
305 std::ostream& out = OH.Externals()
306 << std::setw(8) << lab
307 << " " << pForce->GetVec(1)
308 << " " << pMoms->GetVec(1) << std::endl;
309 if (NodeN > 1) {
310 for (int i=1; i <= NodeN; i++) {
311 out << std::setw(8) << pLabList[i-1]
312 << " " << pForce->GetVec(i+1)
313 << " " << pMoms->GetVec(i+1) << std::endl;
314 }
315 }
316 }
317 }
318
319 Elem *
ReadAerodynamicExternal(DataManager * pDM,MBDynParser & HP,const DofOwner * pDO,unsigned int uLabel)320 ReadAerodynamicExternal(DataManager* pDM, MBDynParser& HP,
321 const DofOwner *pDO, unsigned int uLabel)
322 {
323 /* legge i dati d'ingresso e li passa al costruttore dell'elemento */
324 AerodynamicExternal* pEl = NULL;
325
326 unsigned int NodeN = HP.GetInt();
327
328 std::cout << "Aerodynamic External Block made of " << NodeN << " Nodes." << std::endl;
329
330 const StructNode** ppNodeList = NULL;
331 SAFENEWARR(ppNodeList, const StructNode*, NodeN);
332
333 for (unsigned int iN = 0; iN < NodeN; iN++) {
334 unsigned int NodeL = HP.GetInt();
335 const StructNode* pTmpNode = pDM->pFindNode<const StructNode, Node::STRUCTURAL>(NodeL);
336 if (pTmpNode == 0) {
337 std::cerr << "Aerodynamic External(" << uLabel
338 << "): structural node " << NodeL
339 << " at line " << HP.GetLineData()
340 << " not defined" << std::endl;
341 throw DataManager::ErrGeneric(MBDYN_EXCEPT_ARGS);
342 } else {
343 ppNodeList[iN] = pTmpNode;
344 }
345 }
346 doublereal* pRefLenght = NULL;
347 SAFENEWARR(pRefLenght, doublereal, NodeN);
348 if (HP.IsKeyWord("list")) {
349 for (unsigned int iN=0; iN < NodeN; iN++) {
350 pRefLenght[iN] = HP.GetReal();
351 }
352 } else {
353 pRefLenght[0] = HP.GetReal();
354 for (unsigned int iN=1; iN < NodeN; iN++) {
355 pRefLenght[iN] = pRefLenght[iN-1];
356 }
357 }
358
359 int OffN = 0;
360 Mat3xN* pOffVec = NULL;
361 if (HP.IsKeyWord("offset")) {
362 OffN = HP.GetInt();
363 SAFENEWWITHCONSTRUCTOR(pOffVec, Mat3xN, Mat3xN(OffN));
364 for (int iOf = 1; iOf <= OffN; iOf++) {
365 for (int iCnt = 1; iCnt <= 3; iCnt++) {
366 pOffVec->Put(iCnt, iOf, HP.GetReal());
367 }
368 }
369 }
370
371
372 int comm = 0;
373 if (HP.IsKeyWord("comm")) {
374 comm = HP.GetInt();
375 }
376 std::list<MPI::Intercomm>::iterator iComm = InterfaceComms.begin();
377 std::list<MPI::Intercomm>::iterator iComm_end = InterfaceComms.end();
378 if (iComm == iComm_end) {
379 std::cerr << "Aerodynamic External(" << uLabel
380 << "): there are no external communicators "
381 << "Aborting." << std::endl;
382 throw DataManager::ErrGeneric(MBDYN_EXCEPT_ARGS);
383 }
384
385 for (int i = 0; i++ < comm; ++iComm) {
386 if (iComm == iComm_end) {
387 std::cerr << "Aerodynamic External(" << uLabel
388 << "): external communicator " << comm
389 << "does not exist. Last communicator is " << i-1
390 << ". Aborting." << std::endl;
391 throw DataManager::ErrGeneric(MBDYN_EXCEPT_ARGS);
392 }
393 }
394 MPI::Intercomm* pInterC = &(*(iComm));
395
396 bool VelFlag = false;
397 bool MomFlag = false;
398
399 if (HP.IsKeyWord("velocity")) VelFlag = true;
400 if (HP.IsKeyWord("moment")) MomFlag = true;
401
402 flag fOut = pDM->fReadOutput(HP, Elem::AERODYNAMIC);
403 DriveCaller* pDC = NULL;
404 SAFENEWWITHCONSTRUCTOR(pDC, TimeDriveCaller, TimeDriveCaller(pDM->pGetDrvHdl()));
405
406 if (OffN) {
407 SAFENEWWITHCONSTRUCTOR(pEl,
408 AerodynamicExternal,
409 AerodynamicExternal(uLabel,
410 pDO,
411 pDC,
412 NodeN,
413 ppNodeList,
414 pRefLenght,
415 pInterC,
416 OffN,
417 pOffVec,
418 fOut,
419 VelFlag,
420 MomFlag));
421 } else {
422 SAFENEWWITHCONSTRUCTOR(pEl,
423 AerodynamicExternal,
424 AerodynamicExternal(uLabel,
425 pDO,
426 pDC,
427 NodeN,
428 ppNodeList,
429 pRefLenght,
430 pInterC,
431 fOut,
432 VelFlag,
433 MomFlag));
434 }
435
436 return pEl;
437 }
438
439
440 /* Aerodynamic External - end */
441
442
443 /* Aerodynamic External Modal - begin */
444
AerodynamicExternalModal(unsigned int uLabel,const DofOwner * pDO,const DriveCaller * pDC,Modal * pM,MPI::Intercomm * IC,flag fOut,bool VF,bool MF)445 AerodynamicExternalModal::AerodynamicExternalModal(unsigned int uLabel,
446 const DofOwner *pDO,
447 const DriveCaller* pDC,
448 Modal* pM,
449 MPI::Intercomm* IC,
450 flag fOut,
451 bool VF,
452 bool MF)
453 :Elem(uLabel, fOut),
454 AerodynamicElem(uLabel, pDO, fOut),
455 DriveOwner(pDC),
456 pdBuffer(NULL),
457 pdBufferVel(NULL),
458 pModal(pM),
459 ModalNodes(pM->uGetNFemNodes()),
460 pInterfComm(IC),
461 VelFlag(VF),
462 MomFlag(MF),
463 SentFlag(false),
464 pForce(NULL)
465 {
466 SAFENEWWITHCONSTRUCTOR(pdBuffer, MyVectorHandler,
467 MyVectorHandler(8+3*ModalNodes));
468
469 if (bToBeOutput()) {
470 integer NModes = pModal->uGetNModes();
471 SAFENEWWITHCONSTRUCTOR(pForce, MyVectorHandler,
472 MyVectorHandler(NModes));
473 }
474
475 if (VelFlag || MomFlag) SAFENEWWITHCONSTRUCTOR(pdBufferVel, MyVectorHandler,
476 MyVectorHandler(3*ModalNodes));
477
478 unsigned int* pLabList = NULL;
479 SAFENEWARR(pLabList, unsigned int, 2);
480 pLabList[0] = pModal->GetLabel();
481 pLabList[1] = ModalNodes;
482 /* il rank del processo con cui si comunica e' zero perche' il comuncatore del
483 codice di interfaccia e' sempre fatto da una macchina sola */
484 pInterfComm->Send(pLabList, 2, MPI::UNSIGNED, 0,(this->GetLabel())*10+1);
485
486 Mat3xN* pNodePos = pModal->GetCurrFemNodesPosition();
487 for (int i=0; i < ModalNodes; i++) {
488 pdBuffer->Put(i*3, pNodePos->GetVec(i+1));
489 }
490 pInterfComm->Send(pdBuffer->pdGetVec(), 3*ModalNodes, MPI::DOUBLE, 0,(this->GetLabel())*10+2);
491 }
492
~AerodynamicExternalModal(void)493 AerodynamicExternalModal::~AerodynamicExternalModal(void)
494 {
495 if (pdBuffer != NULL) SAFEDELETE(pdBuffer);
496 if (pdBufferVel != NULL) SAFEDELETE(pdBufferVel);
497 if (pForce != NULL) SAFEDELETE(pForce);
498 }
499
500 /* invia posizione e velocita' predetti */
501 void
AfterPredict(VectorHandler & X,VectorHandler & XP)502 AerodynamicExternalModal::AfterPredict(VectorHandler& X, VectorHandler& XP)
503 {
504 DEBUGCOUTFNAME("AerodynamicExternalModal::AfterPredict");
505 Send(X, XP);
506 }
507
508 /* invia posizione e velocita' predetti */
Update(const VectorHandler & X,const VectorHandler & XP)509 void AerodynamicExternalModal::Update(const VectorHandler& X ,
510 const VectorHandler& XP )
511 {
512 DEBUGCOUTFNAME("AerodynamicExternalModal::Update");
513 Send(X, XP);
514 }
515
Send(const VectorHandler & X,const VectorHandler & XP)516 void AerodynamicExternalModal::Send(const VectorHandler& X ,
517 const VectorHandler& XP )
518 {
519
520 Vec3 Vinf(0.);
521 doublereal rho, c, p, T;
522 GetAirProps(Vinf, rho, c, p, T);
523 /* current time */
524 pdBuffer->PutCoef(1, DriveOwner::dGet());
525 pdBuffer->Put(2, Vinf);
526 pdBuffer->PutCoef(5, rho);
527 pdBuffer->PutCoef(6,c);
528 pdBuffer->PutCoef(7,p);
529 pdBuffer->PutCoef(8,T);
530
531 Mat3xN* pNodePos = pModal->GetCurrFemNodesPosition();
532 Mat3xN* pNodeVel = NULL;
533 if (VelFlag) pNodeVel = pModal->GetCurrFemNodesVelocity();
534 for (int i=0; i < ModalNodes; i++) {
535 pdBuffer->Put(8+i*3, pNodePos->GetVec(i+1));
536 if (VelFlag) pdBufferVel->Put(i*3, pNodeVel->GetVec(i+1));
537 }
538 pInterfComm->Send(pdBuffer->pdGetVec(), 8+3*ModalNodes, MPI::DOUBLE, 0,(this->GetLabel())*10+2);
539 if (VelFlag) pInterfComm->Send(pdBufferVel->pdGetVec(), 3*ModalNodes, MPI::DOUBLE, 0,(this->GetLabel())*10+3);
540 SentFlag = true;
541
542 }
543
544
545 SubVectorHandler&
AssRes(SubVectorHandler & WorkVec,doublereal dCoef,const VectorHandler & XCurr,const VectorHandler & XPrimeCurr)546 AerodynamicExternalModal::AssRes(SubVectorHandler& WorkVec,
547 doublereal dCoef,
548 const VectorHandler& XCurr,
549 const VectorHandler& XPrimeCurr)
550 {
551 DEBUGCOUTFNAME("AerodynamicExternalModal::AssRes");
552 if (!SentFlag) {
553 Send(XCurr, XPrimeCurr);
554 }
555 integer NModes = pModal->uGetNModes();
556 WorkVec.ResizeReset(NModes);
557
558
559 integer iFirstIndex = pModal->iGetFirstIndex() + NModes;
560 for (int iCnt = 0; iCnt < NModes; iCnt++) {
561 WorkVec.PutRowIndex(iCnt, iFirstIndex + iCnt);
562 }
563
564 pInterfComm->Recv(pdBuffer->pdGetVec(), 3*ModalNodes, MPI::DOUBLE, 0, GetLabel()*10 + 4);
565 if (MomFlag) {
566 pInterfComm->Recv(pdBufferVel->pdGetVec(), 3*ModalNodes, MPI::DOUBLE, 0, GetLabel()*10 + 5);
567 }
568
569 /* calcola le forze su ciascun nodo
570 * le forze e i momenti sono espressi nel sistema di riferimento global
571 * vanno quindi riportati nel sistema di riferimento del nodo modal e
572 * poi proiettati sulla base modale
573 */
574
575 Vec3 F(0.);
576 Vec3 M(0.);
577 Mat3x3 RT(pModal->pGetModalNode()->GetRCurr().Transpose());
578
579 #warning "FIXME: review"
580 for (int iNode = 0; iNode < ModalNodes; iNode++) {
581 F[0] += pdBuffer->operator()(iNode*ModalNodes*3 + 1);
582 F[1] += pdBuffer->operator()(iNode*ModalNodes*3 + 2);
583 F[2] += pdBuffer->operator()(iNode*ModalNodes*3 + 3);
584 if (MomFlag) {
585 M[0] += pdBufferVel->operator()(iNode*ModalNodes*3 + 1);
586 M[1] += pdBufferVel->operator()(iNode*ModalNodes*3 + 2);
587 M[2] += pdBufferVel->operator()(iNode*ModalNodes*3 + 3);
588 }
589 for (int iMode = 0; iMode < NModes; iMode++) {
590 WorkVec.Add(iMode + 1, RT*F*(pModal->pGetPHIt().GetVec(iMode*ModalNodes + iNode + 1)));
591 if (MomFlag) {
592 WorkVec.Add(iMode + 1, RT*M*(pModal->pGetPHIr().GetVec(iMode*ModalNodes + iNode + 1)));
593 }
594 }
595 }
596
597 if (bToBeOutput()) {
598 *pForce = WorkVec;
599 }
600
601 SentFlag = false;
602 return WorkVec;
603 }
604
605
606 void
Output(OutputHandler & OH) const607 AerodynamicExternalModal::Output(OutputHandler& OH) const
608 {
609 if (bToBeOutput()) {
610 std::ostream& out = OH.Externals()
611 << std::setw(8) << -1*GetLabel();
612 int size = pForce->iGetSize();
613 for (int i = 0; i < size; i++) {
614 out << std::setw(8)
615 << " " << pForce->operator()(i+1) << std::endl;
616 }
617 }
618 }
619
620 Elem *
ReadAerodynamicExternalModal(DataManager * pDM,MBDynParser & HP,const DofOwner * pDO,unsigned int uLabel)621 ReadAerodynamicExternalModal(DataManager* pDM, MBDynParser& HP,
622 const DofOwner *pDO, unsigned int uLabel)
623 {
624 /* legge i dati d'ingresso e li passa al costruttore dell'elemento */
625 AerodynamicExternalModal* pEl = NULL;
626
627 /* giunto modale collegato */
628 const Modal *pModalJoint = pDM->ReadElem<const Modal, const Joint, Elem::JOINT>(HP);
629 int comm = 0;
630 if (HP.IsKeyWord("comm")) {
631 comm = HP.GetInt();
632 }
633
634 std::list<MPI::Intercomm>::iterator iComm = InterfaceComms.begin();
635 for (int i = 0; i++ < comm; ++iComm) {}
636 MPI::Intercomm* pInterC = &(*(iComm));
637
638 bool VelFlag = false;
639 bool MomFlag = false;
640
641 if (HP.IsKeyWord("velocity")) VelFlag = true;
642 if (HP.IsKeyWord("moment")) MomFlag = true;
643
644 flag fOut = pDM->fReadOutput(HP, Elem::AERODYNAMIC);
645 DriveCaller* pDC = NULL;
646 SAFENEWWITHCONSTRUCTOR(pDC, TimeDriveCaller, TimeDriveCaller(pDM->pGetDrvHdl()));
647
648 SAFENEWWITHCONSTRUCTOR(pEl,
649 AerodynamicExternalModal,
650 AerodynamicExternalModal(uLabel,
651 pDO,
652 pDC,
653 pModalJoint,
654 pInterC,
655 fOut,
656 VelFlag,
657 MomFlag));
658
659 return pEl;
660 }
661
662 /* Aerodynamic External Modal - end */
663
664 #endif /* USE_EXTERNAL */
665
666