KratosMultiphysics
KRATOS Multiphysics (Kratos) is a framework for building parallel, multi-disciplinary simulation software, aiming at modularity, extensibility, and high performance. Kratos is written in C++, and counts with an extensive Python interface.
hexahedra_3d_20.h
Go to the documentation of this file.
1 // | / |
2 // ' / __| _` | __| _ \ __|
3 // . \ | ( | | ( |\__ `
4 // _|\_\_| \__,_|\__|\___/ ____/
5 // Multi-Physics
6 //
7 // License: BSD License
8 // Kratos default license: kratos/license.txt
9 //
10 // Main authors: Riccardo Rossi
11 // Janosch Stascheit
12 // Felix Nagel
13 // contributors: Hoang Giang Bui
14 // Josep Maria Carbonell
15 //
16 
17 #pragma once
18 
19 // System includes
20 
21 // External includes
22 
23 // Project includes
27 
28 namespace Kratos
29 {
50 template<class TPointType> class Hexahedra3D20 : public Geometry<TPointType>
51 {
52 public:
53 
62 
68 
73 
78 
84 
88  typedef TPointType PointType;
89 
95  typedef typename BaseType::IndexType IndexType;
96 
102  typedef typename BaseType::SizeType SizeType;
103 
109 
116 
124 
131 
138 
145 
152 
159 
164 
169 
174 
175 
179  Hexahedra3D20( const PointType& Point1, const PointType& Point2,
180  const PointType& Point3, const PointType& Point4,
181  const PointType& Point5, const PointType& Point6,
182  const PointType& Point7, const PointType& Point8,
183  const PointType& Point9, const PointType& Point10,
184  const PointType& Point11, const PointType& Point12,
185  const PointType& Point13, const PointType& Point14,
186  const PointType& Point15, const PointType& Point16,
187  const PointType& Point17, const PointType& Point18,
188  const PointType& Point19, const PointType& Point20
189  )
190  : BaseType( PointsArrayType(), &msGeometryData )
191  {
192  this->Points().push_back( typename PointType::Pointer( new PointType( Point1 ) ) );
193  this->Points().push_back( typename PointType::Pointer( new PointType( Point2 ) ) );
194  this->Points().push_back( typename PointType::Pointer( new PointType( Point3 ) ) );
195  this->Points().push_back( typename PointType::Pointer( new PointType( Point4 ) ) );
196  this->Points().push_back( typename PointType::Pointer( new PointType( Point5 ) ) );
197  this->Points().push_back( typename PointType::Pointer( new PointType( Point6 ) ) );
198  this->Points().push_back( typename PointType::Pointer( new PointType( Point7 ) ) );
199  this->Points().push_back( typename PointType::Pointer( new PointType( Point8 ) ) );
200  this->Points().push_back( typename PointType::Pointer( new PointType( Point9 ) ) );
201  this->Points().push_back( typename PointType::Pointer( new PointType( Point10 ) ) );
202  this->Points().push_back( typename PointType::Pointer( new PointType( Point11 ) ) );
203  this->Points().push_back( typename PointType::Pointer( new PointType( Point12 ) ) );
204  this->Points().push_back( typename PointType::Pointer( new PointType( Point13 ) ) );
205  this->Points().push_back( typename PointType::Pointer( new PointType( Point14 ) ) );
206  this->Points().push_back( typename PointType::Pointer( new PointType( Point15 ) ) );
207  this->Points().push_back( typename PointType::Pointer( new PointType( Point16 ) ) );
208  this->Points().push_back( typename PointType::Pointer( new PointType( Point17 ) ) );
209  this->Points().push_back( typename PointType::Pointer( new PointType( Point18 ) ) );
210  this->Points().push_back( typename PointType::Pointer( new PointType( Point19 ) ) );
211  this->Points().push_back( typename PointType::Pointer( new PointType( Point20 ) ) );
212  }
213 
214  Hexahedra3D20( typename PointType::Pointer pPoint1,
215  typename PointType::Pointer pPoint2,
216  typename PointType::Pointer pPoint3,
217  typename PointType::Pointer pPoint4,
218  typename PointType::Pointer pPoint5,
219  typename PointType::Pointer pPoint6,
220  typename PointType::Pointer pPoint7,
221  typename PointType::Pointer pPoint8,
222  typename PointType::Pointer pPoint9,
223  typename PointType::Pointer pPoint10,
224  typename PointType::Pointer pPoint11,
225  typename PointType::Pointer pPoint12,
226  typename PointType::Pointer pPoint13,
227  typename PointType::Pointer pPoint14,
228  typename PointType::Pointer pPoint15,
229  typename PointType::Pointer pPoint16,
230  typename PointType::Pointer pPoint17,
231  typename PointType::Pointer pPoint18,
232  typename PointType::Pointer pPoint19,
233  typename PointType::Pointer pPoint20 )
234  : BaseType( PointsArrayType(), &msGeometryData )
235  {
236  this->Points().push_back( pPoint1 );
237  this->Points().push_back( pPoint2 );
238  this->Points().push_back( pPoint3 );
239  this->Points().push_back( pPoint4 );
240  this->Points().push_back( pPoint5 );
241  this->Points().push_back( pPoint6 );
242  this->Points().push_back( pPoint7 );
243  this->Points().push_back( pPoint8 );
244  this->Points().push_back( pPoint9 );
245  this->Points().push_back( pPoint10 );
246  this->Points().push_back( pPoint11 );
247  this->Points().push_back( pPoint12 );
248  this->Points().push_back( pPoint13 );
249  this->Points().push_back( pPoint14 );
250  this->Points().push_back( pPoint15 );
251  this->Points().push_back( pPoint16 );
252  this->Points().push_back( pPoint17 );
253  this->Points().push_back( pPoint18 );
254  this->Points().push_back( pPoint19 );
255  this->Points().push_back( pPoint20 );
256  }
257 
258  Hexahedra3D20( const PointsArrayType& ThisPoints )
259  : BaseType( ThisPoints, &msGeometryData )
260  {
261  KRATOS_ERROR_IF(this->PointsNumber() != 20) << "Invalid points number. Expected 20, given " << this->PointsNumber() << std::endl;
262  }
263 
265  explicit Hexahedra3D20(
266  const IndexType GeometryId,
267  const PointsArrayType& rThisPoints
268  ) : BaseType(GeometryId, rThisPoints, &msGeometryData)
269  {
270  KRATOS_ERROR_IF( this->PointsNumber() != 20 ) << "Invalid points number. Expected 20, given " << this->PointsNumber() << std::endl;
271  }
272 
274  explicit Hexahedra3D20(
275  const std::string& GeometryName,
276  const PointsArrayType& rThisPoints
277  ) : BaseType( GeometryName, rThisPoints, &msGeometryData)
278  {
279  KRATOS_ERROR_IF(this->PointsNumber() != 20) << "Invalid points number. Expected 20, given " << this->PointsNumber() << std::endl;
280  }
281 
291  Hexahedra3D20( Hexahedra3D20 const& rOther )
292  : BaseType( rOther )
293  {
294  }
295 
307  template<class TOtherPointType> Hexahedra3D20( Hexahedra3D20<TOtherPointType> const& rOther )
308  : BaseType( rOther )
309  {
310  }
311 
315  ~Hexahedra3D20() override {}
316 
317 
319  {
321  }
322 
324  {
326  }
327 
343  {
344  BaseType::operator=( rOther );
345  return *this;
346  }
347 
359  template<class TOtherPointType>
361  {
362  BaseType::operator=( rOther );
363 
364  return *this;
365  }
366 
370 
377  typename BaseType::Pointer Create(
378  const IndexType NewGeometryId,
379  PointsArrayType const& rThisPoints
380  ) const override
381  {
382  return typename BaseType::Pointer( new Hexahedra3D20( NewGeometryId, rThisPoints ) );
383  }
384 
391  typename BaseType::Pointer Create(
392  const IndexType NewGeometryId,
393  const BaseType& rGeometry
394  ) const override
395  {
396  auto p_geometry = typename BaseType::Pointer( new Hexahedra3D20( NewGeometryId, rGeometry.Points() ) );
397  p_geometry->SetData(rGeometry.GetData());
398  return p_geometry;
399  }
400 
420  double Length() const override
421  {
422  return sqrt( fabs( this->DeterminantOfJacobian( PointType() ) ) );
423  }
424 
438  double Area() const override
439  {
440  return Volume();
441 
442  }
443 
452  double Volume() const override
453  {
455  }
456 
470  double DomainSize() const override
471  {
472  return Volume();
473  }
474 
483  bool IsInside(
484  const CoordinatesArrayType& rPoint,
485  CoordinatesArrayType& rResult,
486  const double Tolerance = std::numeric_limits<double>::epsilon()
487  ) const override
488  {
489  this->PointLocalCoordinates( rResult, rPoint );
490 
491  if ( std::abs( rResult[0] ) <= (1.0 + Tolerance) )
492  {
493  if ( std::abs( rResult[1] ) <= (1.0 + Tolerance) )
494  {
495  if ( std::abs( rResult[2] ) <= (1.0 + Tolerance) )
496  {
497  return true;
498  }
499  }
500  }
501 
502  return false;
503  }
504 
505 
513 
525  SizeType EdgesNumber() const override
526  {
527  return 12;
528  }
529 
539  {
541  typedef typename Geometry<TPointType>::Pointer EdgePointerType;
542  //0
543  edges.push_back( EdgePointerType( new EdgeType(
544  this->pGetPoint( 0 ),
545  this->pGetPoint( 1 ),
546  this->pGetPoint( 8 ) ) ) );
547  edges.push_back( EdgePointerType( new EdgeType(
548  this->pGetPoint( 1 ),
549  this->pGetPoint( 2 ),
550  this->pGetPoint( 9 ) ) ) );
551  edges.push_back( EdgePointerType( new EdgeType(
552  this->pGetPoint( 2 ),
553  this->pGetPoint( 3 ),
554  this->pGetPoint( 10 ) ) ) );
555  edges.push_back( EdgePointerType( new EdgeType(
556  this->pGetPoint( 3 ),
557  this->pGetPoint( 0 ),
558  this->pGetPoint( 11 ) ) ) );
559 
560  //4
561  edges.push_back( EdgePointerType( new EdgeType(
562  this->pGetPoint( 4 ),
563  this->pGetPoint( 5 ),
564  this->pGetPoint( 16 ) ) ) );
565  edges.push_back( EdgePointerType( new EdgeType(
566  this->pGetPoint( 5 ),
567  this->pGetPoint( 6 ),
568  this->pGetPoint( 17 ) ) ) );
569  edges.push_back( EdgePointerType( new EdgeType(
570  this->pGetPoint( 6 ),
571  this->pGetPoint( 7 ),
572  this->pGetPoint( 18 ) ) ) );
573  edges.push_back( EdgePointerType( new EdgeType(
574  this->pGetPoint( 7 ),
575  this->pGetPoint( 4 ),
576  this->pGetPoint( 19 ) ) ) );
577 
578  //8
579  edges.push_back( EdgePointerType( new EdgeType(
580  this->pGetPoint( 0 ),
581  this->pGetPoint( 4 ),
582  this->pGetPoint( 12 ) ) ) );
583  edges.push_back( EdgePointerType( new EdgeType(
584  this->pGetPoint( 1 ),
585  this->pGetPoint( 5 ),
586  this->pGetPoint( 13 ) ) ) );
587  edges.push_back( EdgePointerType( new EdgeType(
588  this->pGetPoint( 2 ),
589  this->pGetPoint( 6 ),
590  this->pGetPoint( 14 ) ) ) );
591  edges.push_back( EdgePointerType( new EdgeType(
592  this->pGetPoint( 3 ),
593  this->pGetPoint( 7 ),
594  this->pGetPoint( 15 ) ) ) );
595  return edges;
596  }
597 
601 
609  SizeType FacesNumber() const override
610  {
611  return 6;
612  }
613 
623  {
625  typedef typename Geometry<TPointType>::Pointer FacePointerType;
626 
627  faces.push_back( FacePointerType( new FaceType(
628  this->pGetPoint( 3 ),
629  this->pGetPoint( 2 ),
630  this->pGetPoint( 1 ),
631  this->pGetPoint( 0 ),
632  this->pGetPoint( 10 ),
633  this->pGetPoint( 9 ),
634  this->pGetPoint( 8 ),
635  this->pGetPoint( 11 ) ) ) );
636  faces.push_back( FacePointerType( new FaceType(
637  this->pGetPoint( 0 ),
638  this->pGetPoint( 1 ),
639  this->pGetPoint( 5 ),
640  this->pGetPoint( 4 ),
641  this->pGetPoint( 8 ),
642  this->pGetPoint( 13 ),
643  this->pGetPoint( 16 ),
644  this->pGetPoint( 12 ) ) ) );
645  faces.push_back( FacePointerType( new FaceType(
646  this->pGetPoint( 2 ),
647  this->pGetPoint( 6 ),
648  this->pGetPoint( 5 ),
649  this->pGetPoint( 1 ),
650  this->pGetPoint( 14 ),
651  this->pGetPoint( 17 ),
652  this->pGetPoint( 13 ),
653  this->pGetPoint( 9 ) ) ) );
654  faces.push_back( FacePointerType( new FaceType(
655  this->pGetPoint( 7 ),
656  this->pGetPoint( 6 ),
657  this->pGetPoint( 2 ),
658  this->pGetPoint( 3 ),
659  this->pGetPoint( 14 ),
660  this->pGetPoint( 18 ),
661  this->pGetPoint( 10 ),
662  this->pGetPoint( 15 ) ) ) );
663  faces.push_back( FacePointerType( new FaceType(
664  this->pGetPoint( 7 ),
665  this->pGetPoint( 3 ),
666  this->pGetPoint( 0 ),
667  this->pGetPoint( 4 ),
668  this->pGetPoint( 15 ),
669  this->pGetPoint( 11 ),
670  this->pGetPoint( 12 ),
671  this->pGetPoint( 19 ) ) ) );
672  faces.push_back( FacePointerType( new FaceType(
673  this->pGetPoint( 4 ),
674  this->pGetPoint( 5 ),
675  this->pGetPoint( 6 ),
676  this->pGetPoint( 7 ),
677  this->pGetPoint( 16 ),
678  this->pGetPoint( 17 ),
679  this->pGetPoint( 18 ),
680  this->pGetPoint( 19 ) ) ) );
681  return faces;
682  }
683 
698  double ShapeFunctionValue( IndexType ShapeFunctionIndex,
699  const CoordinatesArrayType& rPoint ) const override
700  {
701  switch ( ShapeFunctionIndex )
702  {
703  case 0 :
704  return -(( 1.0 + rPoint[0] )
705  *( 1.0 - rPoint[1] )*( 2.0
706  - rPoint[0] + rPoint[1] - rPoint[2] )*( 1.0 + rPoint[2] ) ) / 8.0;
707  case 1 :
708  return -(( 1.0 + rPoint[0] )
709  *( 1.0 + rPoint[1] )*( 2.0
710  - rPoint[0] - rPoint[1] - rPoint[2] )*( 1.0 + rPoint[2] ) ) / 8.0;
711  case 2 :
712  return -(( 1.0 + rPoint[0] )
713  *( 1.0 + rPoint[1] )*( 1.0 - rPoint[2] )*( 2.0
714  - rPoint[0] - rPoint[1] + rPoint[2] ) ) / 8.0;
715  case 3:
716  return -(( 1.0 + rPoint[0] )
717  *( 1.0 - rPoint[1] )*( 1.0 - rPoint[2] )*( 2.0
718  - rPoint[0] + rPoint[1] + rPoint[2] ) ) / 8.0;
719  case 4 :
720  return -(( 1.0 - rPoint[0] )
721  *( 1.0 - rPoint[1] )*( 2.0
722  + rPoint[0] + rPoint[1] - rPoint[2] )*( 1.0 + rPoint[2] ) ) / 8.0;
723  case 5 :
724  return -(( 1.0 - rPoint[0] )
725  *( 1.0 + rPoint[1] )*( 2.0
726  + rPoint[0] - rPoint[1] - rPoint[2] )*( 1.0 + rPoint[2] ) ) / 8.0;
727  case 6 :
728  return -(( 1.0 - rPoint[0] )*( 1.0
729  + rPoint[1] )*( 1.0 - rPoint[2] )*( 2.0
730  + rPoint[0] - rPoint[1] + rPoint[2] ) ) / 8.0;
731  case 7 :
732  return -(( 1.0 - rPoint[0] )*( 1.0 - rPoint[1] )
733  *( 1.0 - rPoint[2] )*( 2.0 + rPoint[0]
734  + rPoint[1] + rPoint[2] ) ) / 8.0;
735  case 8 :
736  return (( 1.0 + rPoint[0] )
737  *( 1.0 - rPoint[1]*rPoint[1] )*( 1.0 + rPoint[2] ) ) / 4.0 ;
738  case 9 :
739  return (( 1.0 + rPoint[0] )*( 1.0 + rPoint[1] )
740  *( 1.0 - rPoint[2]*rPoint[2] ) ) / 4.0 ;
741  case 10 :
742  return (( 1.0 + rPoint[0] )
743  *( 1.0 - rPoint[1]*rPoint[1] )*( 1.0 - rPoint[2] ) ) / 4.0 ;
744  case 11 :
745  return (( 1.0 + rPoint[0] )*( 1.0 - rPoint[1] )
746  *( 1.0 - rPoint[2]*rPoint[2] ) ) / 4.0 ;
747  case 12 :
748  return (( 1.0 -rPoint[0]
749  *rPoint[0] )*( 1.0 - rPoint[1] )*( 1.0 + rPoint[2] ) ) / 4.0 ;
750  case 13 :
751  return (( 1.0 -rPoint[0]
752  *rPoint[0] )*( 1.0 + rPoint[1] )*( 1.0 + rPoint[2] ) ) / 4.0 ;
753  case 14 :
754  return (( 1.0 -rPoint[0]
755  *rPoint[0] )*( 1.0 + rPoint[1] )*( 1.0 - rPoint[2] ) ) / 4.0 ;
756  case 15 :
757  return (( 1.0 -rPoint[0]*rPoint[0] )
758  *( 1.0 - rPoint[1] )*( 1.0 - rPoint[2] ) ) / 4.0;
759  case 16 :
760  return (( 1.0 -rPoint[0] )
761  *( 1.0 - rPoint[1]*rPoint[1] )*( 1.0 + rPoint[2] ) ) / 4.0 ;
762  case 17 :
763  return (( 1.0 -rPoint[0] )*( 1.0 + rPoint[1] )
764  *( 1.0 - rPoint[2]*rPoint[2] ) ) / 4.0 ;
765  case 18 :
766  return (( 1.0 -rPoint[0] )
767  *( 1.0 - rPoint[1]*rPoint[1] )*( 1.0 - rPoint[2] ) ) / 4.0 ;
768  case 19 :
769  return (( 1.0 -rPoint[0] )
770  *( 1.0 - rPoint[1] )*( 1.0 - rPoint[2]*rPoint[2] ) ) / 4.0 ;
771 
772  default:
773  KRATOS_ERROR << "Wrong index of shape function!" << *this << std::endl;
774  }
775 
776  return 0;
777  }
778 
790  Vector& ShapeFunctionsValues (Vector &rResult, const CoordinatesArrayType& rCoordinates) const override
791  {
792  if(rResult.size() != 20) rResult.resize(20,false);
793  rResult[0] = -(( 1.0 + rCoordinates[0] )*( 1.0 - rCoordinates[1] )*( 2.0
794  - rCoordinates[0] + rCoordinates[1] - rCoordinates[2] )*( 1.0 + rCoordinates[2] ) ) / 8.0;
795  rResult[1] = -(( 1.0 + rCoordinates[0] )
796  *( 1.0 + rCoordinates[1] )*( 2.0
797  - rCoordinates[0] - rCoordinates[1] - rCoordinates[2] )*( 1.0 + rCoordinates[2] ) ) / 8.0;
798  rResult[2] = -(( 1.0 + rCoordinates[0] )
799  *( 1.0 + rCoordinates[1] )*( 1.0 - rCoordinates[2] )*( 2.0
800  - rCoordinates[0] - rCoordinates[1] + rCoordinates[2] ) ) / 8.0;
801  rResult[3] = -(( 1.0 + rCoordinates[0] )
802  *( 1.0 - rCoordinates[1] )*( 1.0 - rCoordinates[2] )*( 2.0
803  - rCoordinates[0] + rCoordinates[1] + rCoordinates[2] ) ) / 8.0;
804  rResult[4] = -(( 1.0 - rCoordinates[0] )
805  *( 1.0 - rCoordinates[1] )*( 2.0
806  + rCoordinates[0] + rCoordinates[1] - rCoordinates[2] )*( 1.0 + rCoordinates[2] ) ) / 8.0;
807  rResult[5] = -(( 1.0 - rCoordinates[0] )
808  *( 1.0 + rCoordinates[1] )*( 2.0
809  + rCoordinates[0] - rCoordinates[1] - rCoordinates[2] )*( 1.0 + rCoordinates[2] ) ) / 8.0;
810  rResult[6] = -(( 1.0 - rCoordinates[0] )*( 1.0
811  + rCoordinates[1] )*( 1.0 - rCoordinates[2] )*( 2.0
812  + rCoordinates[0] - rCoordinates[1] + rCoordinates[2] ) ) / 8.0;
813  rResult[7] = -(( 1.0 - rCoordinates[0] )*( 1.0 - rCoordinates[1] )
814  *( 1.0 - rCoordinates[2] )*( 2.0 + rCoordinates[0]
815  + rCoordinates[1] + rCoordinates[2] ) ) / 8.0;
816  rResult[8] = (( 1.0 + rCoordinates[0] )
817  *( 1.0 - rCoordinates[1]*rCoordinates[1] )*( 1.0 + rCoordinates[2] ) ) / 4.0 ;
818  rResult[9] = (( 1.0 + rCoordinates[0] )*( 1.0 + rCoordinates[1] )
819  *( 1.0 - rCoordinates[2]*rCoordinates[2] ) ) / 4.0 ;
820  rResult[10] = (( 1.0 + rCoordinates[0] )
821  *( 1.0 - rCoordinates[1]*rCoordinates[1] )*( 1.0 - rCoordinates[2] ) ) / 4.0 ;
822  rResult[11] = (( 1.0 + rCoordinates[0] )*( 1.0 - rCoordinates[1] )
823  *( 1.0 - rCoordinates[2]*rCoordinates[2] ) ) / 4.0 ;
824  rResult[12] = (( 1.0 -rCoordinates[0]
825  *rCoordinates[0] )*( 1.0 - rCoordinates[1] )*( 1.0 + rCoordinates[2] ) ) / 4.0 ;
826  rResult[13] = (( 1.0 -rCoordinates[0]
827  *rCoordinates[0] )*( 1.0 + rCoordinates[1] )*( 1.0 + rCoordinates[2] ) ) / 4.0 ;
828  rResult[14] = (( 1.0 -rCoordinates[0]
829  *rCoordinates[0] )*( 1.0 + rCoordinates[1] )*( 1.0 - rCoordinates[2] ) ) / 4.0 ;
830  rResult[15] = (( 1.0 -rCoordinates[0]*rCoordinates[0] )
831  *( 1.0 - rCoordinates[1] )*( 1.0 - rCoordinates[2] ) ) / 4.0;
832  rResult[16] = (( 1.0 -rCoordinates[0] )
833  *( 1.0 - rCoordinates[1]*rCoordinates[1] )*( 1.0 + rCoordinates[2] ) ) / 4.0 ;
834  rResult[17] = (( 1.0 -rCoordinates[0] )*( 1.0 + rCoordinates[1] )
835  *( 1.0 - rCoordinates[2]*rCoordinates[2] ) ) / 4.0 ;
836  rResult[18] = (( 1.0 -rCoordinates[0] )
837  *( 1.0 - rCoordinates[1]*rCoordinates[1] )*( 1.0 - rCoordinates[2] ) ) / 4.0 ;
838  rResult[19] = (( 1.0 -rCoordinates[0] )
839  *( 1.0 - rCoordinates[1] )*( 1.0 - rCoordinates[2]*rCoordinates[2] ) ) / 4.0 ;
840  return rResult;
841  }
842 
853  std::string Info() const override
854  {
855  return "3 dimensional hexahedra with 20 nodes and quadratic shape functions in 3D space";
856  }
857 
865  void PrintInfo( std::ostream& rOStream ) const override
866  {
867  rOStream << "3 dimensional hexahedra with 20 nodes and quadratic shape functions in 3D space";
868  }
869 
879  void PrintData( std::ostream& rOStream ) const override
880  {
881  // Base Geometry class PrintData call
882  BaseType::PrintData( rOStream );
883  std::cout << std::endl;
884 
885  // If the geometry has valid points, calculate and output its data
886  if (this->AllPointsAreValid()) {
887  Matrix jacobian;
888  this->Jacobian( jacobian, PointType() );
889  rOStream << " Jacobian in the origin\t : " << jacobian;
890  }
891  }
892 
905  const CoordinatesArrayType& rPoint ) const override
906  {
907  //setting up result matrix
908  if ( result.size1() != 20 || result.size2() != 3 )
909  result.resize( 20, 3, false );
910 
911  result( 0, 0 ) = (( -1.0 + rPoint[1] ) * ( 1.0 - 2.0 * rPoint[0] + rPoint[1] - rPoint[2] ) * ( 1.0
912  + rPoint[2] ) ) / 8.0;
913 
914  result( 0, 1 ) = -(( 1.0 + rPoint[0] ) * ( 1.0 + rPoint[2] ) * ( -1.0 + rPoint[0] - 2.0
915  * rPoint[1] + rPoint[2] ) ) / 8.0;
916 
917  result( 0, 2 ) = -(( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] ) * ( -1.0 + rPoint[0] - rPoint[1] + 2.0
918  * rPoint[2] ) ) / 8.0;
919 
920  result( 1, 0 ) = (( 1.0 + rPoint[1] ) * ( 1.0 + rPoint[2] ) * ( -1.0 + 2.0
921  * rPoint[0] + rPoint[1] + rPoint[2] ) ) / 8.0;
922 
923  result( 1, 1 ) = (( 1.0 + rPoint[0] ) * ( 1.0 + rPoint[2] ) * ( -1.0 + rPoint[0] + 2.0
924  * rPoint[1] + rPoint[2] ) ) / 8.0;
925 
926  result( 1, 2 ) = (( 1.0 + rPoint[0] ) * ( 1.0 + rPoint[1] ) * ( -1.0 + rPoint[0] + rPoint[1] + 2.0
927  * rPoint[2] ) ) / 8.0;
928 
929  result( 2, 0 ) = -(( 1.0 + rPoint[1] ) * ( -1.0 + 2.0 * rPoint[0] + rPoint[1] - rPoint[2] ) * ( -1.0
930  + rPoint[2] ) ) / 8.0;
931 
932  result( 2, 1 ) = -(( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[0] + 2.0 * rPoint[1] - rPoint[2] ) * ( -1.0
933  + rPoint[2] ) ) / 8.0;
934 
935  result( 2, 2 ) = -(( 1.0 + rPoint[0] ) * ( 1.0 + rPoint[1] ) * ( -1.0 + rPoint[0] + rPoint[1] - 2.0
936  * rPoint[2] ) ) / 8.0;
937 
938  result( 3, 0 ) = -(( -1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] ) * ( 1.0 - 2.0
939  * rPoint[0] + rPoint[1] + rPoint[2] ) ) / 8.0;
940 
941  result( 3, 1 ) = (( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[0] - 2.0 * rPoint[1] - rPoint[2] ) * ( -1.0
942  + rPoint[2] ) ) / 8.0;
943 
944  result( 3, 2 ) = (( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] ) * ( -1.0 + rPoint[0] - rPoint[1] - 2.0
945  * rPoint[2] ) ) / 8.0;
946 
947  result( 4, 0 ) = -(( -1.0 + rPoint[1] ) * ( 1.0 + 2.0 * rPoint[0] + rPoint[1] - rPoint[2] ) * ( 1.0
948  + rPoint[2] ) ) / 8.0;
949 
950  result( 4, 1 ) = -(( -1.0 + rPoint[0] ) * ( 1.0 + rPoint[0] + 2.0 * rPoint[1] - rPoint[2] ) * ( 1.0
951  + rPoint[2] ) ) / 8.0;
952 
953  result( 4, 2 ) = -(( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] ) * ( 1.0 + rPoint[0] + rPoint[1] - 2.0
954  * rPoint[2] ) ) / 8.0;
955 
956  result( 5, 0 ) = -(( 1.0 + rPoint[1] ) * ( 1.0 + rPoint[2] ) * ( -1.0 - 2.0
957  * rPoint[0] + rPoint[1] + rPoint[2] ) ) / 8.0;
958 
959  result( 5, 1 ) = (( -1.0 + rPoint[0] ) * ( 1.0 + rPoint[0] - 2.0 * rPoint[1] - rPoint[2] ) * ( 1.0
960  + rPoint[2] ) ) / 8.0;
961 
962  result( 5, 2 ) = (( -1.0 + rPoint[0] ) * ( 1.0 + rPoint[1] ) * ( 1.0 + rPoint[0] - rPoint[1] - 2.0
963  * rPoint[2] ) ) / 8.0;
964 
965  result( 6, 0 ) = (( 1.0 + rPoint[1] ) * ( -1.0 - 2.0 * rPoint[0] + rPoint[1] - rPoint[2] ) * ( -1.0
966  + rPoint[2] ) ) / 8.0;
967 
968  result( 6, 1 ) = -(( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[2] ) * ( 1.0 + rPoint[0] - 2.0
969  * rPoint[1] + rPoint[2] ) ) / 8.0;
970 
971  result( 6, 2 ) = -(( -1.0 + rPoint[0] ) * ( 1.0 + rPoint[1] ) * ( 1.0 + rPoint[0] - rPoint[1] + 2.0
972  * rPoint[2] ) ) / 8.0;
973 
974  result( 7, 0 ) = (( -1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] ) * ( 1.0 + 2.0
975  * rPoint[0] + rPoint[1] + rPoint[2] ) ) / 8.0;
976 
977  result( 7, 1 ) = (( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[2] ) * ( 1.0 + rPoint[0] + 2.0
978  * rPoint[1] + rPoint[2] ) ) / 8.0;
979 
980  result( 7, 2 ) = (( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] ) * ( 1.0 + rPoint[0] + rPoint[1] + 2.0
981  * rPoint[2] ) ) / 8.0;
982 
983  result( 8, 0 ) = -(( -1.0 + rPoint[1] * rPoint[1] ) * ( 1.0 + rPoint[2] ) ) / 4.0;
984 
985  result( 8, 1 ) = -(( 1.0 + rPoint[0] ) * rPoint[1] * ( 1.0 + rPoint[2] ) ) / 2.0;
986 
987  result( 8, 2 ) = -(( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] * rPoint[1] ) ) / 4.0;
988 
989  result( 9, 0 ) = -(( 1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
990 
991  result( 9, 1 ) = -(( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
992 
993  result( 9, 2 ) = -(( 1.0 + rPoint[0] ) * ( 1.0 + rPoint[1] ) * rPoint[2] ) / 2.0;
994 
995  result( 10, 0 ) = (( -1.0 + rPoint[1] * rPoint[1] ) * ( -1.0 + rPoint[2] ) ) / 4.0;
996 
997  result( 10, 1 ) = (( 1.0 + rPoint[0] ) * rPoint[1] * ( -1.0 + rPoint[2] ) ) / 2.0;
998 
999  result( 10, 2 ) = (( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] * rPoint[1] ) ) / 4.0;
1000 
1001  result( 11, 0 ) = (( -1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
1002 
1003  result( 11, 1 ) = (( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
1004 
1005  result( 11, 2 ) = (( 1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] ) * rPoint[2] ) / 2.0;
1006 
1007  result( 12, 0 ) = ( rPoint[0] * ( -1.0 + rPoint[1] ) * ( 1.0 + rPoint[2] ) ) / 2.0;
1008 
1009  result( 12, 1 ) = (( -1.0 + rPoint[0] * rPoint[0] ) * ( 1.0 + rPoint[2] ) ) / 4.0;
1010 
1011  result( 12, 2 ) = (( -1.0 + rPoint[0] * rPoint[0] ) * ( -1.0 + rPoint[1] ) ) / 4.0;
1012 
1013  result( 13, 0 ) = -( rPoint[0] * ( 1.0 + rPoint[1] ) * ( 1.0 + rPoint[2] ) ) / 2.0;
1014 
1015  result( 13, 1 ) = -(( -1.0 + rPoint[0] * rPoint[0] ) * ( 1.0 + rPoint[2] ) ) / 4.0;
1016 
1017  result( 13, 2 ) = -(( -1.0 + rPoint[0] * rPoint[0] ) * ( 1.0 + rPoint[1] ) ) / 4.0;
1018 
1019  result( 14, 0 ) = ( rPoint[0] * ( 1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] ) ) / 2.0;
1020 
1021  result( 14, 1 ) = (( -1.0 + rPoint[0] * rPoint[0] ) * ( -1.0 + rPoint[2] ) ) / 4.0;
1022 
1023  result( 14, 2 ) = (( -1.0 + rPoint[0] * rPoint[0] ) * ( 1.0 + rPoint[1] ) ) / 4.0;
1024 
1025  result( 15, 0 ) = -( rPoint[0] * ( -1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] ) ) / 2.0;
1026 
1027  result( 15, 1 ) = -(( -1.0 + rPoint[0] * rPoint[0] ) * ( -1.0 + rPoint[2] ) ) / 4.0;
1028 
1029  result( 15, 2 ) = -(( -1.0 + rPoint[0] * rPoint[0] ) * ( -1.0 + rPoint[1] ) ) / 4.0;
1030 
1031  result( 16, 0 ) = (( -1.0 + rPoint[1] * rPoint[1] ) * ( 1.0 + rPoint[2] ) ) / 4.0;
1032 
1033  result( 16, 1 ) = (( -1.0 + rPoint[0] ) * rPoint[1] * ( 1.0 + rPoint[2] ) ) / 2.0;
1034 
1035  result( 16, 2 ) = (( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] * rPoint[1] ) ) / 4.0;
1036 
1037  result( 17, 0 ) = (( 1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
1038 
1039  result( 17, 1 ) = (( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
1040 
1041  result( 17, 2 ) = (( -1.0 + rPoint[0] ) * ( 1.0 + rPoint[1] ) * rPoint[2] ) / 2.0;
1042 
1043  result( 18, 0 ) = -(( -1.0 + rPoint[1] * rPoint[1] ) * ( -1.0 + rPoint[2] ) ) / 4.0;
1044 
1045  result( 18, 1 ) = -(( -1.0 + rPoint[0] ) * rPoint[1] * ( -1.0 + rPoint[2] ) ) / 2.0;
1046 
1047  result( 18, 2 ) = -(( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] * rPoint[1] ) ) / 4.0;
1048 
1049  result( 19, 0 ) = -(( -1.0 + rPoint[1] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
1050 
1051  result( 19, 1 ) = -(( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[2] * rPoint[2] ) ) / 4.0;
1052 
1053  result( 19, 2 ) = -(( -1.0 + rPoint[0] ) * ( -1.0 + rPoint[1] ) * rPoint[2] ) / 2.0;
1054 
1055  return( result );
1056  }
1057 
1058 
1059 
1060 protected:
1061 
1062 
1067 private:
1068 
1073  static const GeometryData msGeometryData;
1074 
1075  static const GeometryDimension msGeometryDimension;
1076 
1080 
1081  friend class Serializer;
1082 
1083  void save( Serializer& rSerializer ) const override
1084  {
1085 
1087  }
1088 
1089  void load( Serializer& rSerializer ) override
1090  {
1092  }
1093 
1094  Hexahedra3D20(): BaseType( PointsArrayType(), &msGeometryData ) {}
1095 
1113  static Matrix CalculateShapeFunctionsIntegrationPointsValues(
1114  typename BaseType::IntegrationMethod ThisMethod )
1115  {
1116  const IntegrationPointsContainerType all_integration_points =
1117  AllIntegrationPoints();
1118  const IntegrationPointsArrayType integration_points = all_integration_points[static_cast<int>(ThisMethod)];
1119  //number of integration points
1120  const int integration_points_number = integration_points.size();
1121  //number of nodes in current geometry
1122  const int points_number = 20;
1123  //setting up return matrix
1124  Matrix shape_function_values( integration_points_number, points_number );
1125  //loop over all integration points
1126 
1127  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
1128  {
1129  shape_function_values(pnt, 0 ) =
1130  -(( 1.0 + integration_points[pnt].X() )
1131  * ( 1.0 - integration_points[pnt].Y() ) * ( 2.0
1132  - integration_points[pnt].X() + integration_points[pnt].Y()
1133  - integration_points[pnt].Z() ) * ( 1.0 + integration_points[pnt].Z() ) ) / 8.0;
1134  shape_function_values(pnt, 1 ) =
1135  -(( 1.0 + integration_points[pnt].X() )
1136  * ( 1.0 + integration_points[pnt].Y() ) * ( 2.0
1137  - integration_points[pnt].X() - integration_points[pnt].Y()
1138  - integration_points[pnt].Z() ) * ( 1.0 + integration_points[pnt].Z() ) ) / 8.0;
1139  shape_function_values(pnt, 2 ) =
1140  -(( 1.0 + integration_points[pnt].X() )
1141  * ( 1.0 + integration_points[pnt].Y() ) * ( 1.0 - integration_points[pnt].Z() ) * ( 2.0
1142  - integration_points[pnt].X() - integration_points[pnt].Y()
1143  + integration_points[pnt].Z() ) ) / 8.0;
1144  shape_function_values(pnt, 3 ) =
1145  -(( 1.0 + integration_points[pnt].X() )
1146  * ( 1.0 - integration_points[pnt].Y() ) * ( 1.0 - integration_points[pnt].Z() )
1147  * ( 2.0 - integration_points[pnt].X()
1148  + integration_points[pnt].Y() + integration_points[pnt].Z() ) ) / 8.0;
1149  shape_function_values(pnt, 4 ) =
1150  -(( 1.0 - integration_points[pnt].X() )
1151  * ( 1.0 - integration_points[pnt].Y() ) * ( 2.0
1152  + integration_points[pnt].X() + integration_points[pnt].Y()
1153  - integration_points[pnt].Z() ) * ( 1.0 + integration_points[pnt].Z() ) ) / 8.0;
1154  shape_function_values(pnt, 5 ) =
1155  -(( 1.0 - integration_points[pnt].X() )
1156  * ( 1.0 + integration_points[pnt].Y() ) * ( 2.0
1157  + integration_points[pnt].X() - integration_points[pnt].Y()
1158  - integration_points[pnt].Z() ) * ( 1.0 + integration_points[pnt].Z() ) ) / 8.0;
1159  shape_function_values(pnt, 6 ) =
1160  -(( 1.0 - integration_points[pnt].X() ) * ( 1.0
1161  + integration_points[pnt].Y() ) * ( 1.0 - integration_points[pnt].Z() ) * ( 2.0
1162  + integration_points[pnt].X() - integration_points[pnt].Y()
1163  + integration_points[pnt].Z() ) ) / 8.0;
1164  shape_function_values(pnt, 7 ) =
1165  -(( 1.0 - integration_points[pnt].X() ) * ( 1.0 - integration_points[pnt].Y() )
1166  * ( 1.0 - integration_points[pnt].Z() ) * ( 2.0 + integration_points[pnt].X()
1167  + integration_points[pnt].Y() + integration_points[pnt].Z() ) ) / 8.0;
1168  shape_function_values(pnt, 8 ) =
1169  (( 1.0 + integration_points[pnt].X() )
1170  * ( 1.0 - integration_points[pnt].Y() * integration_points[pnt].Y() )
1171  * ( 1.0 + integration_points[pnt].Z() ) ) / 4.0 ;
1172  shape_function_values(pnt, 9 ) =
1173  (( 1.0 + integration_points[pnt].X() ) * ( 1.0 + integration_points[pnt].Y() )
1174  * ( 1.0 - integration_points[pnt].Z() * integration_points[pnt].Z() ) ) / 4.0 ;
1175  shape_function_values(pnt, 10 ) =
1176  (( 1.0 + integration_points[pnt].X() )
1177  * ( 1.0 - integration_points[pnt].Y() * integration_points[pnt].Y() )
1178  * ( 1.0 - integration_points[pnt].Z() ) ) / 4.0 ;
1179  shape_function_values(pnt, 11 ) =
1180  (( 1.0 + integration_points[pnt].X() ) * ( 1.0 - integration_points[pnt].Y() )
1181  * ( 1.0 - integration_points[pnt].Z() * integration_points[pnt].Z() ) ) / 4.0 ;
1182  shape_function_values(pnt, 12 ) =
1183  (( 1.0 - integration_points[pnt].X()
1184  * integration_points[pnt].X() ) * ( 1.0 - integration_points[pnt].Y() )
1185  * ( 1.0 + integration_points[pnt].Z() ) ) / 4.0 ;
1186  shape_function_values(pnt, 13 ) =
1187  (( 1.0 - integration_points[pnt].X()
1188  * integration_points[pnt].X() ) * ( 1.0 + integration_points[pnt].Y() )
1189  * ( 1.0 + integration_points[pnt].Z() ) ) / 4.0 ;
1190  shape_function_values(pnt, 14 ) =
1191  (( 1.0 - integration_points[pnt].X()
1192  * integration_points[pnt].X() ) * ( 1.0 + integration_points[pnt].Y() )
1193  * ( 1.0 - integration_points[pnt].Z() ) ) / 4.0 ;
1194  shape_function_values(pnt, 15 ) =
1195  (( 1.0 - integration_points[pnt].X() * integration_points[pnt].X() )
1196  * ( 1.0 - integration_points[pnt].Y() ) * ( 1.0 - integration_points[pnt].Z() ) ) / 4.0;
1197  shape_function_values(pnt, 16 ) =
1198  (( 1.0 - integration_points[pnt].X() )
1199  * ( 1.0 - integration_points[pnt].Y() * integration_points[pnt].Y() )
1200  * ( 1.0 + integration_points[pnt].Z() ) ) / 4.0 ;
1201  shape_function_values(pnt, 17 ) =
1202  (( 1.0 - integration_points[pnt].X() ) * ( 1.0 + integration_points[pnt].Y() )
1203  * ( 1.0 - integration_points[pnt].Z() * integration_points[pnt].Z() ) ) / 4.0 ;
1204  shape_function_values(pnt, 18 ) =
1205  (( 1.0 - integration_points[pnt].X() )
1206  * ( 1.0 - integration_points[pnt].Y() * integration_points[pnt].Y() )
1207  * ( 1.0 - integration_points[pnt].Z() ) ) / 4.0 ;
1208  shape_function_values(pnt, 19 ) =
1209  (( 1.0 - integration_points[pnt].X() )
1210  * ( 1.0 - integration_points[pnt].Y() ) * ( 1.0
1211  - integration_points[pnt].Z() * integration_points[pnt].Z() ) ) / 4.0 ;
1212  }
1213 
1214  return shape_function_values;
1215  }
1216 
1231  CalculateShapeFunctionsIntegrationPointsLocalGradients(
1232  typename BaseType::IntegrationMethod ThisMethod )
1233  {
1234  const IntegrationPointsContainerType all_integration_points =
1235  AllIntegrationPoints();
1236  const IntegrationPointsArrayType integration_points = all_integration_points[static_cast<int>(ThisMethod)];
1237  //number of integration points
1238  const int integration_points_number = integration_points.size();
1239  ShapeFunctionsGradientsType d_shape_f_values( integration_points_number );
1240  //initialising container
1241  //std::fill(d_shape_f_values.begin(), d_shape_f_values.end(), Matrix(8,3));
1242  //loop over all integration points
1243 
1244  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
1245  {
1246  Matrix result = ZeroMatrix( 20, 3 );
1247 
1248  result( 0, 0 ) = (( -1.0 + integration_points[pnt].Y() )
1249  * ( 1.0 - 2.0 * integration_points[pnt].X()
1250  + integration_points[pnt].Y()
1251  - integration_points[pnt].Z() ) * ( 1.0
1252  + integration_points[pnt].Z() ) ) / 8.0;
1253  result( 0, 1 ) = -(( 1.0 + integration_points[pnt].X() )
1254  * ( 1.0 + integration_points[pnt].Z() ) * ( -1.0
1255  + integration_points[pnt].X() - 2.0
1256  * integration_points[pnt].Y()
1257  + integration_points[pnt].Z() ) ) / 8.0;
1258  result( 0, 2 ) = -(( 1.0 + integration_points[pnt].X() )
1259  * ( -1.0 + integration_points[pnt].Y() ) * ( -1.0
1260  + integration_points[pnt].X()
1261  - integration_points[pnt].Y() + 2.0
1262  * integration_points[pnt].Z() ) ) / 8.0;
1263 
1264  result( 1, 0 ) = (( 1.0 + integration_points[pnt].Y() )
1265  * ( 1.0 + integration_points[pnt].Z() )
1266  * ( -1.0 + 2.0 * integration_points[pnt].X()
1267  + integration_points[pnt].Y()
1268  + integration_points[pnt].Z() ) ) / 8.0;
1269  result( 1, 1 ) = (( 1.0 + integration_points[pnt].X() )
1270  * ( 1.0 + integration_points[pnt].Z() ) * ( -1.0
1271  + integration_points[pnt].X() + 2.0
1272  * integration_points[pnt].Y()
1273  + integration_points[pnt].Z() ) ) / 8.0;
1274  result( 1, 2 ) = (( 1.0 + integration_points[pnt].X() )
1275  * ( 1.0 + integration_points[pnt].Y() ) * ( -1.0
1276  + integration_points[pnt].X()
1277  + integration_points[pnt].Y() + 2.0
1278  * integration_points[pnt].Z() ) ) / 8.0;
1279 
1280  result( 2, 0 ) = -(( 1.0 + integration_points[pnt].Y() )
1281  * ( -1.0 + 2.0 * integration_points[pnt].X()
1282  + integration_points[pnt].Y()
1283  - integration_points[pnt].Z() ) * ( -1.0 + integration_points[pnt].Z() ) ) / 8.0;
1284  result( 2, 1 ) = -(( 1.0 + integration_points[pnt].X() )
1285  * ( -1.0 + integration_points[pnt].X() + 2.0
1286  * integration_points[pnt].Y()
1287  - integration_points[pnt].Z() ) * ( -1.0
1288  + integration_points[pnt].Z() ) ) / 8.0;
1289  result( 2, 2 ) = -(( 1.0 + integration_points[pnt].X() )
1290  * ( 1.0 + integration_points[pnt].Y() )
1291  * ( -1.0 + integration_points[pnt].X()
1292  + integration_points[pnt].Y() - 2.0
1293  * integration_points[pnt].Z() ) ) / 8.0;
1294 
1295  result( 3, 0 ) = -(( -1.0 + integration_points[pnt].Y() )
1296  * ( -1.0 + integration_points[pnt].Z() ) * ( 1.0 - 2.0 * integration_points[pnt].X()
1297  + integration_points[pnt].Y()
1298  + integration_points[pnt].Z() ) ) / 8.0;
1299  result( 3, 1 ) = (( 1.0 + integration_points[pnt].X() )
1300  * ( -1.0 + integration_points[pnt].X() - 2.0
1301  * integration_points[pnt].Y()
1302  - integration_points[pnt].Z() ) * ( -1.0
1303  + integration_points[pnt].Z() ) ) / 8.0;
1304  result( 3, 2 ) = (( 1.0 + integration_points[pnt].X() )
1305  * ( -1.0 + integration_points[pnt].Y() )
1306  * ( -1.0 + integration_points[pnt].X()
1307  - integration_points[pnt].Y() - 2.0
1308  * integration_points[pnt].Z() ) ) / 8.0;
1309 
1310  result( 4, 0 ) = -(( -1.0 + integration_points[pnt].Y() )
1311  * ( 1.0 + 2.0 * integration_points[pnt].X()
1312  + integration_points[pnt].Y()
1313  - integration_points[pnt].Z() ) * ( 1.0 + integration_points[pnt].Z() ) ) / 8.0;
1314  result( 4, 1 ) = -(( -1.0 + integration_points[pnt].X() )
1315  * ( 1.0 + integration_points[pnt].X() + 2.0
1316  * integration_points[pnt].Y()
1317  - integration_points[pnt].Z() ) * ( 1.0
1318  + integration_points[pnt].Z() ) ) / 8.0;
1319  result( 4, 2 ) = -(( -1.0 + integration_points[pnt].X() )
1320  * ( -1.0 + integration_points[pnt].Y() ) * ( 1.0
1321  + integration_points[pnt].X()
1322  + integration_points[pnt].Y() - 2.0
1323  * integration_points[pnt].Z() ) ) / 8.0;
1324 
1325  result( 5, 0 ) = -(( 1.0 + integration_points[pnt].Y() )
1326  * ( 1.0 + integration_points[pnt].Z() ) * ( -1.0 - 2.0
1327  * integration_points[pnt].X()
1328  + integration_points[pnt].Y()
1329  + integration_points[pnt].Z() ) ) / 8.0;
1330  result( 5, 1 ) = (( -1.0 + integration_points[pnt].X() )
1331  * ( 1.0 + integration_points[pnt].X() - 2.0
1332  * integration_points[pnt].Y()
1333  - integration_points[pnt].Z() ) * ( 1.0
1334  + integration_points[pnt].Z() ) ) / 8.0;
1335  result( 5, 2 ) = (( -1.0 + integration_points[pnt].X() )
1336  * ( 1.0 + integration_points[pnt].Y() ) * ( 1.0
1337  + integration_points[pnt].X()
1338  - integration_points[pnt].Y() - 2.0
1339  * integration_points[pnt].Z() ) ) / 8.0;
1340 
1341  result( 6, 0 ) = (( 1.0 + integration_points[pnt].Y() )
1342  * ( -1.0 - 2.0 * integration_points[pnt].X()
1343  + integration_points[pnt].Y()
1344  - integration_points[pnt].Z() ) * ( -1.0
1345  + integration_points[pnt].Z() ) ) / 8.0;
1346  result( 6, 1 ) = -(( -1.0 + integration_points[pnt].X() )
1347  * ( -1.0 + integration_points[pnt].Z() )
1348  * ( 1.0 + integration_points[pnt].X()
1349  - 2.0 * integration_points[pnt].Y()
1350  + integration_points[pnt].Z() ) ) / 8.0;
1351  result( 6, 2 ) = -(( -1.0 + integration_points[pnt].X() )
1352  * ( 1.0 + integration_points[pnt].Y() ) * ( 1.0
1353  + integration_points[pnt].X()
1354  - integration_points[pnt].Y() + 2.0
1355  * integration_points[pnt].Z() ) ) / 8.0;
1356 
1357  result( 7, 0 ) = (( -1.0 + integration_points[pnt].Y() )
1358  * ( -1.0 + integration_points[pnt].Z() ) * ( 1.0 + 2.0
1359  * integration_points[pnt].X() + integration_points[pnt].Y()
1360  + integration_points[pnt].Z() ) ) / 8.0;
1361  result( 7, 1 ) = (( -1.0 + integration_points[pnt].X() )
1362  * ( -1.0 + integration_points[pnt].Z() )
1363  * ( 1.0 + integration_points[pnt].X() + 2.0 * integration_points[pnt].Y()
1364  + integration_points[pnt].Z() ) ) / 8.0;
1365  result( 7, 2 ) = (( -1.0 + integration_points[pnt].X() )
1366  * ( -1.0 + integration_points[pnt].Y() ) * ( 1.0 + integration_points[pnt].X()
1367  + integration_points[pnt].Y() + 2.0
1368  * integration_points[pnt].Z() ) ) / 8.0;
1369 
1370  result( 8, 0 ) = -(( -1.0 + integration_points[pnt].Y()
1371  * integration_points[pnt].Y() )
1372  * ( 1.0 + integration_points[pnt].Z() ) ) / 4.0;
1373  result( 8, 1 ) = -(( 1.0 + integration_points[pnt].X() )
1374  * integration_points[pnt].Y()
1375  * ( 1.0 + integration_points[pnt].Z() ) ) / 2.0;
1376  result( 8, 2 ) = -(( 1.0 + integration_points[pnt].X() )
1377  * ( -1.0 + integration_points[pnt].Y()
1378  * integration_points[pnt].Y() ) ) / 4.0;
1379 
1380  result( 9, 0 ) = -(( 1.0 + integration_points[pnt].Y() )
1381  * ( -1.0 + integration_points[pnt].Z()
1382  * integration_points[pnt].Z() ) ) / 4.0;
1383  result( 9, 1 ) = -(( 1.0 + integration_points[pnt].X() )
1384  * ( -1.0 + integration_points[pnt].Z()
1385  * integration_points[pnt].Z() ) ) / 4.0;
1386  result( 9, 2 ) = -(( 1.0 + integration_points[pnt].X() )
1387  * ( 1.0 + integration_points[pnt].Y() )
1388  * integration_points[pnt].Z() ) / 2.0;
1389 
1390  result( 10, 0 ) = (( -1.0 + integration_points[pnt].Y()
1391  * integration_points[pnt].Y() )
1392  * ( -1.0 + integration_points[pnt].Z() ) ) / 4.0;
1393  result( 10, 1 ) = (( 1.0 + integration_points[pnt].X() )
1394  * integration_points[pnt].Y()
1395  * ( -1.0 + integration_points[pnt].Z() ) ) / 2.0;
1396  result( 10, 2 ) = (( 1.0 + integration_points[pnt].X() )
1397  * ( -1.0 + integration_points[pnt].Y()
1398  * integration_points[pnt].Y() ) ) / 4.0;
1399 
1400  result( 11, 0 ) = (( -1.0 + integration_points[pnt].Y() )
1401  * ( -1.0 + integration_points[pnt].Z()
1402  * integration_points[pnt].Z() ) ) / 4.0;
1403  result( 11, 1 ) = (( 1.0 + integration_points[pnt].X() )
1404  * ( -1.0 + integration_points[pnt].Z()
1405  * integration_points[pnt].Z() ) ) / 4.0;
1406  result( 11, 2 ) = (( 1.0 + integration_points[pnt].X() )
1407  * ( -1.0 + integration_points[pnt].Y() )
1408  * integration_points[pnt].Z() ) / 2.0;
1409 
1410  result( 12, 0 ) = ( integration_points[pnt].X()
1411  * ( -1.0 + integration_points[pnt].Y() )
1412  * ( 1.0 + integration_points[pnt].Z() ) ) / 2.0;
1413  result( 12, 1 ) = (( -1.0 + integration_points[pnt].X()
1414  * integration_points[pnt].X() )
1415  * ( 1.0 + integration_points[pnt].Z() ) ) / 4.0;
1416  result( 12, 2 ) = (( -1.0 + integration_points[pnt].X()
1417  * integration_points[pnt].X() )
1418  * ( -1.0 + integration_points[pnt].Y() ) ) / 4.0;
1419 
1420  result( 13, 0 ) = -( integration_points[pnt].X()
1421  * ( 1.0 + integration_points[pnt].Y() )
1422  * ( 1.0 + integration_points[pnt].Z() ) ) / 2.0;
1423  result( 13, 1 ) = -(( -1.0 + integration_points[pnt].X()
1424  * integration_points[pnt].X() )
1425  * ( 1.0 + integration_points[pnt].Z() ) ) / 4.0;
1426  result( 13, 2 ) = -(( -1.0 + integration_points[pnt].X()
1427  * integration_points[pnt].X() )
1428  * ( 1.0 + integration_points[pnt].Y() ) ) / 4.0;
1429 
1430  result( 14, 0 ) = ( integration_points[pnt].X()
1431  * ( 1.0 + integration_points[pnt].Y() )
1432  * ( -1.0 + integration_points[pnt].Z() ) ) / 2.0;
1433  result( 14, 1 ) = (( -1.0 + integration_points[pnt].X()
1434  * integration_points[pnt].X() )
1435  * ( -1.0 + integration_points[pnt].Z() ) ) / 4.0;
1436  result( 14, 2 ) = (( -1.0 + integration_points[pnt].X()
1437  * integration_points[pnt].X() )
1438  * ( 1.0 + integration_points[pnt].Y() ) ) / 4.0;
1439 
1440  result( 15, 0 ) = -( integration_points[pnt].X() * ( -1.0
1441  + integration_points[pnt].Y() ) * ( -1.0
1442  + integration_points[pnt].Z() ) ) / 2.0;
1443  result( 15, 1 ) = -(( -1.0 + integration_points[pnt].X()
1444  * integration_points[pnt].X() ) * ( -1.0
1445  + integration_points[pnt].Z() ) ) / 4.0;
1446  result( 15, 2 ) = -(( -1.0 + integration_points[pnt].X()
1447  * integration_points[pnt].X() )
1448  * ( -1.0 + integration_points[pnt].Y() ) ) / 4.0;
1449 
1450  result( 16, 0 ) = (( -1.0 + integration_points[pnt].Y()
1451  * integration_points[pnt].Y() )
1452  * ( 1.0 + integration_points[pnt].Z() ) ) / 4.0;
1453  result( 16, 1 ) = (( -1.0 + integration_points[pnt].X() )
1454  * integration_points[pnt].Y()
1455  * ( 1.0 + integration_points[pnt].Z() ) ) / 2.0;
1456  result( 16, 2 ) = (( -1.0 + integration_points[pnt].X() )
1457  * ( -1.0 + integration_points[pnt].Y()
1458  * integration_points[pnt].Y() ) ) / 4.0;
1459 
1460  result( 17, 0 ) = (( 1.0 + integration_points[pnt].Y() )
1461  * ( -1.0 + integration_points[pnt].Z()
1462  * integration_points[pnt].Z() ) ) / 4.0;
1463  result( 17, 1 ) = (( -1.0 + integration_points[pnt].X() )
1464  * ( -1.0 + integration_points[pnt].Z()
1465  * integration_points[pnt].Z() ) ) / 4.0;
1466  result( 17, 2 ) = (( -1.0 + integration_points[pnt].X() )
1467  * ( 1.0 + integration_points[pnt].Y() )
1468  * integration_points[pnt].Z() ) / 2.0;
1469 
1470  result( 18, 0 ) = -(( -1.0 + integration_points[pnt].Y()
1471  * integration_points[pnt].Y() )
1472  * ( -1.0 + integration_points[pnt].Z() ) ) / 4.0;
1473  result( 18, 1 ) = -(( -1.0 + integration_points[pnt].X() )
1474  * integration_points[pnt].Y()
1475  * ( -1.0 + integration_points[pnt].Z() ) ) / 2.0;
1476  result( 18, 2 ) = -(( -1.0 + integration_points[pnt].X() )
1477  * ( -1.0 + integration_points[pnt].Y()
1478  * integration_points[pnt].Y() ) ) / 4.0;
1479 
1480  result( 19, 0 ) = -(( -1.0 + integration_points[pnt].Y() )
1481  * ( -1.0 + integration_points[pnt].Z()
1482  * integration_points[pnt].Z() ) ) / 4.0;
1483  result( 19, 1 ) = -(( -1.0 + integration_points[pnt].X() )
1484  * ( -1.0 + integration_points[pnt].Z()
1485  * integration_points[pnt].Z() ) ) / 4.0;
1486  result( 19, 2 ) = -(( -1.0 + integration_points[pnt].X() )
1487  * ( -1.0 + integration_points[pnt].Y() )
1488  * integration_points[pnt].Z() ) / 2.0;
1489 
1490  d_shape_f_values[pnt] = result;
1491  }
1492 
1493  return d_shape_f_values;
1494  }
1495 
1499  static const IntegrationPointsContainerType AllIntegrationPoints()
1500  {
1501  IntegrationPointsContainerType integration_points =
1502  {
1503  {
1504  Quadrature < HexahedronGaussLegendreIntegrationPoints1,
1505  3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1506  Quadrature < HexahedronGaussLegendreIntegrationPoints2,
1507  3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1508  Quadrature < HexahedronGaussLegendreIntegrationPoints3,
1509  3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1510  Quadrature < HexahedronGaussLegendreIntegrationPoints4,
1511  3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1512  Quadrature < HexahedronGaussLegendreIntegrationPoints5,
1513  3, IntegrationPoint<3> >::GenerateIntegrationPoints()
1514  }
1515  };
1516  return integration_points;
1517  }
1518 
1522  static const ShapeFunctionsValuesContainerType AllShapeFunctionsValues()
1523  {
1524  ShapeFunctionsValuesContainerType shape_functions_values =
1525  {
1526  {
1527  Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1529  Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1531  Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1533  Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1535  Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1537  }
1538  };
1539  return shape_functions_values;
1540  }
1541 
1546  AllShapeFunctionsLocalGradients()
1547  {
1548  ShapeFunctionsLocalGradientsContainerType shape_functions_local_gradients =
1549  {
1550  {
1551  Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1553  Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1555  Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1557  Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1559  Hexahedra3D20<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients(
1561  }
1562  };
1563  return shape_functions_local_gradients;
1564  }
1565 
1566 
1571  template<class TOtherPointType> friend class Hexahedra3D20;
1572 
1573 
1578 };// Class Hexahedra3D20
1579 
1580 
1588 template<class TPointType> inline std::istream& operator >> (
1589  std::istream& rIStream, Hexahedra3D20<TPointType>& rThis );
1590 
1594 template<class TPointType> inline std::ostream& operator << (
1595  std::ostream& rOStream, const Hexahedra3D20<TPointType>& rThis )
1596 {
1597  rThis.PrintInfo( rOStream );
1598  rOStream << std::endl;
1599  rThis.PrintData( rOStream );
1600 
1601  return rOStream;
1602 }
1603 
1604 template<class TPointType> const
1605 GeometryData Hexahedra3D20<TPointType>::msGeometryData(
1608  Hexahedra3D20<TPointType>::AllIntegrationPoints(),
1609  Hexahedra3D20<TPointType>::AllShapeFunctionsValues(),
1610  AllShapeFunctionsLocalGradients()
1611 );
1612 
1613 template<class TPointType> const
1614 GeometryDimension Hexahedra3D20<TPointType>::msGeometryDimension(3, 3);
1615 
1616 }// namespace Kratos.
Definition: geometry_data.h:60
KratosGeometryType
Definition: geometry_data.h:110
IntegrationMethod
Definition: geometry_data.h:76
KratosGeometryFamily
Definition: geometry_data.h:91
Definition: geometry_dimension.h:42
Geometry base class.
Definition: geometry.h:71
SizeType PointsNumber() const
Definition: geometry.h:528
virtual void PrintData(std::ostream &rOStream) const
Print object's data.
Definition: geometry.h:3834
Geometry & operator=(const Geometry &rOther)
Definition: geometry.h:400
std::vector< IntegrationPointType > IntegrationPointsArrayType
Definition: geometry.h:161
DataValueContainer & GetData()
Definition: geometry.h:591
std::size_t SizeType
Definition: geometry.h:144
const TPointType::Pointer pGetPoint(const int Index) const
Definition: geometry.h:1790
std::size_t IndexType
Definition: geometry.h:137
std::array< Matrix, static_cast< int >GeometryData::IntegrationMethod::NumberOfIntegrationMethods)> ShapeFunctionsValuesContainerType
Definition: geometry.h:172
Vector & DeterminantOfJacobian(Vector &rResult) const
Definition: geometry.h:3154
virtual CoordinatesArrayType & PointLocalCoordinates(CoordinatesArrayType &rResult, const CoordinatesArrayType &rPoint) const
Returns the local coordinates of a given arbitrary point.
Definition: geometry.h:1854
const PointsArrayType & Points() const
Definition: geometry.h:1768
bool AllPointsAreValid() const
Checks if the geometry points are valid Checks if the geometry points are valid from the pointer valu...
Definition: geometry.h:4093
std::array< IntegrationPointsArrayType, static_cast< int >GeometryData::IntegrationMethod::NumberOfIntegrationMethods)> IntegrationPointsContainerType
Definition: geometry.h:167
JacobiansType & Jacobian(JacobiansType &rResult) const
Definition: geometry.h:2901
GeometryData::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: geometry.h:177
GeometryData::IntegrationMethod IntegrationMethod
Definition: geometry.h:122
A twenty node hexahedra geometry with serendipity shape functions.
Definition: hexahedra_3d_20.h:51
double ShapeFunctionValue(IndexType ShapeFunctionIndex, const CoordinatesArrayType &rPoint) const override
Definition: hexahedra_3d_20.h:698
BaseType::IntegrationPointsArrayType IntegrationPointsArrayType
Definition: hexahedra_3d_20.h:123
BaseType::CoordinatesArrayType CoordinatesArrayType
Definition: hexahedra_3d_20.h:168
SizeType EdgesNumber() const override
This method gives you number of all edges of this geometry.
Definition: hexahedra_3d_20.h:525
Hexahedra3D20 & operator=(const Hexahedra3D20 &rOther)
Definition: hexahedra_3d_20.h:342
friend class Hexahedra3D20
Definition: hexahedra_3d_20.h:1571
Quadrilateral3D8< TPointType > FaceType
Definition: hexahedra_3d_20.h:67
Vector & ShapeFunctionsValues(Vector &rResult, const CoordinatesArrayType &rCoordinates) const override
Definition: hexahedra_3d_20.h:790
BaseType::Pointer Create(const IndexType NewGeometryId, const BaseType &rGeometry) const override
Creates a new geometry pointer.
Definition: hexahedra_3d_20.h:391
Hexahedra3D20(const std::string &GeometryName, const PointsArrayType &rThisPoints)
Constructor with Geometry Name.
Definition: hexahedra_3d_20.h:274
BaseType::IntegrationPointsContainerType IntegrationPointsContainerType
Definition: hexahedra_3d_20.h:130
Hexahedra3D20(const PointsArrayType &ThisPoints)
Definition: hexahedra_3d_20.h:258
GeometryData::KratosGeometryFamily GetGeometryFamily() const override
Definition: hexahedra_3d_20.h:318
Matrix MatrixType
Definition: hexahedra_3d_20.h:173
GeometryData::IntegrationMethod IntegrationMethod
Definition: hexahedra_3d_20.h:77
BaseType::GeometriesArrayType GeometriesArrayType
Definition: hexahedra_3d_20.h:83
SizeType FacesNumber() const override
Returns the number of faces of the current geometry.
Definition: hexahedra_3d_20.h:609
double DomainSize() const override
Definition: hexahedra_3d_20.h:470
TPointType PointType
Definition: hexahedra_3d_20.h:88
KRATOS_CLASS_POINTER_DEFINITION(Hexahedra3D20)
void PrintData(std::ostream &rOStream) const override
Definition: hexahedra_3d_20.h:879
Matrix & ShapeFunctionsLocalGradients(Matrix &result, const CoordinatesArrayType &rPoint) const override
Definition: hexahedra_3d_20.h:904
Hexahedra3D20(const PointType &Point1, const PointType &Point2, const PointType &Point3, const PointType &Point4, const PointType &Point5, const PointType &Point6, const PointType &Point7, const PointType &Point8, const PointType &Point9, const PointType &Point10, const PointType &Point11, const PointType &Point12, const PointType &Point13, const PointType &Point14, const PointType &Point15, const PointType &Point16, const PointType &Point17, const PointType &Point18, const PointType &Point19, const PointType &Point20)
Definition: hexahedra_3d_20.h:179
double Length() const override
Definition: hexahedra_3d_20.h:420
BaseType::NormalType NormalType
Definition: hexahedra_3d_20.h:163
Hexahedra3D20 & operator=(Hexahedra3D20< TOtherPointType > const &rOther)
Definition: hexahedra_3d_20.h:360
BaseType::ShapeFunctionsGradientsType ShapeFunctionsGradientsType
Definition: hexahedra_3d_20.h:158
BaseType::IndexType IndexType
Definition: hexahedra_3d_20.h:95
Line3D3< TPointType > EdgeType
Definition: hexahedra_3d_20.h:66
GeometriesArrayType GenerateFaces() const override
Returns all faces of the current geometry.
Definition: hexahedra_3d_20.h:622
void PrintInfo(std::ostream &rOStream) const override
Definition: hexahedra_3d_20.h:865
double Area() const override
Definition: hexahedra_3d_20.h:438
Hexahedra3D20(const IndexType GeometryId, const PointsArrayType &rThisPoints)
Constructor with Geometry Id.
Definition: hexahedra_3d_20.h:265
std::string Info() const override
Definition: hexahedra_3d_20.h:853
BaseType::JacobiansType JacobiansType
Definition: hexahedra_3d_20.h:151
double Volume() const override
This method calculate and return volume of this geometry.
Definition: hexahedra_3d_20.h:452
BaseType::Pointer Create(const IndexType NewGeometryId, PointsArrayType const &rThisPoints) const override
Creates a new geometry pointer.
Definition: hexahedra_3d_20.h:377
Hexahedra3D20(typename PointType::Pointer pPoint1, typename PointType::Pointer pPoint2, typename PointType::Pointer pPoint3, typename PointType::Pointer pPoint4, typename PointType::Pointer pPoint5, typename PointType::Pointer pPoint6, typename PointType::Pointer pPoint7, typename PointType::Pointer pPoint8, typename PointType::Pointer pPoint9, typename PointType::Pointer pPoint10, typename PointType::Pointer pPoint11, typename PointType::Pointer pPoint12, typename PointType::Pointer pPoint13, typename PointType::Pointer pPoint14, typename PointType::Pointer pPoint15, typename PointType::Pointer pPoint16, typename PointType::Pointer pPoint17, typename PointType::Pointer pPoint18, typename PointType::Pointer pPoint19, typename PointType::Pointer pPoint20)
Definition: hexahedra_3d_20.h:214
GeometryData::KratosGeometryType GetGeometryType() const override
Definition: hexahedra_3d_20.h:323
Geometry< TPointType > BaseType
Definition: hexahedra_3d_20.h:61
Hexahedra3D20(Hexahedra3D20< TOtherPointType > const &rOther)
Definition: hexahedra_3d_20.h:307
GeometriesArrayType GenerateEdges() const override
This method gives you all edges of this geometry.
Definition: hexahedra_3d_20.h:538
Hexahedra3D20(Hexahedra3D20 const &rOther)
Definition: hexahedra_3d_20.h:291
BaseType::ShapeFunctionsValuesContainerType ShapeFunctionsValuesContainerType
Definition: hexahedra_3d_20.h:137
BaseType::SizeType SizeType
Definition: hexahedra_3d_20.h:102
bool IsInside(const CoordinatesArrayType &rPoint, CoordinatesArrayType &rResult, const double Tolerance=std::numeric_limits< double >::epsilon()) const override
Returns whether given arbitrary point is inside the Geometry and the respective local point for the g...
Definition: hexahedra_3d_20.h:483
BaseType::PointsArrayType PointsArrayType
Definition: hexahedra_3d_20.h:108
BaseType::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: hexahedra_3d_20.h:144
~Hexahedra3D20() override
Definition: hexahedra_3d_20.h:315
BaseType::IntegrationPointType IntegrationPointType
Definition: hexahedra_3d_20.h:115
Short class definition.
Definition: integration_point.h:52
static double ComputeVolume3DGeometry(const Geometry< TPointType > &rGeometry)
This method calculates and returns the volume of the geometry from a 3D geometry.
Definition: integration_utilities.h:168
Definition: amatrix_interface.h:41
void resize(std::size_t NewSize1, std::size_t NewSize2, bool preserve=0)
Definition: amatrix_interface.h:224
An three node 3D line geometry with quadratic shape functions.
Definition: line_3d_3.h:66
This class defines the node.
Definition: node.h:65
PointerVector is a container like stl vector but using a vector to store pointers to its data.
Definition: pointer_vector.h:72
void push_back(const TPointerType &x)
Definition: pointer_vector.h:270
A eight node 3D quadrilateral geometry with quadratic shape functions.
Definition: quadrilateral_3d_8.h:49
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
Short class definition.
Definition: array_1d.h:61
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
#define KRATOS_ERROR
Definition: exception.h:161
#define KRATOS_ERROR_IF(conditional)
Definition: exception.h:162
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
KratosZeroMatrix< double > ZeroMatrix
Definition: amatrix_interface.h:559
Internals::Matrix< double, AMatrix::dynamic, AMatrix::dynamic > Matrix
Definition: amatrix_interface.h:470
const GeometryData Hexahedra3D20< TPointType >::msGeometryData & msGeometryDimension(), Hexahedra3D20< TPointType >::AllShapeFunctionsValues(), AllShapeFunctionsLocalGradients()
Definition: brep_curve.h:483
std::istream & operator>>(std::istream &rIStream, LinearMasterSlaveConstraint &rThis)
input stream function
const GeometryData BrepCurve< TContainerPointType, TContainerPointEmbeddedType >::msGeometryData & msGeometryDimension
Definition: brep_curve.h:483
std::ostream & operator<<(std::ostream &rOStream, const LinearMasterSlaveConstraint &rThis)
output stream function
Definition: linear_master_slave_constraint.h:432
def load(f)
Definition: ode_solve.py:307