1 // Copyright (C) 2007-2023 CEA, EDF, OPEN CASCADE
3 // Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
6 // This library is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License, or (at your option) any later version.
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 // Lesser General Public License for more details.
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
23 #include "MED_GaussUtils.hxx"
24 #include "MED_Utilities.hxx"
28 //---------------------------------------------------------------
31 TModeSwitchInfo(eFULL_INTERLACE),
41 ::Init(TInt theNbElem,
46 myModeSwitch = theMode;
49 myNbGauss = theNbGauss;
52 myGaussStep = myNbGauss*myDim;
54 myGaussCoord.resize(theNbElem*myGaussStep);
82 return (unsigned char*)&(myGaussCoord[0]);
87 ::GetCoordSliceArr(TInt theElemId) const
89 TCCoordSliceArr aCoordSliceArr(myNbGauss);
90 if(GetModeSwitch() == eFULL_INTERLACE){
91 TInt anId = theElemId*myGaussStep;
92 for(TInt anGaussId = 0; anGaussId < myNbGauss; anGaussId++){
93 aCoordSliceArr[anGaussId] =
94 TCCoordSlice(myGaussCoord,std::slice(anId,myDim,1));
99 for(TInt anGaussId = 0; anGaussId < myNbGauss; anGaussId++){
100 aCoordSliceArr[anGaussId] =
101 TCCoordSlice(myGaussCoord,std::slice(theElemId,myDim,myGaussStep));
104 return aCoordSliceArr;
109 ::GetCoordSliceArr(TInt theElemId)
111 TCoordSliceArr aCoordSliceArr(myNbGauss);
112 if(GetModeSwitch() == eFULL_INTERLACE){
113 TInt anId = theElemId*myGaussStep;
114 for(TInt anGaussId = 0; anGaussId < myNbGauss; anGaussId++){
115 aCoordSliceArr[anGaussId] =
116 TCoordSlice(myGaussCoord,std::slice(anId,myDim,1));
121 for(TInt anGaussId = 0; anGaussId < myNbGauss; anGaussId++){
122 aCoordSliceArr[anGaussId] =
123 TCoordSlice(myGaussCoord,std::slice(theElemId,myDim,myGaussStep));
126 return aCoordSliceArr;
129 //---------------------------------------------------------------
132 IsEqual(TFloat theLeft, TFloat theRight)
134 static TFloat EPS = 1.0E-3;
135 if(fabs(theLeft) + fabs(theRight) > EPS)
136 return fabs(theLeft-theRight)/(fabs(theLeft)+fabs(theRight)) < EPS;
140 //---------------------------------------------------------------
141 class TShapeFun::TFun
149 Init(TInt theNbGauss,
152 myFun.resize(theNbGauss*theNbRef);
157 GetFunSlice(TInt theGaussId) const
159 return TCFloatVecSlice(myFun,std::slice(theGaussId*myNbRef,myNbRef,1));
163 GetFunSlice(TInt theGaussId)
165 return TFloatVecSlice(myFun,std::slice(theGaussId*myNbRef,myNbRef,1));
169 //---------------------------------------------------------------
171 TShapeFun::TShapeFun(TInt theDim, TInt theNbRef):
172 myRefCoord(theNbRef*theDim),
178 TShapeFun::GetCoord(TInt theRefId) const
180 return TCCoordSlice(myRefCoord,std::slice(theRefId*myDim,myDim,1));
184 TShapeFun::GetCoord(TInt theRefId)
186 return TCoordSlice(myRefCoord,std::slice(theRefId*myDim,myDim,1));
190 TShapeFun::GetFun(const TCCoordSliceArr& theRef,
191 const TCCoordSliceArr& theGauss,
194 TInt aNbRef = theRef.size();
195 TInt aNbGauss = theGauss.size();
196 theFun.Init(aNbGauss,aNbRef);
200 TShapeFun::IsSatisfy(const TCCoordSliceArr& theRefCoord) const
202 TInt aNbRef = theRefCoord.size();
203 TInt aNbRef2 = GetNbRef();
204 INITMSG("TShapeFun::IsSatisfy "<<
205 "- aNbRef("<<aNbRef<<")"<<
206 "; aNbRef2("<<aNbRef2<<")\n");
207 bool anIsSatisfy = (aNbRef == aNbRef2);
209 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
210 const TCCoordSlice& aCoord2 = theRefCoord[aRefId];
211 TCCoordSlice aCoord = GetCoord(aRefId);
212 TInt aDim = aCoord.size();
213 bool anIsEqual = false;
214 for(TInt anId = 0; anId < aDim; anId++){
215 anIsEqual = IsEqual(aCoord[anId],aCoord2[anId]);
222 if(SALOME::VerbosityActivated()){
223 TCCoordSlice aCoord = GetCoord(aRefId);
224 INITMSG(aRefId + 1<<": aCoord = {");
225 TInt aDim = aCoord.size();
226 for(TInt anId = 0; anId < aDim; anId++)
227 ADDMSG("\t"<<aCoord[anId]);
228 const TCCoordSlice& aCoord2 = theRefCoord[aRefId];
229 ADDMSG("}\t!=\taCoord2 = {");
230 for(TInt anId = 0; anId < aDim; anId++)
231 ADDMSG("\t"<<aCoord2[anId]);
235 BEGMSG("anIsSatisfy = "<<anIsSatisfy<<"\n");
241 BEGMSG("anIsSatisfy = "<<anIsSatisfy<<"\n");
246 TShapeFun::Eval(const TCellInfo& theCellInfo,
247 const TNodeInfo& theNodeInfo,
248 const TElemNum& theElemNum,
249 const TCCoordSliceArr& theRef,
250 const TCCoordSliceArr& theGauss,
251 TGaussCoord& theGaussCoord,
254 INITMSG("TShapeFun::Eval"<<std::endl);
256 if(IsSatisfy(theRef)){
257 const PMeshInfo& aMeshInfo = theCellInfo.GetMeshInfo();
258 TInt aDim = aMeshInfo->GetDim();
259 TInt aNbGauss = theGauss.size();
261 bool anIsSubMesh = !theElemNum.empty();
264 aNbElem = theElemNum.size();
266 aNbElem = theCellInfo.GetNbElem();
268 theGaussCoord.Init(aNbElem,aNbGauss,aDim,theMode);
271 InitFun(theRef,theGauss,aFun);
272 TInt aConnDim = theCellInfo.GetConnDim();
274 INITMSG("aDim = "<<aDim<<
275 "; aNbGauss = "<<aNbGauss<<
276 "; aNbElem = "<<aNbElem<<
277 "; aNbNodes = "<<theNodeInfo.GetNbElem()<<
280 for(TInt anElemId = 0; anElemId < aNbElem; anElemId++){
281 TInt aCellId = anIsSubMesh? theElemNum[anElemId]-1: anElemId;
282 TCConnSlice aConnSlice = theCellInfo.GetConnSlice(aCellId);
283 TCoordSliceArr aCoordSliceArr = theGaussCoord.GetCoordSliceArr(anElemId);
285 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
286 TCoordSlice& aGaussCoordSlice = aCoordSliceArr[aGaussId];
287 TCFloatVecSlice aFunSlice = aFun.GetFunSlice(aGaussId);
289 for(TInt aConnId = 0; aConnId < aConnDim; aConnId++){
290 TInt aNodeId = aConnSlice[aConnId] - 1;
291 TCCoordSlice aNodeCoordSlice = theNodeInfo.GetCoordSlice(aNodeId);
293 for(TInt aDimId = 0; aDimId < aDim; aDimId++){
294 aGaussCoordSlice[aDimId] += aNodeCoordSlice[aDimId]*aFunSlice[aConnId];
300 if (SALOME::VerbosityActivated())
302 INITMSG("theGauss: ");
303 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
304 TCCoordSlice aCoordSlice = theGauss[aGaussId];
306 for(TInt aDimId = 0; aDimId < aDim; aDimId++){
307 ADDMSG(aCoordSlice[aDimId]<<" ");
313 for(TInt anElemId = 0; anElemId < aNbElem; anElemId++){
314 TCCoordSliceArr aCoordSliceArr = theGaussCoord.GetCoordSliceArr(anElemId);
316 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
317 TCCoordSlice aCoordSlice = aCoordSliceArr[aGaussId];
319 for(TInt aDimId = 0; aDimId < aDim; aDimId++){
320 ADDMSG(aCoordSlice[aDimId]<<" ");
334 //---------------------------------------------------------------
335 TSeg2a::TSeg2a():TShapeFun(1,2)
337 TInt aNbRef = GetNbRef();
338 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
339 TCoordSlice aCoord = GetCoord(aRefId);
341 case 0: aCoord[0] = -1.0; break;
342 case 1: aCoord[0] = 1.0; break;
348 TSeg2a::InitFun(const TCCoordSliceArr& theRef,
349 const TCCoordSliceArr& theGauss,
352 GetFun(theRef,theGauss,theFun);
354 TInt aNbGauss = theGauss.size();
355 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
356 const TCCoordSlice& aCoord = theGauss[aGaussId];
357 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
359 aSlice[0] = 0.5*(1.0 - aCoord[0]);
360 aSlice[1] = 0.5*(1.0 + aCoord[0]);
364 //---------------------------------------------------------------
365 TSeg3a::TSeg3a():TShapeFun(1,3)
367 TInt aNbRef = GetNbRef();
368 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
369 TCoordSlice aCoord = GetCoord(aRefId);
371 case 0: aCoord[0] = -1.0; break;
372 case 1: aCoord[0] = 1.0; break;
373 case 2: aCoord[0] = 0.0; break;
379 TSeg3a::InitFun(const TCCoordSliceArr& theRef,
380 const TCCoordSliceArr& theGauss,
383 GetFun(theRef,theGauss,theFun);
385 TInt aNbGauss = theGauss.size();
386 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
387 const TCCoordSlice& aCoord = theGauss[aGaussId];
388 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
390 aSlice[0] = 0.5*(1.0 - aCoord[0])*aCoord[0];
391 aSlice[1] = 0.5*(1.0 + aCoord[0])*aCoord[0];
392 aSlice[2] = (1.0 + aCoord[0])*(1.0 - aCoord[0]);
396 //---------------------------------------------------------------
400 TInt aNbRef = GetNbRef();
401 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
402 TCoordSlice aCoord = GetCoord(aRefId);
404 case 0: aCoord[0] = -1.0; aCoord[1] = 1.0; break;
405 case 1: aCoord[0] = -1.0; aCoord[1] = -1.0; break;
406 case 2: aCoord[0] = 1.0; aCoord[1] = -1.0; break;
412 TTria3a::InitFun(const TCCoordSliceArr& theRef,
413 const TCCoordSliceArr& theGauss,
416 GetFun(theRef,theGauss,theFun);
418 TInt aNbGauss = theGauss.size();
419 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
420 const TCCoordSlice& aCoord = theGauss[aGaussId];
421 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
423 aSlice[0] = 0.5*(1.0 + aCoord[1]);
424 aSlice[1] = -0.5*(aCoord[0] + aCoord[1]);
425 aSlice[2] = 0.5*(1.0 + aCoord[0]);
429 //---------------------------------------------------------------
430 TTria6a::TTria6a():TShapeFun(2,6)
432 TInt aNbRef = GetNbRef();
433 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
434 TCoordSlice aCoord = GetCoord(aRefId);
436 case 0: aCoord[0] = -1.0; aCoord[1] = 1.0; break;
437 case 1: aCoord[0] = -1.0; aCoord[1] = -1.0; break;
438 case 2: aCoord[0] = 1.0; aCoord[1] = -1.0; break;
440 case 3: aCoord[0] = -1.0; aCoord[1] = 1.0; break;
441 case 4: aCoord[0] = 0.0; aCoord[1] = -1.0; break;
442 case 5: aCoord[0] = 0.0; aCoord[1] = 0.0; break;
448 TTria6a::InitFun(const TCCoordSliceArr& theRef,
449 const TCCoordSliceArr& theGauss,
452 GetFun(theRef,theGauss,theFun);
454 TInt aNbGauss = theGauss.size();
455 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
456 const TCCoordSlice& aCoord = theGauss[aGaussId];
457 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
459 aSlice[0] = 0.5*(1.0 + aCoord[1])*aCoord[1];
460 aSlice[1] = 0.5*(aCoord[0] + aCoord[1])*(aCoord[0] + aCoord[1] + 1);
461 aSlice[2] = 0.5*(1.0 + aCoord[0])*aCoord[0];
463 aSlice[3] = -1.0*(1.0 + aCoord[1])*(aCoord[0] + aCoord[1]);
464 aSlice[4] = -1.0*(1.0 + aCoord[0])*(aCoord[0] + aCoord[1]);
465 aSlice[5] = (1.0 + aCoord[1])*(1.0 + aCoord[1]);
469 //---------------------------------------------------------------
473 TInt aNbRef = GetNbRef();
474 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
475 TCoordSlice aCoord = GetCoord(aRefId);
477 case 0: aCoord[0] = 0.0; aCoord[1] = 0.0; break;
478 case 1: aCoord[0] = 1.0; aCoord[1] = 0.0; break;
479 case 2: aCoord[0] = 0.0; aCoord[1] = 1.0; break;
485 TTria3b::InitFun(const TCCoordSliceArr& theRef,
486 const TCCoordSliceArr& theGauss,
489 GetFun(theRef,theGauss,theFun);
491 TInt aNbGauss = theGauss.size();
492 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
493 const TCCoordSlice& aCoord = theGauss[aGaussId];
494 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
496 aSlice[0] = 1.0 - aCoord[0] - aCoord[1];
497 aSlice[1] = aCoord[0];
498 aSlice[2] = aCoord[1];
502 //---------------------------------------------------------------
506 TInt aNbRef = GetNbRef();
507 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
508 TCoordSlice aCoord = GetCoord(aRefId);
510 case 0: aCoord[0] = 0.0; aCoord[1] = 0.0; break;
511 case 1: aCoord[0] = 1.0; aCoord[1] = 0.0; break;
512 case 2: aCoord[0] = 0.0; aCoord[1] = 1.0; break;
514 case 3: aCoord[0] = 0.5; aCoord[1] = 0.0; break;
515 case 4: aCoord[0] = 0.5; aCoord[1] = 0.5; break;
516 case 5: aCoord[0] = 0.0; aCoord[1] = 0.5; break;
522 TTria6b::InitFun(const TCCoordSliceArr& theRef,
523 const TCCoordSliceArr& theGauss,
526 GetFun(theRef,theGauss,theFun);
528 TInt aNbGauss = theGauss.size();
529 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
530 const TCCoordSlice& aCoord = theGauss[aGaussId];
531 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
533 aSlice[0] = (1.0 - aCoord[0] - aCoord[1])*(1.0 - 2.0*aCoord[0] - 2.0*aCoord[1]);
534 aSlice[1] = aCoord[0]*(2.0*aCoord[0] - 1.0);
535 aSlice[2] = aCoord[1]*(2.0*aCoord[1] - 1.0);
537 aSlice[3] = 4.0*aCoord[0]*(1.0 - aCoord[0] - aCoord[1]);
538 aSlice[4] = 4.0*aCoord[0]*aCoord[1];
539 aSlice[5] = 4.0*aCoord[1]*(1.0 - aCoord[0] - aCoord[1]);
543 //---------------------------------------------------------------
547 TInt aNbRef = GetNbRef();
548 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
549 TCoordSlice aCoord = GetCoord(aRefId);
551 case 0: aCoord[0] = -1.0; aCoord[1] = 1.0; break;
552 case 1: aCoord[0] = -1.0; aCoord[1] = -1.0; break;
553 case 2: aCoord[0] = 1.0; aCoord[1] = -1.0; break;
554 case 3: aCoord[0] = 1.0; aCoord[1] = 1.0; break;
560 TQuad4a::InitFun(const TCCoordSliceArr& theRef,
561 const TCCoordSliceArr& theGauss,
564 GetFun(theRef,theGauss,theFun);
566 TInt aNbGauss = theGauss.size();
567 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
568 const TCCoordSlice& aCoord = theGauss[aGaussId];
569 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
571 aSlice[0] = 0.25*(1.0 + aCoord[1])*(1.0 - aCoord[0]);
572 aSlice[1] = 0.25*(1.0 - aCoord[1])*(1.0 - aCoord[0]);
573 aSlice[2] = 0.25*(1.0 - aCoord[1])*(1.0 + aCoord[0]);
574 aSlice[3] = 0.25*(1.0 + aCoord[0])*(1.0 + aCoord[1]);
578 //---------------------------------------------------------------
582 TInt aNbRef = GetNbRef();
583 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
584 TCoordSlice aCoord = GetCoord(aRefId);
586 case 0: aCoord[0] = -1.0; aCoord[1] = 1.0; break;
587 case 1: aCoord[0] = -1.0; aCoord[1] = -1.0; break;
588 case 2: aCoord[0] = 1.0; aCoord[1] = -1.0; break;
589 case 3: aCoord[0] = 1.0; aCoord[1] = 1.0; break;
591 case 4: aCoord[0] = -1.0; aCoord[1] = 0.0; break;
592 case 5: aCoord[0] = 0.0; aCoord[1] = -1.0; break;
593 case 6: aCoord[0] = 1.0; aCoord[1] = 0.0; break;
594 case 7: aCoord[0] = 0.0; aCoord[1] = 1.0; break;
600 TQuad8a::InitFun(const TCCoordSliceArr& theRef,
601 const TCCoordSliceArr& theGauss,
604 GetFun(theRef,theGauss,theFun);
606 TInt aNbGauss = theGauss.size();
607 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
608 const TCCoordSlice& aCoord = theGauss[aGaussId];
609 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
611 aSlice[0] = 0.25*(1.0 + aCoord[1])*(1.0 - aCoord[0])*(aCoord[1] - aCoord[0] - 1.0);
612 aSlice[1] = 0.25*(1.0 - aCoord[1])*(1.0 - aCoord[0])*(-aCoord[1] - aCoord[0] - 1.0);
613 aSlice[2] = 0.25*(1.0 - aCoord[1])*(1.0 + aCoord[0])*(-aCoord[1] + aCoord[0] - 1.0);
614 aSlice[3] = 0.25*(1.0 + aCoord[1])*(1.0 + aCoord[0])*(aCoord[1] + aCoord[0] - 1.0);
616 aSlice[4] = 0.5*(1.0 - aCoord[0])*(1.0 - aCoord[1])*(1.0 + aCoord[1]);
617 aSlice[5] = 0.5*(1.0 - aCoord[1])*(1.0 - aCoord[0])*(1.0 + aCoord[0]);
618 aSlice[6] = 0.5*(1.0 + aCoord[0])*(1.0 - aCoord[1])*(1.0 + aCoord[1]);
619 aSlice[7] = 0.5*(1.0 + aCoord[1])*(1.0 - aCoord[0])*(1.0 + aCoord[0]);
623 //---------------------------------------------------------------
627 TInt aNbRef = GetNbRef();
628 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
629 TCoordSlice aCoord = GetCoord(aRefId);
631 case 0: aCoord[0] = -1.0; aCoord[1] = 1.0; break;
632 case 1: aCoord[0] = -1.0; aCoord[1] = -1.0; break;
633 case 2: aCoord[0] = 1.0; aCoord[1] = -1.0; break;
634 case 3: aCoord[0] = 1.0; aCoord[1] = 1.0; break;
636 case 4: aCoord[0] = -1.0; aCoord[1] = 0.0; break;
637 case 5: aCoord[0] = 0.0; aCoord[1] = -1.0; break;
638 case 6: aCoord[0] = 1.0; aCoord[1] = 0.0; break;
639 case 7: aCoord[0] = 0.0; aCoord[1] = 1.0; break;
641 case 8: aCoord[0] = 0.0; aCoord[1] = 0.0; break;
647 TQuad9a::InitFun(const TCCoordSliceArr& theRef,
648 const TCCoordSliceArr& theGauss,
651 GetFun(theRef,theGauss,theFun);
653 TInt aNbGauss = theGauss.size();
654 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
655 const TCCoordSlice& aCoord = theGauss[aGaussId];
656 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
658 aSlice[0] = 0.25*aCoord[0]*aCoord[1]*(aCoord[0] + 1.0)*(aCoord[1] - 1.0);
659 aSlice[1] = 0.25*aCoord[0]*aCoord[1]*(aCoord[0] + 1.0)*(aCoord[1] + 1.0);
660 aSlice[2] = 0.25*aCoord[0]*aCoord[1]*(aCoord[0] - 1.0)*(aCoord[1] + 1.0);
661 aSlice[3] = 0.25*aCoord[0]*aCoord[1]*(aCoord[0] - 1.0)*(aCoord[1] - 1.0);
663 aSlice[4] = 0.5*aCoord[0]*(aCoord[0] + 1.0)*(1.0 - aCoord[1]*aCoord[1]);
664 aSlice[5] = 0.5*(1.0 - aCoord[0]*aCoord[0])*aCoord[1]*(aCoord[1] + 1.0);
665 aSlice[6] = 0.5*aCoord[0]*(aCoord[0] - 1.0)*(1.0 - aCoord[1]*aCoord[1]);
666 aSlice[7] = 0.5*(1.0 - aCoord[0]*aCoord[0])*aCoord[1]*(aCoord[1] - 1.0);
668 aSlice[8] = (1.0 - aCoord[0]*aCoord[0])*(1.0 - aCoord[1]*aCoord[1]);
672 //---------------------------------------------------------------
676 TInt aNbRef = GetNbRef();
677 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
678 TCoordSlice aCoord = GetCoord(aRefId);
680 case 0: aCoord[0] = -1.0; aCoord[1] = -1.0; break;
681 case 1: aCoord[0] = 1.0; aCoord[1] = -1.0; break;
682 case 2: aCoord[0] = 1.0; aCoord[1] = 1.0; break;
683 case 3: aCoord[0] = -1.0; aCoord[1] = 1.0; break;
689 TQuad4b::InitFun(const TCCoordSliceArr& theRef,
690 const TCCoordSliceArr& theGauss,
693 GetFun(theRef,theGauss,theFun);
695 TInt aNbGauss = theGauss.size();
696 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
697 const TCCoordSlice& aCoord = theGauss[aGaussId];
698 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
700 aSlice[0] = 0.25*(1.0 - aCoord[0])*(1.0 - aCoord[1]);
701 aSlice[1] = 0.25*(1.0 + aCoord[0])*(1.0 - aCoord[1]);
702 aSlice[2] = 0.25*(1.0 + aCoord[0])*(1.0 + aCoord[1]);
703 aSlice[3] = 0.25*(1.0 - aCoord[0])*(1.0 + aCoord[1]);
707 //---------------------------------------------------------------
711 TInt aNbRef = GetNbRef();
712 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
713 TCoordSlice aCoord = GetCoord(aRefId);
715 case 0: aCoord[0] = -1.0; aCoord[1] = -1.0; break;
716 case 1: aCoord[0] = 1.0; aCoord[1] = -1.0; break;
717 case 2: aCoord[0] = 1.0; aCoord[1] = 1.0; break;
718 case 3: aCoord[0] = -1.0; aCoord[1] = 1.0; break;
720 case 4: aCoord[0] = 0.0; aCoord[1] = -1.0; break;
721 case 5: aCoord[0] = 1.0; aCoord[1] = 0.0; break;
722 case 6: aCoord[0] = 0.0; aCoord[1] = 1.0; break;
723 case 7: aCoord[0] = -1.0; aCoord[1] = 0.0; break;
729 TQuad8b::InitFun(const TCCoordSliceArr& theRef,
730 const TCCoordSliceArr& theGauss,
733 GetFun(theRef,theGauss,theFun);
735 TInt aNbGauss = theGauss.size();
736 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
737 const TCCoordSlice& aCoord = theGauss[aGaussId];
738 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
740 aSlice[0] = 0.25*(1.0 - aCoord[0])*(1.0 - aCoord[1])*(-1.0 - aCoord[0] - aCoord[1]);
741 aSlice[1] = 0.25*(1.0 + aCoord[0])*(1.0 - aCoord[1])*(-1.0 + aCoord[0] - aCoord[1]);
742 aSlice[2] = 0.25*(1.0 + aCoord[0])*(1.0 + aCoord[1])*(-1.0 + aCoord[0] + aCoord[1]);
743 aSlice[3] = 0.25*(1.0 - aCoord[0])*(1.0 + aCoord[1])*(-1.0 - aCoord[0] + aCoord[1]);
745 aSlice[4] = 0.5*(1.0 - aCoord[0]*aCoord[0])*(1.0 - aCoord[1]);
746 aSlice[5] = 0.5*(1.0 - aCoord[1]*aCoord[1])*(1.0 + aCoord[0]);
747 aSlice[6] = 0.5*(1.0 - aCoord[0]*aCoord[0])*(1.0 + aCoord[1]);
748 aSlice[7] = 0.5*(1.0 - aCoord[1]*aCoord[1])*(1.0 - aCoord[0]);
750 //aSlice[4] = 0.5*(1.0 - aCoord[0])*(1.0 - aCoord[0])*(1.0 - aCoord[1]);
751 //aSlice[5] = 0.5*(1.0 + aCoord[0])*(1.0 - aCoord[1])*(1.0 - aCoord[1]);
752 //aSlice[6] = 0.5*(1.0 - aCoord[0])*(1.0 - aCoord[0])*(1.0 + aCoord[1]);
753 //aSlice[7] = 0.5*(1.0 - aCoord[0])*(1.0 - aCoord[1])*(1.0 - aCoord[1]);
757 //---------------------------------------------------------------
761 TInt aNbRef = GetNbRef();
762 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
763 TCoordSlice aCoord = GetCoord(aRefId);
765 case 0: aCoord[0] = -1.0; aCoord[1] = -1.0; break;
766 case 1: aCoord[0] = 1.0; aCoord[1] = -1.0; break;
767 case 2: aCoord[0] = 1.0; aCoord[1] = 1.0; break;
768 case 3: aCoord[0] = -1.0; aCoord[1] = 1.0; break;
770 case 4: aCoord[0] = 0.0; aCoord[1] = -1.0; break;
771 case 5: aCoord[0] = 1.0; aCoord[1] = 0.0; break;
772 case 6: aCoord[0] = 0.0; aCoord[1] = 1.0; break;
773 case 7: aCoord[0] = -1.0; aCoord[1] = 0.0; break;
775 case 8: aCoord[0] = 0.0; aCoord[1] = 0.0; break;
781 TQuad9b::InitFun(const TCCoordSliceArr& theRef,
782 const TCCoordSliceArr& theGauss,
785 GetFun(theRef,theGauss,theFun);
787 TInt aNbGauss = theGauss.size();
788 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
789 const TCCoordSlice& aCoord = theGauss[aGaussId];
790 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
792 aSlice[0] = 0.25*aCoord[0]*aCoord[1]*(aCoord[0] - 1.0)*(aCoord[1] - 1.0);
793 aSlice[1] = 0.25*aCoord[0]*aCoord[1]*(aCoord[0] + 1.0)*(aCoord[1] - 1.0);
794 aSlice[2] = 0.25*aCoord[0]*aCoord[1]*(aCoord[0] + 1.0)*(aCoord[1] + 1.0);
795 aSlice[3] = 0.25*aCoord[0]*aCoord[1]*(aCoord[0] - 1.0)*(aCoord[1] + 1.0);
797 aSlice[4] = 0.5*(1.0 - aCoord[0]*aCoord[0])*aCoord[1]*(aCoord[1] - 1.0);
798 aSlice[5] = 0.5*aCoord[0]*(aCoord[0] + 1.0)*(1.0 - aCoord[1]*aCoord[1]);
799 aSlice[6] = 0.5*(1.0 - aCoord[0]*aCoord[0])*aCoord[1]*(aCoord[1] + 1.0);
800 aSlice[7] = 0.5*aCoord[0]*(aCoord[0] - 1.0)*(1.0 - aCoord[1]*aCoord[1]);
802 aSlice[8] = (1.0 - aCoord[0]*aCoord[0])*(1.0 - aCoord[1]*aCoord[1]);
806 //---------------------------------------------------------------
807 TTetra4a::TTetra4a():
810 TInt aNbRef = GetNbRef();
811 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
812 TCoordSlice aCoord = GetCoord(aRefId);
814 case 0: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
815 case 1: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
816 case 2: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
817 case 3: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
823 TTetra4a::InitFun(const TCCoordSliceArr& theRef,
824 const TCCoordSliceArr& theGauss,
827 GetFun(theRef,theGauss,theFun);
829 TInt aNbGauss = theGauss.size();
830 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
831 const TCCoordSlice& aCoord = theGauss[aGaussId];
832 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
834 aSlice[0] = aCoord[1];
835 aSlice[1] = aCoord[2];
836 aSlice[2] = 1.0 - aCoord[0] - aCoord[1] - aCoord[2];
837 aSlice[3] = aCoord[0];
841 //---------------------------------------------------------------
842 TTetra10a::TTetra10a():
845 TInt aNbRef = GetNbRef();
846 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
847 TCoordSlice aCoord = GetCoord(aRefId);
849 case 0: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
850 case 1: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
851 case 2: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
852 case 3: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
854 case 4: aCoord[0] = 0.0; aCoord[1] = 0.5; aCoord[2] = 0.5; break;
855 case 5: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 0.5; break;
856 case 6: aCoord[0] = 0.0; aCoord[1] = 0.5; aCoord[2] = 0.0; break;
858 case 7: aCoord[0] = 0.5; aCoord[1] = 0.5; aCoord[2] = 0.0; break;
859 case 8: aCoord[0] = 0.5; aCoord[1] = 0.0; aCoord[2] = 0.5; break;
860 case 9: aCoord[0] = 0.5; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
866 TTetra10a::InitFun(const TCCoordSliceArr& theRef,
867 const TCCoordSliceArr& theGauss,
870 GetFun(theRef,theGauss,theFun);
872 TInt aNbGauss = theGauss.size();
873 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
874 const TCCoordSlice& aCoord = theGauss[aGaussId];
875 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
877 aSlice[0] = aCoord[1]*(2.0*aCoord[1] - 1.0);
878 aSlice[1] = aCoord[2]*(2.0*aCoord[2] - 1.0);
879 aSlice[2] = (1.0 - aCoord[0] - aCoord[1] - aCoord[2])*(1.0 - 2.0*aCoord[0] - 2.0*aCoord[1] - 2.0*aCoord[2]);
880 aSlice[3] = aCoord[0]*(2.0*aCoord[0] - 1.0);
882 aSlice[4] = 4.0*aCoord[1]*aCoord[2];
883 aSlice[5] = 4.0*aCoord[2]*(1.0 - aCoord[0] - aCoord[1] - aCoord[2]);
884 aSlice[6] = 4.0*aCoord[1]*(1.0 - aCoord[0] - aCoord[1] - aCoord[2]);
886 aSlice[7] = 4.0*aCoord[0]*aCoord[1];
887 aSlice[8] = 4.0*aCoord[0]*aCoord[2];
888 aSlice[9] = 4.0*aCoord[0]*(1.0 - aCoord[0] - aCoord[1] - aCoord[2]);
892 //---------------------------------------------------------------
894 TTetra4b::TTetra4b():
897 TInt aNbRef = GetNbRef();
898 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
899 TCoordSlice aCoord = GetCoord(aRefId);
901 case 0: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
902 case 2: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
903 case 1: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
904 case 3: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
910 TTetra4b::InitFun(const TCCoordSliceArr& theRef,
911 const TCCoordSliceArr& theGauss,
914 GetFun(theRef,theGauss,theFun);
916 TInt aNbGauss = theGauss.size();
917 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
918 const TCCoordSlice& aCoord = theGauss[aGaussId];
919 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
921 aSlice[0] = aCoord[1];
922 aSlice[2] = aCoord[2];
923 aSlice[1] = 1.0 - aCoord[0] - aCoord[1] - aCoord[2];
924 aSlice[3] = aCoord[0];
928 //---------------------------------------------------------------
929 TTetra10b::TTetra10b():
932 TInt aNbRef = GetNbRef();
933 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
934 TCoordSlice aCoord = GetCoord(aRefId);
936 case 0: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
937 case 2: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
938 case 1: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
939 case 3: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
941 case 6: aCoord[0] = 0.0; aCoord[1] = 0.5; aCoord[2] = 0.5; break;
942 case 5: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 0.5; break;
943 case 4: aCoord[0] = 0.0; aCoord[1] = 0.5; aCoord[2] = 0.0; break;
945 case 7: aCoord[0] = 0.5; aCoord[1] = 0.5; aCoord[2] = 0.0; break;
946 case 9: aCoord[0] = 0.5; aCoord[1] = 0.0; aCoord[2] = 0.5; break;
947 case 8: aCoord[0] = 0.5; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
953 TTetra10b::InitFun(const TCCoordSliceArr& theRef,
954 const TCCoordSliceArr& theGauss,
957 GetFun(theRef,theGauss,theFun);
959 TInt aNbGauss = theGauss.size();
960 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
961 const TCCoordSlice& aCoord = theGauss[aGaussId];
962 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
964 aSlice[0] = aCoord[1]*(2.0*aCoord[1] - 1.0);
965 aSlice[2] = aCoord[2]*(2.0*aCoord[2] - 1.0);
966 aSlice[1] = (1.0 - aCoord[0] - aCoord[1] - aCoord[2])*(1.0 - 2.0*aCoord[0] - 2.0*aCoord[1] - 2.0*aCoord[2]);
967 aSlice[3] = aCoord[0]*(2.0*aCoord[0] - 1.0);
969 aSlice[6] = 4.0*aCoord[1]*aCoord[2];
970 aSlice[5] = 4.0*aCoord[2]*(1.0 - aCoord[0] - aCoord[1] - aCoord[2]);
971 aSlice[4] = 4.0*aCoord[1]*(1.0 - aCoord[0] - aCoord[1] - aCoord[2]);
973 aSlice[7] = 4.0*aCoord[0]*aCoord[1];
974 aSlice[9] = 4.0*aCoord[0]*aCoord[2];
975 aSlice[8] = 4.0*aCoord[0]*(1.0 - aCoord[0] - aCoord[1] - aCoord[2]);
979 //---------------------------------------------------------------
983 TInt aNbRef = GetNbRef();
984 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
985 TCoordSlice aCoord = GetCoord(aRefId);
987 case 0: aCoord[0] = -1.0; aCoord[1] = -1.0; aCoord[2] = -1.0; break;
988 case 1: aCoord[0] = 1.0; aCoord[1] = -1.0; aCoord[2] = -1.0; break;
989 case 2: aCoord[0] = 1.0; aCoord[1] = 1.0; aCoord[2] = -1.0; break;
990 case 3: aCoord[0] = -1.0; aCoord[1] = 1.0; aCoord[2] = -1.0; break;
991 case 4: aCoord[0] = -1.0; aCoord[1] = -1.0; aCoord[2] = 1.0; break;
992 case 5: aCoord[0] = 1.0; aCoord[1] = -1.0; aCoord[2] = 1.0; break;
993 case 6: aCoord[0] = 1.0; aCoord[1] = 1.0; aCoord[2] = 1.0; break;
994 case 7: aCoord[0] = -1.0; aCoord[1] = 1.0; aCoord[2] = 1.0; break;
1000 THexa8a::InitFun(const TCCoordSliceArr& theRef,
1001 const TCCoordSliceArr& theGauss,
1004 GetFun(theRef,theGauss,theFun);
1006 TInt aNbGauss = theGauss.size();
1007 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1008 const TCCoordSlice& aCoord = theGauss[aGaussId];
1009 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
1011 aSlice[0] = 0.125*(1.0 - aCoord[0])*(1.0 - aCoord[1])*(1.0 - aCoord[2]);
1012 aSlice[1] = 0.125*(1.0 + aCoord[0])*(1.0 - aCoord[1])*(1.0 - aCoord[2]);
1013 aSlice[2] = 0.125*(1.0 + aCoord[0])*(1.0 + aCoord[1])*(1.0 - aCoord[2]);
1014 aSlice[3] = 0.125*(1.0 - aCoord[0])*(1.0 + aCoord[1])*(1.0 - aCoord[2]);
1016 aSlice[4] = 0.125*(1.0 - aCoord[0])*(1.0 - aCoord[1])*(1.0 + aCoord[2]);
1017 aSlice[5] = 0.125*(1.0 + aCoord[0])*(1.0 - aCoord[1])*(1.0 + aCoord[2]);
1018 aSlice[6] = 0.125*(1.0 + aCoord[0])*(1.0 + aCoord[1])*(1.0 + aCoord[2]);
1019 aSlice[7] = 0.125*(1.0 - aCoord[0])*(1.0 + aCoord[1])*(1.0 + aCoord[2]);
1023 //---------------------------------------------------------------
1024 THexa20a::THexa20a(TInt theDim, TInt theNbRef):
1025 TShapeFun(theDim,theNbRef)
1027 TInt aNbRef = myRefCoord.size();
1028 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
1029 TCoordSlice aCoord = GetCoord(aRefId);
1031 case 0: aCoord[0] = -1.0; aCoord[1] = -1.0; aCoord[2] = -1.0; break;
1032 case 1: aCoord[0] = 1.0; aCoord[1] = -1.0; aCoord[2] = -1.0; break;
1033 case 2: aCoord[0] = 1.0; aCoord[1] = 1.0; aCoord[2] = -1.0; break;
1034 case 3: aCoord[0] = -1.0; aCoord[1] = 1.0; aCoord[2] = -1.0; break;
1035 case 4: aCoord[0] = -1.0; aCoord[1] = -1.0; aCoord[2] = 1.0; break;
1036 case 5: aCoord[0] = 1.0; aCoord[1] = -1.0; aCoord[2] = 1.0; break;
1037 case 6: aCoord[0] = 1.0; aCoord[1] = 1.0; aCoord[2] = 1.0; break;
1038 case 7: aCoord[0] = -1.0; aCoord[1] = 1.0; aCoord[2] = 1.0; break;
1040 case 8: aCoord[0] = 0.0; aCoord[1] = -1.0; aCoord[2] = -1.0; break;
1041 case 9: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = -1.0; break;
1042 case 10: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = -1.0; break;
1043 case 11: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = -1.0; break;
1044 case 12: aCoord[0] = -1.0; aCoord[1] = -1.0; aCoord[2] = 0.0; break;
1045 case 13: aCoord[0] = 1.0; aCoord[1] = -1.0; aCoord[2] = 0.0; break;
1046 case 14: aCoord[0] = 1.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1047 case 15: aCoord[0] = -1.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1048 case 16: aCoord[0] = 0.0; aCoord[1] = -1.0; aCoord[2] = 1.0; break;
1049 case 17: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1050 case 18: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = 1.0; break;
1051 case 19: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1057 THexa20a::InitFun(const TCCoordSliceArr& theRef,
1058 const TCCoordSliceArr& theGauss,
1061 GetFun(theRef,theGauss,theFun);
1063 TInt aNbGauss = theGauss.size();
1064 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1065 const TCCoordSlice& aCoord = theGauss[aGaussId];
1066 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
1068 aSlice[0] = 0.125*(1.0 - aCoord[0])*(1.0 - aCoord[1])*(1.0 - aCoord[2])*
1069 (-2.0 - aCoord[0] - aCoord[1] - aCoord[2]);
1070 aSlice[1] = 0.125*(1.0 + aCoord[0])*(1.0 - aCoord[1])*(1.0 - aCoord[2])*
1071 (-2.0 + aCoord[0] - aCoord[1] - aCoord[2]);
1072 aSlice[2] = 0.125*(1.0 + aCoord[0])*(1.0 + aCoord[1])*(1.0 - aCoord[2])*
1073 (-2.0 + aCoord[0] + aCoord[1] - aCoord[2]);
1074 aSlice[3] = 0.125*(1.0 - aCoord[0])*(1.0 + aCoord[1])*(1.0 - aCoord[2])*
1075 (-2.0 - aCoord[0] + aCoord[1] - aCoord[2]);
1076 aSlice[4] = 0.125*(1.0 - aCoord[0])*(1.0 - aCoord[1])*(1.0 + aCoord[2])*
1077 (-2.0 - aCoord[0] - aCoord[1] + aCoord[2]);
1078 aSlice[5] = 0.125*(1.0 + aCoord[0])*(1.0 - aCoord[1])*(1.0 + aCoord[2])*
1079 (-2.0 + aCoord[0] - aCoord[1] + aCoord[2]);
1080 aSlice[6] = 0.125*(1.0 + aCoord[0])*(1.0 + aCoord[1])*(1.0 + aCoord[2])*
1081 (-2.0 + aCoord[0] + aCoord[1] + aCoord[2]);
1082 aSlice[7] = 0.125*(1.0 - aCoord[0])*(1.0 + aCoord[1])*(1.0 + aCoord[2])*
1083 (-2.0 - aCoord[0] + aCoord[1] + aCoord[2]);
1085 aSlice[8] = 0.25*(1.0 - aCoord[0]*aCoord[0])*(1.0 - aCoord[1])*(1.0 - aCoord[2]);
1086 aSlice[9] = 0.25*(1.0 - aCoord[1]*aCoord[1])*(1.0 + aCoord[0])*(1.0 - aCoord[2]);
1087 aSlice[10] = 0.25*(1.0 - aCoord[0]*aCoord[0])*(1.0 + aCoord[1])*(1.0 - aCoord[2]);
1088 aSlice[11] = 0.25*(1.0 - aCoord[1]*aCoord[1])*(1.0 - aCoord[0])*(1.0 - aCoord[2]);
1089 aSlice[12] = 0.25*(1.0 - aCoord[2]*aCoord[2])*(1.0 - aCoord[0])*(1.0 - aCoord[1]);
1090 aSlice[13] = 0.25*(1.0 - aCoord[2]*aCoord[2])*(1.0 + aCoord[0])*(1.0 - aCoord[1]);
1091 aSlice[14] = 0.25*(1.0 - aCoord[2]*aCoord[2])*(1.0 + aCoord[0])*(1.0 + aCoord[1]);
1092 aSlice[15] = 0.25*(1.0 - aCoord[2]*aCoord[2])*(1.0 - aCoord[0])*(1.0 + aCoord[1]);
1093 aSlice[16] = 0.25*(1.0 - aCoord[0]*aCoord[0])*(1.0 - aCoord[1])*(1.0 + aCoord[2]);
1094 aSlice[17] = 0.25*(1.0 - aCoord[1]*aCoord[1])*(1.0 + aCoord[0])*(1.0 + aCoord[2]);
1095 aSlice[18] = 0.25*(1.0 - aCoord[0]*aCoord[0])*(1.0 + aCoord[1])*(1.0 + aCoord[2]);
1096 aSlice[19] = 0.25*(1.0 - aCoord[1]*aCoord[1])*(1.0 - aCoord[0])*(1.0 + aCoord[2]);
1100 //---------------------------------------------------------------
1101 THexa27a::THexa27a():
1104 TInt aNbRef = myRefCoord.size();
1105 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
1106 TCoordSlice aCoord = GetCoord(aRefId);
1108 case 20: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = -1.0; break;
1109 case 21: aCoord[0] = 0.0; aCoord[1] = -1.0; aCoord[2] = 0.0; break;
1110 case 22: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1111 case 23: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1112 case 24: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1113 case 25: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1114 case 26: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1120 THexa27a::InitFun(const TCCoordSliceArr& theRef,
1121 const TCCoordSliceArr& theGauss,
1124 GetFun(theRef,theGauss,theFun);
1126 TInt aNbGauss = theGauss.size();
1127 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1128 const TCCoordSlice& aCoord = theGauss[aGaussId];
1129 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
1131 aSlice[0] = 0.125*aCoord[0]*(aCoord[0] - 1.0)*aCoord[1]*(aCoord[1] - 1.0)*aCoord[2]*(aCoord[2] - 1.0);
1132 aSlice[1] = 0.125*aCoord[0]*(aCoord[0] + 1.0)*aCoord[1]*(aCoord[1] - 1.0)*aCoord[2]*(aCoord[2] - 1.0);
1133 aSlice[2] = 0.125*aCoord[0]*(aCoord[0] + 1.0)*aCoord[1]*(aCoord[1] + 1.0)*aCoord[2]*(aCoord[2] - 1.0);
1134 aSlice[3] = 0.125*aCoord[0]*(aCoord[0] - 1.0)*aCoord[1]*(aCoord[1] + 1.0)*aCoord[2]*(aCoord[2] - 1.0);
1135 aSlice[4] = 0.125*aCoord[0]*(aCoord[0] + 1.0)*aCoord[1]*(aCoord[1] - 1.0)*aCoord[2]*(aCoord[2] + 1.0);
1136 aSlice[5] = 0.125*aCoord[0]*(aCoord[0] + 1.0)*aCoord[1]*(aCoord[1] - 1.0)*aCoord[2]*(aCoord[2] + 1.0);
1137 aSlice[6] = 0.125*aCoord[0]*(aCoord[0] + 1.0)*aCoord[1]*(aCoord[1] + 1.0)*aCoord[2]*(aCoord[2] + 1.0);
1138 aSlice[7] = 0.125*aCoord[0]*(aCoord[0] - 1.0)*aCoord[1]*(aCoord[1] + 1.0)*aCoord[2]*(aCoord[2] + 1.0);
1140 aSlice[8] = 0.25*(1.0 - aCoord[0]*aCoord[0])*aCoord[1]*(aCoord[1] - 1.0)*aCoord[2]*(aCoord[2] - 1.0);
1141 aSlice[9] = 0.25*aCoord[0]*(aCoord[0] + 1.0)*(1.0 - aCoord[1]*aCoord[1])*aCoord[2]*(aCoord[2] - 1.0);
1142 aSlice[10]= 0.25*aCoord[0]*(aCoord[0] + 1.0)*(1.0 - aCoord[1]*aCoord[1])*aCoord[2]*(aCoord[2] - 1.0);
1143 aSlice[11]= 0.25*aCoord[0]*(aCoord[0] - 1.0)*(1.0 - aCoord[1]*aCoord[1])*aCoord[2]*(aCoord[2] - 1.0);
1144 aSlice[12]= 0.25*aCoord[0]*(aCoord[0] - 1.0)*aCoord[1]*(aCoord[1] - 1.0)*(1.0 - aCoord[2]*aCoord[2]);
1145 aSlice[13]= 0.25*aCoord[0]*(aCoord[0] + 1.0)*aCoord[1]*(aCoord[1] - 1.0)*(1.0 - aCoord[2]*aCoord[2]);
1146 aSlice[14]= 0.25*aCoord[0]*(aCoord[0] + 1.0)*aCoord[1]*(aCoord[1] + 1.0)*(1.0 - aCoord[2]*aCoord[2]);
1147 aSlice[15]= 0.25*aCoord[0]*(aCoord[0] - 1.0)*aCoord[1]*(aCoord[1] + 1.0)*(1.0 - aCoord[2]*aCoord[2]);
1148 aSlice[16]= 0.25*(1.0 - aCoord[0]*aCoord[0])*aCoord[1]*(aCoord[1] - 1.0)*aCoord[2]*(aCoord[2] + 1.0);
1149 aSlice[17]= 0.25*aCoord[0]*(aCoord[0] + 1.0)*(1.0 - aCoord[1]*aCoord[1])*aCoord[2]*(aCoord[2] + 1.0);
1150 aSlice[18]= 0.25*(1.0 - aCoord[0]*aCoord[0])*aCoord[1]*(aCoord[1] + 1.0)*aCoord[2]*(aCoord[2] + 1.0);
1151 aSlice[19]= 0.25*aCoord[0]*(aCoord[0] - 1.0)*(1.0 - aCoord[1]*aCoord[1])*aCoord[2]*(aCoord[2] + 1.0);
1152 aSlice[20]= 0.25*(1.0 - aCoord[0]*aCoord[0])*(1.0 - aCoord[1]*aCoord[1])*aCoord[2]*(aCoord[2] - 1.0);
1153 aSlice[21]= 0.25*(1.0 - aCoord[0]*aCoord[0])*aCoord[1]*(aCoord[1] - 1.0)*(1.0 - aCoord[2]*aCoord[2]);
1154 aSlice[22]= 0.25*aCoord[0]*(aCoord[0] + 1.0)*(1.0 - aCoord[1]*aCoord[1])*(1.0 - aCoord[2]*aCoord[2]);
1155 aSlice[23]= 0.25*(1.0 - aCoord[0]*aCoord[0])*aCoord[1]*(aCoord[1] + 1.0)*(1.0 - aCoord[2]*aCoord[2]);
1156 aSlice[24]= 0.25*aCoord[0]*(aCoord[0] - 1.0)*(1.0 - aCoord[1]*aCoord[1])*(1.0 - aCoord[2]*aCoord[2]);
1157 aSlice[25]= 0.25*(1.0 - aCoord[0]*aCoord[0])*(1.0 - aCoord[1]*aCoord[1])*aCoord[2]*(aCoord[2] + 1.0);
1158 aSlice[26]= 0.25*(1.0 - aCoord[0]*aCoord[0])*(1.0 - aCoord[1]*aCoord[1])*(1.0 - aCoord[1]*aCoord[1]);
1162 //---------------------------------------------------------------
1166 TInt aNbRef = GetNbRef();
1167 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
1168 TCoordSlice aCoord = GetCoord(aRefId);
1170 case 0: aCoord[0] = -1.0; aCoord[1] = -1.0; aCoord[2] = -1.0; break;
1171 case 3: aCoord[0] = 1.0; aCoord[1] = -1.0; aCoord[2] = -1.0; break;
1172 case 2: aCoord[0] = 1.0; aCoord[1] = 1.0; aCoord[2] = -1.0; break;
1173 case 1: aCoord[0] = -1.0; aCoord[1] = 1.0; aCoord[2] = -1.0; break;
1174 case 4: aCoord[0] = -1.0; aCoord[1] = -1.0; aCoord[2] = 1.0; break;
1175 case 7: aCoord[0] = 1.0; aCoord[1] = -1.0; aCoord[2] = 1.0; break;
1176 case 6: aCoord[0] = 1.0; aCoord[1] = 1.0; aCoord[2] = 1.0; break;
1177 case 5: aCoord[0] = -1.0; aCoord[1] = 1.0; aCoord[2] = 1.0; break;
1183 THexa8b::InitFun(const TCCoordSliceArr& theRef,
1184 const TCCoordSliceArr& theGauss,
1187 GetFun(theRef,theGauss,theFun);
1189 TInt aNbGauss = theGauss.size();
1190 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1191 const TCCoordSlice& aCoord = theGauss[aGaussId];
1192 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
1194 aSlice[0] = 0.125*(1.0 - aCoord[0])*(1.0 - aCoord[1])*(1.0 - aCoord[2]);
1195 aSlice[3] = 0.125*(1.0 + aCoord[0])*(1.0 - aCoord[1])*(1.0 - aCoord[2]);
1196 aSlice[2] = 0.125*(1.0 + aCoord[0])*(1.0 + aCoord[1])*(1.0 - aCoord[2]);
1197 aSlice[1] = 0.125*(1.0 - aCoord[0])*(1.0 + aCoord[1])*(1.0 - aCoord[2]);
1199 aSlice[4] = 0.125*(1.0 - aCoord[0])*(1.0 - aCoord[1])*(1.0 + aCoord[2]);
1200 aSlice[7] = 0.125*(1.0 + aCoord[0])*(1.0 - aCoord[1])*(1.0 + aCoord[2]);
1201 aSlice[6] = 0.125*(1.0 + aCoord[0])*(1.0 + aCoord[1])*(1.0 + aCoord[2]);
1202 aSlice[5] = 0.125*(1.0 - aCoord[0])*(1.0 + aCoord[1])*(1.0 + aCoord[2]);
1206 //---------------------------------------------------------------
1207 THexa20b::THexa20b(TInt theDim, TInt theNbRef):
1208 TShapeFun(theDim,theNbRef)
1210 TInt aNbRef = myRefCoord.size();
1211 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
1212 TCoordSlice aCoord = GetCoord(aRefId);
1214 case 0: aCoord[0] = -1.0; aCoord[1] = -1.0; aCoord[2] = -1.0; break;
1215 case 3: aCoord[0] = 1.0; aCoord[1] = -1.0; aCoord[2] = -1.0; break;
1216 case 2: aCoord[0] = 1.0; aCoord[1] = 1.0; aCoord[2] = -1.0; break;
1217 case 1: aCoord[0] = -1.0; aCoord[1] = 1.0; aCoord[2] = -1.0; break;
1218 case 4: aCoord[0] = -1.0; aCoord[1] = -1.0; aCoord[2] = 1.0; break;
1219 case 7: aCoord[0] = 1.0; aCoord[1] = -1.0; aCoord[2] = 1.0; break;
1220 case 6: aCoord[0] = 1.0; aCoord[1] = 1.0; aCoord[2] = 1.0; break;
1221 case 5: aCoord[0] = -1.0; aCoord[1] = 1.0; aCoord[2] = 1.0; break;
1223 case 11: aCoord[0] = 0.0; aCoord[1] = -1.0; aCoord[2] = -1.0; break;
1224 case 10: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = -1.0; break;
1225 case 9: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = -1.0; break;
1226 case 8: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = -1.0; break;
1227 case 16: aCoord[0] = -1.0; aCoord[1] = -1.0; aCoord[2] = 0.0; break;
1228 case 19: aCoord[0] = 1.0; aCoord[1] = -1.0; aCoord[2] = 0.0; break;
1229 case 18: aCoord[0] = 1.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1230 case 17: aCoord[0] = -1.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1231 case 15: aCoord[0] = 0.0; aCoord[1] = -1.0; aCoord[2] = 1.0; break;
1232 case 14: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1233 case 13: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = 1.0; break;
1234 case 12: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1240 THexa20b::InitFun(const TCCoordSliceArr& theRef,
1241 const TCCoordSliceArr& theGauss,
1244 GetFun(theRef,theGauss,theFun);
1246 TInt aNbGauss = theGauss.size();
1247 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1248 const TCCoordSlice& aCoord = theGauss[aGaussId];
1249 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
1251 aSlice[0] = 0.125*(1.0 - aCoord[0])*(1.0 - aCoord[1])*(1.0 - aCoord[2])*
1252 (-2.0 - aCoord[0] - aCoord[1] - aCoord[2]);
1253 aSlice[3] = 0.125*(1.0 + aCoord[0])*(1.0 - aCoord[1])*(1.0 - aCoord[2])*
1254 (-2.0 + aCoord[0] - aCoord[1] - aCoord[2]);
1255 aSlice[2] = 0.125*(1.0 + aCoord[0])*(1.0 + aCoord[1])*(1.0 - aCoord[2])*
1256 (-2.0 + aCoord[0] + aCoord[1] - aCoord[2]);
1257 aSlice[1] = 0.125*(1.0 - aCoord[0])*(1.0 + aCoord[1])*(1.0 - aCoord[2])*
1258 (-2.0 - aCoord[0] + aCoord[1] - aCoord[2]);
1259 aSlice[4] = 0.125*(1.0 - aCoord[0])*(1.0 - aCoord[1])*(1.0 + aCoord[2])*
1260 (-2.0 - aCoord[0] - aCoord[1] + aCoord[2]);
1261 aSlice[7] = 0.125*(1.0 + aCoord[0])*(1.0 - aCoord[1])*(1.0 + aCoord[2])*
1262 (-2.0 + aCoord[0] - aCoord[1] + aCoord[2]);
1263 aSlice[6] = 0.125*(1.0 + aCoord[0])*(1.0 + aCoord[1])*(1.0 + aCoord[2])*
1264 (-2.0 + aCoord[0] + aCoord[1] + aCoord[2]);
1265 aSlice[5] = 0.125*(1.0 - aCoord[0])*(1.0 + aCoord[1])*(1.0 + aCoord[2])*
1266 (-2.0 - aCoord[0] + aCoord[1] + aCoord[2]);
1268 aSlice[11] = 0.25*(1.0 - aCoord[0]*aCoord[0])*(1.0 - aCoord[1])*(1.0 - aCoord[2]);
1269 aSlice[10] = 0.25*(1.0 - aCoord[1]*aCoord[1])*(1.0 + aCoord[0])*(1.0 - aCoord[2]);
1270 aSlice[9] = 0.25*(1.0 - aCoord[0]*aCoord[0])*(1.0 + aCoord[1])*(1.0 - aCoord[2]);
1271 aSlice[8] = 0.25*(1.0 - aCoord[1]*aCoord[1])*(1.0 - aCoord[0])*(1.0 - aCoord[2]);
1272 aSlice[16] = 0.25*(1.0 - aCoord[2]*aCoord[2])*(1.0 - aCoord[0])*(1.0 - aCoord[1]);
1273 aSlice[19] = 0.25*(1.0 - aCoord[2]*aCoord[2])*(1.0 + aCoord[0])*(1.0 - aCoord[1]);
1274 aSlice[18] = 0.25*(1.0 - aCoord[2]*aCoord[2])*(1.0 + aCoord[0])*(1.0 + aCoord[1]);
1275 aSlice[17] = 0.25*(1.0 - aCoord[2]*aCoord[2])*(1.0 - aCoord[0])*(1.0 + aCoord[1]);
1276 aSlice[15] = 0.25*(1.0 - aCoord[0]*aCoord[0])*(1.0 - aCoord[1])*(1.0 + aCoord[2]);
1277 aSlice[14] = 0.25*(1.0 - aCoord[1]*aCoord[1])*(1.0 + aCoord[0])*(1.0 + aCoord[2]);
1278 aSlice[13] = 0.25*(1.0 - aCoord[0]*aCoord[0])*(1.0 + aCoord[1])*(1.0 + aCoord[2]);
1279 aSlice[12] = 0.25*(1.0 - aCoord[1]*aCoord[1])*(1.0 - aCoord[0])*(1.0 + aCoord[2]);
1283 //---------------------------------------------------------------
1284 TPenta6a::TPenta6a():
1287 TInt aNbRef = myRefCoord.size();
1288 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
1289 TCoordSlice aCoord = GetCoord(aRefId);
1291 case 0: aCoord[0] = -1.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1292 case 1: aCoord[0] = -1.0; aCoord[1] = -0.0; aCoord[2] = 1.0; break;
1293 case 2: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1294 case 3: aCoord[0] = 1.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1295 case 4: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1296 case 5: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1302 TPenta6a::InitFun(const TCCoordSliceArr& theRef,
1303 const TCCoordSliceArr& theGauss,
1306 GetFun(theRef,theGauss,theFun);
1308 TInt aNbGauss = theGauss.size();
1309 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1310 const TCCoordSlice& aCoord = theGauss[aGaussId];
1311 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
1313 aSlice[0] = 0.5*aCoord[1]*(1.0 - aCoord[0]);
1314 aSlice[1] = 0.5*aCoord[2]*(1.0 - aCoord[0]);
1315 aSlice[2] = 0.5*(1.0 - aCoord[1] - aCoord[2])*(1.0 - aCoord[0]);
1317 aSlice[3] = 0.5*aCoord[1]*(aCoord[0] + 1.0);
1318 aSlice[4] = 0.5*aCoord[2]*(aCoord[0] + 1.0);
1319 aSlice[5] = 0.5*(1.0 - aCoord[1] - aCoord[2])*(1.0 + aCoord[0]);
1323 //---------------------------------------------------------------
1324 TPenta6b::TPenta6b():
1327 TInt aNbRef = myRefCoord.size();
1328 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
1329 TCoordSlice aCoord = GetCoord(aRefId);
1331 case 0: aCoord[0] = -1.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1332 case 2: aCoord[0] = -1.0; aCoord[1] = -0.0; aCoord[2] = 1.0; break;
1333 case 1: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1334 case 3: aCoord[0] = 1.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1335 case 5: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1336 case 4: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1342 TPenta6b::InitFun(const TCCoordSliceArr& theRef,
1343 const TCCoordSliceArr& theGauss,
1346 GetFun(theRef,theGauss,theFun);
1348 TInt aNbGauss = theGauss.size();
1349 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1350 const TCCoordSlice& aCoord = theGauss[aGaussId];
1351 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
1353 aSlice[0] = 0.5*aCoord[1]*(1.0 - aCoord[0]);
1354 aSlice[2] = 0.5*aCoord[2]*(1.0 - aCoord[0]);
1355 aSlice[1] = 0.5*(1.0 - aCoord[1] - aCoord[2])*(1.0 - aCoord[0]);
1357 aSlice[3] = 0.5*aCoord[1]*(aCoord[0] + 1.0);
1358 aSlice[5] = 0.5*aCoord[2]*(aCoord[0] + 1.0);
1359 aSlice[4] = 0.5*(1.0 - aCoord[1] - aCoord[2])*(1.0 + aCoord[0]);
1363 //---------------------------------------------------------------
1364 TPenta15a::TPenta15a():
1367 TInt aNbRef = myRefCoord.size();
1368 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
1369 TCoordSlice aCoord = GetCoord(aRefId);
1371 case 0: aCoord[0] = -1.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1372 case 1: aCoord[0] = -1.0; aCoord[1] = -0.0; aCoord[2] = 1.0; break;
1373 case 2: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1374 case 3: aCoord[0] = 1.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1375 case 4: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1376 case 5: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1378 case 6: aCoord[0] = -1.0; aCoord[1] = 0.5; aCoord[2] = 0.5; break;
1379 case 7: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = 0.5; break;
1380 case 8: aCoord[0] = -1.0; aCoord[1] = 0.5; aCoord[2] = 0.0; break;
1381 case 9: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1382 case 10: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1383 case 11: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1384 case 12: aCoord[0] = 1.0; aCoord[1] = 0.5; aCoord[2] = 0.5; break;
1385 case 13: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.5; break;
1386 case 14: aCoord[0] = 1.0; aCoord[1] = 0.5; aCoord[2] = 0.0; break;
1392 TPenta15a::InitFun(const TCCoordSliceArr& theRef,
1393 const TCCoordSliceArr& theGauss,
1396 GetFun(theRef,theGauss,theFun);
1398 TInt aNbGauss = theGauss.size();
1399 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1400 const TCCoordSlice& aCoord = theGauss[aGaussId];
1401 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
1403 aSlice[0] = 0.5*aCoord[1]*(1.0 - aCoord[0])*(2.0*aCoord[1] - 2.0 - aCoord[0]);
1404 aSlice[1] = 0.5*aCoord[2]*(1.0 - aCoord[0])*(2.0*aCoord[2] - 2.0 - aCoord[0]);
1405 aSlice[2] = 0.5*(aCoord[0] - 1.0)*(1.0 - aCoord[1] - aCoord[2])*(aCoord[0] + 2.0*aCoord[1] + 2.0*aCoord[2]);
1407 aSlice[3] = 0.5*aCoord[1]*(1.0 + aCoord[0])*(2.0*aCoord[1] - 2.0 + aCoord[0]);
1408 aSlice[4] = 0.5*aCoord[2]*(1.0 + aCoord[0])*(2.0*aCoord[2] - 2.0 + aCoord[0]);
1409 aSlice[5] = 0.5*(-aCoord[0] - 1.0)*(1.0 - aCoord[1] - aCoord[2])*(-aCoord[0] + 2.0*aCoord[1] + 2.0*aCoord[2]);
1411 aSlice[6] = 2.0*aCoord[1]*aCoord[2]*(1.0 - aCoord[0]);
1412 aSlice[7] = 2.0*aCoord[2]*(1.0 - aCoord[1] - aCoord[2])*(1.0 - aCoord[0]);
1413 aSlice[8] = 2.0*aCoord[1]*(1.0 - aCoord[1] - aCoord[2])*(1.0 - aCoord[0]);
1415 aSlice[9] = aCoord[1]*(1.0 - aCoord[0]*aCoord[0]);
1416 aSlice[10] = aCoord[2]*(1.0 - aCoord[0]*aCoord[0]);
1417 aSlice[11] = (1.0 - aCoord[1] - aCoord[2])*(1.0 - aCoord[0]*aCoord[0]);
1419 aSlice[12] = 2.0*aCoord[1]*aCoord[2]*(1.0 + aCoord[0]);
1420 aSlice[13] = 2.0*aCoord[2]*(1.0 - aCoord[1] - aCoord[2])*(1.0 + aCoord[0]);
1421 aSlice[14] = 2.0*aCoord[1]*(1.0 - aCoord[1] - aCoord[2])*(1.0 + aCoord[0]);
1425 //---------------------------------------------------------------
1426 TPenta15b::TPenta15b():
1429 TInt aNbRef = myRefCoord.size();
1430 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
1431 TCoordSlice aCoord = GetCoord(aRefId);
1433 case 0: aCoord[0] = -1.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1434 case 2: aCoord[0] = -1.0; aCoord[1] = -0.0; aCoord[2] = 1.0; break;
1435 case 1: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1436 case 3: aCoord[0] = 1.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1437 case 5: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1438 case 4: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1440 case 8: aCoord[0] = -1.0; aCoord[1] = 0.5; aCoord[2] = 0.5; break;
1441 case 7: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = 0.5; break;
1442 case 6: aCoord[0] = -1.0; aCoord[1] = 0.5; aCoord[2] = 0.0; break;
1443 case 12: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1444 case 14: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1445 case 13: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1446 case 11: aCoord[0] = 1.0; aCoord[1] = 0.5; aCoord[2] = 0.5; break;
1447 case 10: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.5; break;
1448 case 9: aCoord[0] = 1.0; aCoord[1] = 0.5; aCoord[2] = 0.0; break;
1454 TPenta15b::InitFun(const TCCoordSliceArr& theRef,
1455 const TCCoordSliceArr& theGauss,
1458 GetFun(theRef,theGauss,theFun);
1460 TInt aNbGauss = theGauss.size();
1461 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1462 const TCCoordSlice& aCoord = theGauss[aGaussId];
1463 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
1465 aSlice[0] = 0.5*aCoord[1]*(1.0 - aCoord[0])*(2.0*aCoord[1] - 2.0 - aCoord[0]);
1466 aSlice[2] = 0.5*aCoord[2]*(1.0 - aCoord[0])*(2.0*aCoord[2] - 2.0 - aCoord[0]);
1467 aSlice[1] = 0.5*(aCoord[0] - 1.0)*(1.0 - aCoord[1] - aCoord[2])*(aCoord[0] + 2.0*aCoord[1] + 2.0*aCoord[2]);
1469 aSlice[3] = 0.5*aCoord[1]*(1.0 + aCoord[0])*(2.0*aCoord[1] - 2.0 + aCoord[0]);
1470 aSlice[5] = 0.5*aCoord[2]*(1.0 + aCoord[0])*(2.0*aCoord[2] - 2.0 + aCoord[0]);
1471 aSlice[4] = 0.5*(-aCoord[0] - 1.0)*(1.0 - aCoord[1] - aCoord[2])*(-aCoord[0] + 2.0*aCoord[1] + 2.0*aCoord[2]);
1473 aSlice[8] = 2.0*aCoord[1]*aCoord[2]*(1.0 - aCoord[0]);
1474 aSlice[7] = 2.0*aCoord[2]*(1.0 - aCoord[1] - aCoord[2])*(1.0 - aCoord[0]);
1475 aSlice[6] = 2.0*aCoord[1]*(1.0 - aCoord[1] - aCoord[2])*(1.0 - aCoord[0]);
1477 aSlice[12] = aCoord[1]*(1.0 - aCoord[0]*aCoord[0]);
1478 aSlice[14] = aCoord[2]*(1.0 - aCoord[0]*aCoord[0]);
1479 aSlice[13] = (1.0 - aCoord[1] - aCoord[2])*(1.0 - aCoord[0]*aCoord[0]);
1481 aSlice[11] = 2.0*aCoord[1]*aCoord[2]*(1.0 + aCoord[0]);
1482 aSlice[10] = 2.0*aCoord[2]*(1.0 - aCoord[1] - aCoord[2])*(1.0 + aCoord[0]);
1483 aSlice[9] = 2.0*aCoord[1]*(1.0 - aCoord[1] - aCoord[2])*(1.0 + aCoord[0]);
1487 //---------------------------------------------------------------
1491 TInt aNbRef = myRefCoord.size();
1492 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
1493 TCoordSlice aCoord = GetCoord(aRefId);
1495 case 0: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1496 case 1: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1497 case 2: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1498 case 3: aCoord[0] = 0.0; aCoord[1] = -1.0; aCoord[2] = 0.0; break;
1499 case 4: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1505 TPyra5a::InitFun(const TCCoordSliceArr& theRef,
1506 const TCCoordSliceArr& theGauss,
1509 GetFun(theRef,theGauss,theFun);
1511 TInt aNbGauss = theGauss.size();
1512 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1513 const TCCoordSlice& aCoord = theGauss[aGaussId];
1514 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
1516 // Fix for 0019920: EDF 788 VISU: Bad localisation of gauss points for pyramids
1517 // Seems shape function for ePYRA5 elements is:
1518 // w0 = 0.25*(-X + Y - 1.0)*(-X - Y - 1.0)*(1.0 - Z);
1519 // w1 = 0.25*(-X - Y - 1.0)*(+X - Y - 1.0)*(1.0 - Z);
1520 // w2 = 0.25*(+X + Y - 1.0)*(+X - Y - 1.0)*(1.0 - Z);
1521 // w3 = 0.25*(+X + Y - 1.0)*(-X + Y - 1.0)*(1.0 - Z);
1523 aSlice[0] = 0.25*(-aCoord[0] + aCoord[1] - 1.0)*(-aCoord[0] - aCoord[1] - 1.0)*(1.0 - aCoord[2]);
1524 aSlice[1] = 0.25*(-aCoord[0] - aCoord[1] - 1.0)*(+aCoord[0] - aCoord[1] - 1.0)*(1.0 - aCoord[2]);
1525 aSlice[2] = 0.25*(+aCoord[0] + aCoord[1] - 1.0)*(+aCoord[0] - aCoord[1] - 1.0)*(1.0 - aCoord[2]);
1526 aSlice[3] = 0.25*(+aCoord[0] + aCoord[1] - 1.0)*(-aCoord[0] + aCoord[1] - 1.0)*(1.0 - aCoord[2]);
1527 aSlice[4] = aCoord[2];
1531 //---------------------------------------------------------------
1535 TInt aNbRef = myRefCoord.size();
1536 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
1537 TCoordSlice aCoord = GetCoord(aRefId);
1539 case 0: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1540 case 3: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1541 case 2: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1542 case 1: aCoord[0] = 0.0; aCoord[1] = -1.0; aCoord[2] = 0.0; break;
1543 case 4: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1549 TPyra5b::InitFun(const TCCoordSliceArr& theRef,
1550 const TCCoordSliceArr& theGauss,
1553 GetFun(theRef,theGauss,theFun);
1555 TInt aNbGauss = theGauss.size();
1556 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1557 const TCCoordSlice& aCoord = theGauss[aGaussId];
1558 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
1560 // Fix for 0019920: EDF 788 VISU: Bad localisation of gauss points for pyramids
1561 // Seems shape function for ePYRA5 elements is:
1562 // w0 = 0.25*(-X + Y - 1.0)*(-X - Y - 1.0)*(1.0 - Z);
1563 // w1 = 0.25*(-X - Y - 1.0)*(+X - Y - 1.0)*(1.0 - Z);
1564 // w2 = 0.25*(+X + Y - 1.0)*(+X - Y - 1.0)*(1.0 - Z);
1565 // w3 = 0.25*(+X + Y - 1.0)*(-X + Y - 1.0)*(1.0 - Z);
1567 aSlice[0] = 0.25*(-aCoord[0] + aCoord[1] - 1.0)*(-aCoord[0] - aCoord[1] - 1.0)*(1.0 - aCoord[2]);
1568 aSlice[3] = 0.25*(-aCoord[0] - aCoord[1] - 1.0)*(+aCoord[0] - aCoord[1] - 1.0)*(1.0 - aCoord[2]);
1569 aSlice[2] = 0.25*(+aCoord[0] + aCoord[1] - 1.0)*(+aCoord[0] - aCoord[1] - 1.0)*(1.0 - aCoord[2]);
1570 aSlice[1] = 0.25*(+aCoord[0] + aCoord[1] - 1.0)*(-aCoord[0] + aCoord[1] - 1.0)*(1.0 - aCoord[2]);
1571 aSlice[4] = aCoord[2];
1575 //---------------------------------------------------------------
1576 TPyra13a::TPyra13a():
1579 TInt aNbRef = myRefCoord.size();
1580 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
1581 TCoordSlice aCoord = GetCoord(aRefId);
1583 case 0: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1584 case 1: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1585 case 2: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1586 case 3: aCoord[0] = 0.0; aCoord[1] = -1.0; aCoord[2] = 0.0; break;
1587 case 4: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1589 case 5: aCoord[0] = 0.5; aCoord[1] = 0.5; aCoord[2] = 0.0; break;
1590 case 6: aCoord[0] = -0.5; aCoord[1] = 0.5; aCoord[2] = 0.0; break;
1591 case 7: aCoord[0] = -0.5; aCoord[1] = -0.5; aCoord[2] = 0.0; break;
1592 case 8: aCoord[0] = 0.5; aCoord[1] = -0.5; aCoord[2] = 0.0; break;
1593 case 9: aCoord[0] = 0.5; aCoord[1] = 0.0; aCoord[2] = 0.5; break;
1594 case 10: aCoord[0] = 0.0; aCoord[1] = 0.5; aCoord[2] = 0.5; break;
1595 case 11: aCoord[0] = -0.5; aCoord[1] = 0.0; aCoord[2] = 0.5; break;
1596 case 12: aCoord[0] = 0.0; aCoord[1] = -0.5; aCoord[2] = 0.5; break;
1602 TPyra13a::InitFun(const TCCoordSliceArr& theRef,
1603 const TCCoordSliceArr& theGauss,
1606 GetFun(theRef,theGauss,theFun);
1608 TInt aNbGauss = theGauss.size();
1609 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1610 const TCCoordSlice& aCoord = theGauss[aGaussId];
1611 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
1613 aSlice[0] = 0.5*(-aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*(-aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*
1614 (aCoord[0] - 0.5)/(1.0 - aCoord[2]);
1615 aSlice[1] = 0.5*(-aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*(+aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*
1616 (aCoord[1] - 0.5)/(1.0 - aCoord[2]);
1617 aSlice[2] = 0.5*(+aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*(+aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*
1618 (-aCoord[0] - 0.5)/(1.0 - aCoord[2]);
1619 aSlice[3] = 0.5*(+aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*(-aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*
1620 (-aCoord[1] - 0.5)/(1.0 - aCoord[2]);
1622 aSlice[4] = 2.0*aCoord[2]*(aCoord[2] - 0.5);
1624 aSlice[5] = 0.5*(-aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*(-aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*
1625 (aCoord[0] - aCoord[1] + aCoord[2] - 1.0)/(1.0 - aCoord[2]);
1626 aSlice[6] = 0.5*(-aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*(aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*
1627 (aCoord[0] + aCoord[1] + aCoord[2] - 1.0)/(1.0 - aCoord[2]);
1628 aSlice[7] = 0.5*(aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*(aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*
1629 (-aCoord[0] + aCoord[1] + aCoord[2] - 1.0)/(1.0 - aCoord[2]);
1630 aSlice[8] = 0.5*(aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*(-aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*
1631 (-aCoord[0] - aCoord[1] + aCoord[2] - 1.0)/(1.0 - aCoord[2]);
1633 aSlice[9] = 0.5*aCoord[2]*(-aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*(-aCoord[0] - aCoord[1] + aCoord[2] - 1.0)/
1635 aSlice[10] = 0.5*aCoord[2]*(-aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*(aCoord[0] - aCoord[1] + aCoord[2] - 1.0)/
1637 aSlice[11] = 0.5*aCoord[2]*(aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*(aCoord[0] + aCoord[1] + aCoord[2] - 1.0)/
1639 aSlice[12] = 0.5*aCoord[2]*(aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*(-aCoord[0] + aCoord[1] + aCoord[2] - 1.0)/
1644 //---------------------------------------------------------------
1645 TPyra13b::TPyra13b():
1648 TInt aNbRef = myRefCoord.size();
1649 for(TInt aRefId = 0; aRefId < aNbRef; aRefId++){
1650 TCoordSlice aCoord = GetCoord(aRefId);
1652 case 0: aCoord[0] = 1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1653 case 3: aCoord[0] = 0.0; aCoord[1] = 1.0; aCoord[2] = 0.0; break;
1654 case 2: aCoord[0] = -1.0; aCoord[1] = 0.0; aCoord[2] = 0.0; break;
1655 case 1: aCoord[0] = 0.0; aCoord[1] = -1.0; aCoord[2] = 0.0; break;
1656 case 4: aCoord[0] = 0.0; aCoord[1] = 0.0; aCoord[2] = 1.0; break;
1658 case 8: aCoord[0] = 0.5; aCoord[1] = 0.5; aCoord[2] = 0.0; break;
1659 case 7: aCoord[0] = -0.5; aCoord[1] = 0.5; aCoord[2] = 0.0; break;
1660 case 6: aCoord[0] = -0.5; aCoord[1] = -0.5; aCoord[2] = 0.0; break;
1661 case 5: aCoord[0] = 0.5; aCoord[1] = -0.5; aCoord[2] = 0.0; break;
1662 case 9: aCoord[0] = 0.5; aCoord[1] = 0.0; aCoord[2] = 0.5; break;
1663 case 12: aCoord[0] = 0.0; aCoord[1] = 0.5; aCoord[2] = 0.5; break;
1664 case 11: aCoord[0] = -0.5; aCoord[1] = 0.0; aCoord[2] = 0.5; break;
1665 case 10: aCoord[0] = 0.0; aCoord[1] = -0.5; aCoord[2] = 0.5; break;
1671 TPyra13b::InitFun(const TCCoordSliceArr& theRef,
1672 const TCCoordSliceArr& theGauss,
1675 GetFun(theRef,theGauss,theFun);
1677 TInt aNbGauss = theGauss.size();
1678 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1679 const TCCoordSlice& aCoord = theGauss[aGaussId];
1680 TFloatVecSlice aSlice = theFun.GetFunSlice(aGaussId);
1682 aSlice[0] = 0.5*(-aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*(-aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*
1683 (aCoord[0] - 0.5)/(1.0 - aCoord[2]);
1684 aSlice[3] = 0.5*(-aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*(+aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*
1685 (aCoord[1] - 0.5)/(1.0 - aCoord[2]);
1686 aSlice[2] = 0.5*(+aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*(+aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*
1687 (-aCoord[0] - 0.5)/(1.0 - aCoord[2]);
1688 aSlice[1] = 0.5*(+aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*(-aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*
1689 (-aCoord[1] - 0.5)/(1.0 - aCoord[2]);
1691 aSlice[4] = 2.0*aCoord[2]*(aCoord[2] - 0.5);
1693 aSlice[8] = 0.5*(-aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*(-aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*
1694 (aCoord[0] - aCoord[1] + aCoord[2] - 1.0)/(1.0 - aCoord[2]);
1695 aSlice[7] = 0.5*(-aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*(aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*
1696 (aCoord[0] + aCoord[1] + aCoord[2] - 1.0)/(1.0 - aCoord[2]);
1697 aSlice[6] = 0.5*(aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*(aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*
1698 (-aCoord[0] + aCoord[1] + aCoord[2] - 1.0)/(1.0 - aCoord[2]);
1699 aSlice[5] = 0.5*(aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*(-aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*
1700 (-aCoord[0] - aCoord[1] + aCoord[2] - 1.0)/(1.0 - aCoord[2]);
1702 aSlice[9] = 0.5*aCoord[2]*(-aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*(-aCoord[0] - aCoord[1] + aCoord[2] - 1.0)/
1704 aSlice[12] = 0.5*aCoord[2]*(-aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*(aCoord[0] - aCoord[1] + aCoord[2] - 1.0)/
1706 aSlice[11] = 0.5*aCoord[2]*(aCoord[0] - aCoord[1] + aCoord[2] - 1.0)*(aCoord[0] + aCoord[1] + aCoord[2] - 1.0)/
1708 aSlice[10] = 0.5*aCoord[2]*(aCoord[0] + aCoord[1] + aCoord[2] - 1.0)*(-aCoord[0] + aCoord[1] + aCoord[2] - 1.0)/
1713 //---------------------------------------------------------------
1715 GetGaussCoord3D(const TGaussInfo& theGaussInfo,
1716 const TCellInfo& theCellInfo,
1717 const TNodeInfo& theNodeInfo,
1718 TGaussCoord& theGaussCoord,
1719 const TElemNum& theElemNum,
1720 EModeSwitch theMode)
1722 INITMSG("GetGaussCoord3D\n");
1724 if(theGaussInfo.myGeom == theCellInfo.myGeom){
1725 EGeometrieElement aGeom = theGaussInfo.myGeom;
1727 TInt aNbRef = theGaussInfo.GetNbRef();
1728 TCCoordSliceArr aRefSlice(aNbRef);
1729 for(TInt anId = 0; anId < aNbRef; anId++)
1730 aRefSlice[anId] = theGaussInfo.GetRefCoordSlice(anId);
1732 TInt aNbGauss = theGaussInfo.GetNbGauss();
1733 TCCoordSliceArr aGaussSlice(aNbGauss);
1734 for(TInt anId = 0; anId < aNbGauss; anId++)
1735 aGaussSlice[anId] = theGaussInfo.GetGaussCoordSlice(anId);
1739 INITMSG("eSEG2"<<std::endl);
1741 if(TSeg2a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1747 INITMSG("eSEG3"<<std::endl);
1749 if(TSeg3a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1755 INITMSG("eTRIA3"<<std::endl);
1757 if(TTria3a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1760 if(TTria3b().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1766 INITMSG("eTRIA6"<<std::endl);
1768 if(TTria6a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1771 if(TTria6b().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1777 INITMSG("eQUAD4"<<std::endl);
1779 if(TQuad4a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1782 if(TQuad4b().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1788 INITMSG("eQUAD8"<<std::endl);
1790 if(TQuad8a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1793 if(TQuad8b().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1799 INITMSG("eQUAD9"<<std::endl);
1801 if(TQuad9a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1804 if(TQuad9b().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1810 INITMSG("eTETRA4"<<std::endl);
1812 if(TTetra4a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1815 if(TTetra4b().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1821 INITMSG("ePYRA5"<<std::endl);
1823 if(TPyra5a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1826 if(TPyra5b().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1832 INITMSG("ePENTA6"<<std::endl);
1834 if(TPenta6a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1837 if(TPenta6b().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1843 INITMSG("eHEXA8"<<std::endl);
1845 if(THexa8a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1848 if(THexa8b().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1854 INITMSG("eTETRA10"<<std::endl);
1856 if(TTetra10a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1859 if(TTetra10b().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1865 INITMSG("ePYRA13"<<std::endl);
1867 if(TPyra13a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1870 if(TPyra13b().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1876 INITMSG("ePENTA15"<<std::endl);
1878 if(TPenta15a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1881 if(TPenta15b().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1887 INITMSG("eHEXA20"<<std::endl);
1889 if(THexa20a().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1892 if(THexa20b().Eval(theCellInfo,theNodeInfo,theElemNum,aRefSlice,aGaussSlice,theGaussCoord,theMode))
1898 INITMSG("eNONE"<<std::endl);
1906 //---------------------------------------------------------------
1908 GetBaryCenter(const TCellInfo& theCellInfo,
1909 const TNodeInfo& theNodeInfo,
1910 TGaussCoord& theGaussCoord,
1911 const TElemNum& theElemNum,
1912 EModeSwitch theMode)
1914 INITMSG("GetBaryCenter\n");
1915 const PMeshInfo& aMeshInfo = theCellInfo.GetMeshInfo();
1916 TInt aDim = aMeshInfo->GetDim();
1917 static TInt aNbGauss = 1;
1919 bool anIsSubMesh = !theElemNum.empty();
1922 aNbElem = theElemNum.size();
1924 aNbElem = theCellInfo.GetNbElem();
1926 theGaussCoord.Init(aNbElem,aNbGauss,aDim,theMode);
1928 TInt aConnDim = theCellInfo.GetConnDim();
1932 "; aNbGauss = "<<aNbGauss<<
1933 "; aNbElem = "<<aNbElem<<
1934 "; aNbNodes = "<<theNodeInfo.GetNbElem()<<
1937 for(TInt anElemId = 0; anElemId < aNbElem; anElemId++){
1938 TInt aCellId = anIsSubMesh? theElemNum[anElemId]-1: anElemId;
1939 TCConnSlice aConnSlice = theCellInfo.GetConnSlice(aCellId);
1940 TCoordSliceArr aCoordSliceArr = theGaussCoord.GetCoordSliceArr(anElemId);
1942 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1943 TCoordSlice& aGaussCoordSlice = aCoordSliceArr[aGaussId];
1945 for(TInt aConnId = 0; aConnId < aConnDim; aConnId++){
1946 TInt aNodeId = aConnSlice[aConnId] - 1;
1947 TCCoordSlice aNodeCoordSlice = theNodeInfo.GetCoordSlice(aNodeId);
1949 for(TInt aDimId = 0; aDimId < aDim; aDimId++){
1950 aGaussCoordSlice[aDimId] += aNodeCoordSlice[aDimId];
1954 for(TInt aDimId = 0; aDimId < aDim; aDimId++){
1955 aGaussCoordSlice[aDimId] /= aConnDim;
1960 if (SALOME::VerbosityActivated())
1962 for(TInt anElemId = 0; anElemId < aNbElem; anElemId++){
1963 TCoordSliceArr aCoordSliceArr = theGaussCoord.GetCoordSliceArr(anElemId);
1965 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
1966 TCoordSlice& aCoordSlice = aCoordSliceArr[aGaussId];
1968 for(TInt aDimId = 0; aDimId < aDim; aDimId++){
1969 ADDMSG(aCoordSlice[aDimId]<<" ");
1980 //---------------------------------------------------------------
1982 GetBaryCenter(const TPolygoneInfo& thePolygoneInfo,
1983 const TNodeInfo& theNodeInfo,
1984 TGaussCoord& theGaussCoord,
1985 const TElemNum& theElemNum,
1986 EModeSwitch theMode)
1988 INITMSG("GetBaryCenter\n");
1989 const PMeshInfo& aMeshInfo = thePolygoneInfo.GetMeshInfo();
1990 TInt aDim = aMeshInfo->GetDim();
1991 static TInt aNbGauss = 1;
1993 bool anIsSubMesh = !theElemNum.empty();
1996 aNbElem = theElemNum.size();
1998 aNbElem = thePolygoneInfo.GetNbElem();
2000 theGaussCoord.Init(aNbElem,aNbGauss,aDim,theMode);
2004 "; aNbGauss = "<<aNbGauss<<
2005 "; aNbElem = "<<aNbElem<<
2006 "; aNbNodes = "<<theNodeInfo.GetNbElem()<<
2009 for(TInt anElemId = 0; anElemId < aNbElem; anElemId++){
2010 TInt aCellId = anIsSubMesh? theElemNum[anElemId]-1: anElemId;
2012 TCoordSliceArr aCoordSliceArr = theGaussCoord.GetCoordSliceArr(anElemId);
2013 TCConnSlice aConnSlice = thePolygoneInfo.GetConnSlice(aCellId);
2014 TInt aNbConn = thePolygoneInfo.GetNbConn(aCellId);
2015 TInt aNbNodes = thePolygoneInfo.GetNbConn(aCellId);
2017 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
2018 TCoordSlice& aGaussCoordSlice = aCoordSliceArr[aGaussId];
2020 for(TInt aConnId = 0; aConnId < aNbConn; aConnId++){
2021 TInt aNodeId = aConnSlice[aConnId] - 1;
2022 TCCoordSlice aNodeCoordSlice = theNodeInfo.GetCoordSlice(aNodeId);
2024 for(TInt aDimId = 0; aDimId < aDim; aDimId++){
2025 aGaussCoordSlice[aDimId] += aNodeCoordSlice[aDimId];
2029 for(TInt aDimId = 0; aDimId < aDim; aDimId++){
2030 aGaussCoordSlice[aDimId] /= aNbNodes;
2038 //---------------------------------------------------------------
2040 GetBaryCenter(const TPolyedreInfo& thePolyedreInfo,
2041 const TNodeInfo& theNodeInfo,
2042 TGaussCoord& theGaussCoord,
2043 const TElemNum& theElemNum,
2044 EModeSwitch theMode)
2046 INITMSG("GetBaryCenter\n");
2047 const PMeshInfo& aMeshInfo = thePolyedreInfo.GetMeshInfo();
2048 TInt aDim = aMeshInfo->GetDim();
2049 static TInt aNbGauss = 1;
2051 bool anIsSubMesh = !theElemNum.empty();
2054 aNbElem = theElemNum.size();
2056 aNbElem = thePolyedreInfo.GetNbElem();
2058 theGaussCoord.Init(aNbElem,aNbGauss,aDim,theMode);
2062 "; aNbGauss = "<<aNbGauss<<
2063 "; aNbElem = "<<aNbElem<<
2064 "; aNbNodes = "<<theNodeInfo.GetNbElem()<<
2067 for(TInt anElemId = 0; anElemId < aNbElem; anElemId++){
2068 TInt aCellId = anIsSubMesh? theElemNum[anElemId]-1: anElemId;
2070 TCoordSliceArr aCoordSliceArr = theGaussCoord.GetCoordSliceArr(anElemId);
2071 TCConnSliceArr aConnSliceArr = thePolyedreInfo.GetConnSliceArr(aCellId);
2072 TInt aNbFaces = aConnSliceArr.size();
2074 TInt aNbNodes = thePolyedreInfo.GetNbNodes(aCellId);
2076 for(TInt aGaussId = 0; aGaussId < aNbGauss; aGaussId++){
2077 TCoordSlice& aGaussCoordSlice = aCoordSliceArr[aGaussId];
2079 for(TInt aFaceId = 0; aFaceId < aNbFaces; aFaceId++){
2080 TCConnSlice aConnSlice = aConnSliceArr[aFaceId];
2081 TInt aNbConn = aConnSlice.size();
2082 for(TInt aConnId = 0; aConnId < aNbConn; aConnId++){
2083 TInt aNodeId = aConnSlice[aConnId] - 1;
2084 TCCoordSlice aNodeCoordSlice = theNodeInfo.GetCoordSlice(aNodeId);
2086 for(TInt aDimId = 0; aDimId < aDim; aDimId++){
2087 aGaussCoordSlice[aDimId] += aNodeCoordSlice[aDimId];
2091 for(TInt aDimId = 0; aDimId < aDim; aDimId++){
2092 aGaussCoordSlice[aDimId] /= aNbNodes;