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.
quadrilateral_2d_9.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
24 #include "geometries/line_2d_3.h"
27 
28 
29 namespace Kratos
30 {
48 template<class TPointType> class Quadrilateral2D9 : public Geometry<TPointType>
49 {
50 public:
59 
64 
69 
74 
80 
84  typedef TPointType PointType;
85 
90 
91 
98  typedef typename BaseType::IndexType IndexType;
99 
105  typedef typename BaseType::SizeType SizeType;
106 
112 
119 
127 
134 
141 
148 
155 
162 
170 
178 
179 
184 
189  Quadrilateral2D9( const PointType& Point1, const PointType& Point2,
190  const PointType& Point3, const PointType& Point4,
191  const PointType& Point5, const PointType& Point6,
192  const PointType& Point7, const PointType& Point8,
193  const PointType& Point9 )
194  : BaseType( PointsArrayType(), &msGeometryData )
195  {
196  this->Points().push_back( typename PointType::Pointer( new PointType( Point1 ) ) );
197  this->Points().push_back( typename PointType::Pointer( new PointType( Point2 ) ) );
198  this->Points().push_back( typename PointType::Pointer( new PointType( Point3 ) ) );
199  this->Points().push_back( typename PointType::Pointer( new PointType( Point4 ) ) );
200  this->Points().push_back( typename PointType::Pointer( new PointType( Point5 ) ) );
201  this->Points().push_back( typename PointType::Pointer( new PointType( Point6 ) ) );
202  this->Points().push_back( typename PointType::Pointer( new PointType( Point7 ) ) );
203  this->Points().push_back( typename PointType::Pointer( new PointType( Point8 ) ) );
204  this->Points().push_back( typename PointType::Pointer( new PointType( Point9 ) ) );
205  }
206 
207  Quadrilateral2D9( typename PointType::Pointer pPoint1,
208  typename PointType::Pointer pPoint2,
209  typename PointType::Pointer pPoint3,
210  typename PointType::Pointer pPoint4,
211  typename PointType::Pointer pPoint5,
212  typename PointType::Pointer pPoint6,
213  typename PointType::Pointer pPoint7,
214  typename PointType::Pointer pPoint8,
215  typename PointType::Pointer pPoint9 )
216  : BaseType( PointsArrayType(), &msGeometryData )
217  {
218  this->Points().push_back( pPoint1 );
219  this->Points().push_back( pPoint2 );
220  this->Points().push_back( pPoint3 );
221  this->Points().push_back( pPoint4 );
222  this->Points().push_back( pPoint5 );
223  this->Points().push_back( pPoint6 );
224  this->Points().push_back( pPoint7 );
225  this->Points().push_back( pPoint8 );
226  this->Points().push_back( pPoint9 );
227  }
228 
229  Quadrilateral2D9( const PointsArrayType& ThisPoints )
230  : BaseType( ThisPoints, &msGeometryData )
231  {
232  if ( this->PointsNumber() != 9 )
233  KRATOS_ERROR << "Invalid points number. Expected 9, given " << this->PointsNumber() << std::endl;
234  }
235 
238  IndexType GeometryId,
239  const PointsArrayType& rThisPoints
240  ) : BaseType(GeometryId, rThisPoints, &msGeometryData)
241  {
242  KRATOS_ERROR_IF( this->PointsNumber() != 9 ) << "Invalid points number. Expected 9, given " << this->PointsNumber() << std::endl;
243  }
244 
247  const std::string& rGeometryName,
248  const PointsArrayType& rThisPoints
249  ) : BaseType(rGeometryName, rThisPoints, &msGeometryData)
250  {
251  KRATOS_ERROR_IF(this->PointsNumber() != 9) << "Invalid points number. Expected 9, given " << this->PointsNumber() << std::endl;
252  }
253 
264  : BaseType( rOther )
265  {
266  }
267 
280  template<class TOtherPointType> Quadrilateral2D9(
281  Quadrilateral2D9<TOtherPointType> const& rOther )
282  : BaseType( rOther )
283  {
284  }
285 
289  ~Quadrilateral2D9() override {}
290 
292  {
294  }
295 
297  {
299  }
300 
317  {
318  BaseType::operator=( rOther );
319 
320  return *this;
321  }
322 
334  template<class TOtherPointType>
336  {
337  BaseType::operator=( rOther );
338 
339  return *this;
340  }
341 
352  typename BaseType::Pointer Create(
353  const IndexType rNewGeometryId,
354  PointsArrayType const& rThisPoints
355  ) const override
356  {
357  return typename BaseType::Pointer( new Quadrilateral2D9( rNewGeometryId, rThisPoints ) );
358  }
359 
366  typename BaseType::Pointer Create(
367  const IndexType NewGeometryId,
368  const BaseType& rGeometry
369  ) const override
370  {
371  auto p_geometry = typename BaseType::Pointer( new Quadrilateral2D9( NewGeometryId, rGeometry.Points() ) );
372  p_geometry->SetData(rGeometry.GetData());
373  return p_geometry;
374  }
375 
381  SizeType PointsNumberInDirection(IndexType LocalDirectionIndex) const override
382  {
383  if ((LocalDirectionIndex == 0) || (LocalDirectionIndex == 1)) {
384  return 3;
385  }
386  KRATOS_ERROR << "Possible direction index reaches from 0-1. Given direction index: "
387  << LocalDirectionIndex << std::endl;
388  }
389 
406  double Length() const override
407  {
408  return std::sqrt( std::abs( this->DeterminantOfJacobian( PointType() ) ) );
409  }
410 
423  double Area() const override
424  {
426  }
427 
435  double Volume() const override
436  {
437  KRATOS_WARNING("Quadrilateral2D9") << "Method not well defined. Replace with DomainSize() instead. This method preserves current behaviour but will be changed in June 2023 (returning zero instead)" << std::endl;
438  return Area();
439  // TODO: Replace in June 2023
440  // KRATOS_ERROR << "Quadrilateral2D9:: Method not well defined. Replace with DomainSize() instead." << std::endl;
441  // return 0.0;
442  }
443 
455  double DomainSize() const override
456  {
457  return Area();
458  }
459 
468  bool IsInside(
469  const CoordinatesArrayType& rPoint,
470  CoordinatesArrayType& rResult,
471  const double Tolerance = std::numeric_limits<double>::epsilon()
472  ) const override
473  {
474  this->PointLocalCoordinates( rResult, rPoint );
475 
476  if ( (rResult[0] >= (-1.0-Tolerance)) && (rResult[0] <= (1.0+Tolerance)) ) {
477  if ( (rResult[1] >= (-1.0-Tolerance)) && (rResult[1] <= (1.0+Tolerance)) ) {
478  return true;
479  }
480  }
481 
482  return false;
483  }
484 
496  SizeType EdgesNumber() const override
497  {
498  return 4;
499  }
500 
510  {
512  edges.push_back( Kratos::make_shared<EdgeType>( this->pGetPoint( 0 ), this->pGetPoint( 1 ), this->pGetPoint( 4 ) ) );
513  edges.push_back( Kratos::make_shared<EdgeType>( this->pGetPoint( 1 ), this->pGetPoint( 2 ), this->pGetPoint( 5 ) ) );
514  edges.push_back( Kratos::make_shared<EdgeType>( this->pGetPoint( 2 ), this->pGetPoint( 3 ), this->pGetPoint( 6 ) ) );
515  edges.push_back( Kratos::make_shared<EdgeType>( this->pGetPoint( 3 ), this->pGetPoint( 0 ), this->pGetPoint( 7 ) ) );
516  return edges;
517  }
518 
535  double ShapeFunctionValue( IndexType ShapeFunctionIndex,
536  const CoordinatesArrayType& rPoint ) const override
537  {
538  double fx1 = 0.5 * ( rPoint[0] - 1 ) * rPoint[0];
539  double fx2 = 0.5 * ( rPoint[0] + 1 ) * rPoint[0];
540  double fx3 = 1 - rPoint[0] * rPoint[0];
541  double fy1 = 0.5 * ( rPoint[1] - 1 ) * rPoint[1];
542  double fy2 = 0.5 * ( rPoint[1] + 1 ) * rPoint[1];
543  double fy3 = 1 - rPoint[1] * rPoint[1];
544 
545  switch ( ShapeFunctionIndex )
546  {
547  case 0:
548  return( fx1*fy1 );
549  case 1:
550  return( fx2*fy1 );
551  case 2:
552  return( fx2*fy2 );
553  case 3:
554  return( fx1*fy2 );
555  case 4:
556  return( fx3*fy1 );
557  case 5:
558  return( fx2*fy3 );
559  case 6:
560  return( fx3*fy2 );
561  case 7:
562  return( fx1*fy3 );
563  case 8:
564  return( fx3*fy3 );
565  default:
566  KRATOS_ERROR << "Wrong index of shape function!" << *this << std::endl;
567  }
568 
569  return 0;
570  }
571 
583  Vector& ShapeFunctionsValues (Vector &rResult, const CoordinatesArrayType& rCoordinates) const override
584  {
585  if(rResult.size() != 9) rResult.resize(9,false);
586 
587  double fx1 = 0.5 * ( rCoordinates[0] - 1 ) * rCoordinates[0];
588  double fx2 = 0.5 * ( rCoordinates[0] + 1 ) * rCoordinates[0];
589  double fx3 = 1 - rCoordinates[0] * rCoordinates[0];
590  double fy1 = 0.5 * ( rCoordinates[1] - 1 ) * rCoordinates[1];
591  double fy2 = 0.5 * ( rCoordinates[1] + 1 ) * rCoordinates[1];
592  double fy3 = 1 - rCoordinates[1] * rCoordinates[1];
593 
594  rResult[0] = fx1*fy1 ;
595  rResult[1] = fx2*fy1 ;
596  rResult[2] = fx2*fy2 ;
597  rResult[3] = fx1*fy2 ;
598  rResult[4] = fx3*fy1 ;
599  rResult[5] = fx2*fy3 ;
600  rResult[6] = fx3*fy2 ;
601  rResult[7] = fx1*fy3 ;
602  rResult[8] = fx3*fy3 ;
603 
604  return rResult;
605  }
606 
607 
619  std::string Info() const override
620  {
621  return "2 dimensional quadrilateral with nine nodes in 2D space";
622  }
623 
630  void PrintInfo( std::ostream& rOStream ) const override
631  {
632  rOStream << "2 dimensional quadrilateral with nine nodes in 2D space";
633  }
634 
644  void PrintData( std::ostream& rOStream ) const override
645  {
646  // Base Geometry class PrintData call
647  BaseType::PrintData( rOStream );
648  std::cout << std::endl;
649 
650  // If the geometry has valid points, calculate and output its data
651  if (this->AllPointsAreValid()) {
652  Matrix jacobian;
653  this->Jacobian( jacobian, PointType() );
654  rOStream << " Jacobian in the origin\t : " << jacobian;
655  }
656  }
657 
663  IntegrationMethod ThisMethod )
664  {
665  ShapeFunctionsGradientsType localGradients
666  = CalculateShapeFunctionsIntegrationPointsLocalGradients( ThisMethod );
667  const int integration_points_number
668  = msGeometryData.IntegrationPointsNumber( ThisMethod );
669  ShapeFunctionsGradientsType Result( integration_points_number );
670 
671  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
672  {
673  Result[pnt] = localGradients[pnt];
674  }
675 
676  return Result;
677  }
678 
684  {
685  IntegrationMethod ThisMethod = msGeometryData.DefaultIntegrationMethod();
686  ShapeFunctionsGradientsType localGradients
687  = CalculateShapeFunctionsIntegrationPointsLocalGradients( ThisMethod );
688  const int integration_points_number
689  = msGeometryData.IntegrationPointsNumber( ThisMethod );
690  ShapeFunctionsGradientsType Result( integration_points_number );
691 
692  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
693  {
694  Result[pnt] = localGradients[pnt];
695  }
696 
697  return Result;
698  }
699 
709  const CoordinatesArrayType& rPoint ) const override
710  {
711  double fx1 = 0.5 * ( rPoint[0] - 1 ) * rPoint[0];
712  double fx2 = 0.5 * ( rPoint[0] + 1 ) * rPoint[0];
713  double fx3 = 1 - rPoint[0] * rPoint[0];
714  double fy1 = 0.5 * ( rPoint[1] - 1 ) * rPoint[1];
715  double fy2 = 0.5 * ( rPoint[1] + 1 ) * rPoint[1];
716  double fy3 = 1 - rPoint[1] * rPoint[1];
717 
718  double gx1 = 0.5 * ( 2 * rPoint[0] - 1 );
719  double gx2 = 0.5 * ( 2 * rPoint[0] + 1 );
720  double gx3 = -2.0 * rPoint[0];
721  double gy1 = 0.5 * ( 2 * rPoint[1] - 1 );
722  double gy2 = 0.5 * ( 2 * rPoint[1] + 1 );
723  double gy3 = -2.0 * rPoint[1];
724 
725  rResult.resize( 9, 2, false );
726  noalias( rResult ) = ZeroMatrix( 9, 2 );
727  rResult( 0, 0 ) = gx1 * fy1;
728  rResult( 0, 1 ) = fx1 * gy1;
729  rResult( 1, 0 ) = gx2 * fy1;
730  rResult( 1, 1 ) = fx2 * gy1;
731  rResult( 2, 0 ) = gx2 * fy2;
732  rResult( 2, 1 ) = fx2 * gy2;
733  rResult( 3, 0 ) = gx1 * fy2;
734  rResult( 3, 1 ) = fx1 * gy2;
735  rResult( 4, 0 ) = gx3 * fy1;
736  rResult( 4, 1 ) = fx3 * gy1;
737  rResult( 5, 0 ) = gx2 * fy3;
738  rResult( 5, 1 ) = fx2 * gy3;
739  rResult( 6, 0 ) = gx3 * fy2;
740  rResult( 6, 1 ) = fx3 * gy2;
741  rResult( 7, 0 ) = gx1 * fy3;
742  rResult( 7, 1 ) = fx1 * gy3;
743  rResult( 8, 0 ) = gx3 * fy3;
744  rResult( 8, 1 ) = fx3 * gy3;
745 
746  return rResult;
747  }
748 
754  Matrix& PointsLocalCoordinates( Matrix& rResult ) const override
755  {
756  rResult.resize( 9, 2, false );
757  noalias( rResult ) = ZeroMatrix( 9, 2 );
758  rResult( 0, 0 ) = -1.0;
759  rResult( 0, 1 ) = -1.0;
760  rResult( 1, 0 ) = 1.0;
761  rResult( 1, 1 ) = -1.0;
762  rResult( 2, 0 ) = 1.0;
763  rResult( 2, 1 ) = 1.0;
764  rResult( 3, 0 ) = -1.0;
765  rResult( 3, 1 ) = 1.0;
766  rResult( 4, 0 ) = 0.0;
767  rResult( 4, 1 ) = -1.0;
768  rResult( 5, 0 ) = 1.0;
769  rResult( 5, 1 ) = 0.0;
770  rResult( 6, 0 ) = 0.0;
771  rResult( 6, 1 ) = 1.0;
772  rResult( 7, 0 ) = -1.0;
773  rResult( 7, 1 ) = 0.0;
774  rResult( 8, 0 ) = 0.0;
775  rResult( 8, 1 ) = 0.0;
776  return rResult;
777  }
778 
787  virtual Matrix& ShapeFunctionsGradients( Matrix& rResult, PointType& rPoint )
788  {
789  double fx1 = 0.5 * ( rPoint.X() - 1 ) * rPoint.X();
790  double fx2 = 0.5 * ( rPoint.X() + 1 ) * rPoint.X();
791  double fx3 = 1 - rPoint.X() * rPoint.X();
792  double fy1 = 0.5 * ( rPoint.Y() - 1 ) * rPoint.Y();
793  double fy2 = 0.5 * ( rPoint.Y() + 1 ) * rPoint.Y();
794  double fy3 = 1 - rPoint.Y() * rPoint.Y();
795 
796  double gx1 = 0.5 * ( 2 * rPoint.X() - 1 );
797  double gx2 = 0.5 * ( 2 * rPoint.X() + 1 );
798  double gx3 = -2.0 * rPoint.X();
799  double gy1 = 0.5 * ( 2 * rPoint.Y() - 1 );
800  double gy2 = 0.5 * ( 2 * rPoint.Y() + 1 );
801  double gy3 = -2.0 * rPoint.Y();
802 
803  rResult.resize( 9, 2, false );
804  noalias( rResult ) = ZeroMatrix( 9, 2 );
805  rResult( 0, 0 ) = gx1 * fy1;
806  rResult( 0, 1 ) = fx1 * gy1;
807  rResult( 1, 0 ) = gx2 * fy1;
808  rResult( 1, 1 ) = fx2 * gy1;
809  rResult( 2, 0 ) = gx2 * fy2;
810  rResult( 2, 1 ) = fx2 * gy2;
811  rResult( 3, 0 ) = gx1 * fy2;
812  rResult( 3, 1 ) = fx1 * gy2;
813  rResult( 4, 0 ) = gx3 * fy1;
814  rResult( 4, 1 ) = fx3 * gy1;
815  rResult( 5, 0 ) = gx2 * fy3;
816  rResult( 5, 1 ) = fx2 * gy3;
817  rResult( 6, 0 ) = gx3 * fy2;
818  rResult( 6, 1 ) = fx3 * gy2;
819  rResult( 7, 0 ) = gx1 * fy3;
820  rResult( 7, 1 ) = fx1 * gy3;
821  rResult( 8, 0 ) = gx3 * fy3;
822  rResult( 8, 1 ) = fx3 * gy3;
823 
824  return rResult;
825  }
826 
834  {
835  if ( rResult.size() != this->PointsNumber() )
836  {
837  // KLUDGE: While there is a bug in ublas vector resize, I have to put this beside resizing!!
839  rResult.swap( temp );
840  }
841 
842  for ( unsigned int i = 0; i < this->PointsNumber(); i++ )
843  {
844  rResult[i].resize( 2, 2, false );
845  noalias( rResult[i] ) = ZeroMatrix( 2, 2 );
846  }
847 
848  double fx1 = 0.5 * ( rPoint[0] - 1 ) * rPoint[0];
849  double fx2 = 0.5 * ( rPoint[0] + 1 ) * rPoint[0];
850  double fx3 = 1 - rPoint[0] * rPoint[0];
851  double fy1 = 0.5 * ( rPoint[1] - 1 ) * rPoint[1];
852  double fy2 = 0.5 * ( rPoint[1] + 1 ) * rPoint[1];
853  double fy3 = 1 - rPoint[1] * rPoint[1];
854 
855  double gx1 = 0.5 * ( 2 * rPoint[0] - 1 );
856  double gx2 = 0.5 * ( 2 * rPoint[0] + 1 );
857  double gx3 = -2.0 * rPoint[0];
858  double gy1 = 0.5 * ( 2 * rPoint[1] - 1 );
859  double gy2 = 0.5 * ( 2 * rPoint[1] + 1 );
860  double gy3 = -2.0 * rPoint[1];
861 
862  double hx1 = 1.0;
863  double hx2 = 1.0;
864  double hx3 = -2.0;
865  double hy1 = 1.0;
866  double hy2 = 1.0;
867  double hy3 = -2.0;
868 
869  rResult[0]( 0, 0 ) = hx1 * fy1;
870  rResult[0]( 0, 1 ) = gx1 * gy1;
871  rResult[0]( 1, 0 ) = gx1 * gy1;
872  rResult[0]( 1, 1 ) = fx1 * hy1;
873 
874  rResult[1]( 0, 0 ) = hx2 * fy1;
875  rResult[1]( 0, 1 ) = gx2 * gy1;
876  rResult[1]( 1, 0 ) = gx2 * gy1;
877  rResult[1]( 1, 1 ) = fx2 * hy1;
878 
879  rResult[2]( 0, 0 ) = hx2 * fy2;
880  rResult[2]( 0, 1 ) = gx2 * gy2;
881  rResult[2]( 1, 0 ) = gx2 * gy2;
882  rResult[2]( 1, 1 ) = fx2 * hy2;
883 
884  rResult[3]( 0, 0 ) = hx1 * fy2;
885  rResult[3]( 0, 1 ) = gx1 * gy2;
886  rResult[3]( 1, 0 ) = gx1 * gy2;
887  rResult[3]( 1, 1 ) = fx1 * hy2;
888 
889  rResult[4]( 0, 0 ) = hx3 * fy1;
890  rResult[4]( 0, 1 ) = gx3 * gy1;
891  rResult[4]( 1, 0 ) = gx3 * gy1;
892  rResult[4]( 1, 1 ) = fx3 * hy1;
893 
894  rResult[5]( 0, 0 ) = hx2 * fy3;
895  rResult[5]( 0, 1 ) = gx2 * gy3;
896  rResult[5]( 1, 0 ) = gx2 * gy3;
897  rResult[5]( 1, 1 ) = fx2 * hy3;
898 
899  rResult[6]( 0, 0 ) = hx3 * fy2;
900  rResult[6]( 0, 1 ) = gx3 * gy2;
901  rResult[6]( 1, 0 ) = gx3 * gy2;
902  rResult[6]( 1, 1 ) = fx3 * hy2;
903 
904  rResult[7]( 0, 0 ) = hx1 * fy3;
905  rResult[7]( 0, 1 ) = gx1 * gy3;
906  rResult[7]( 1, 0 ) = gx1 * gy3;
907  rResult[7]( 1, 1 ) = fx1 * hy3;
908 
909  rResult[8]( 0, 0 ) = hx3 * fy3;
910  rResult[8]( 0, 1 ) = gx3 * gy3;
911  rResult[8]( 1, 0 ) = gx3 * gy3;
912  rResult[8]( 1, 1 ) = fx3 * hy3;
913 
914  return rResult;
915  }
916 
917 
925  {
926 
927  if ( rResult.size() != this->PointsNumber() )
928  {
929  // KLUDGE: While there is a bug in
930  // ublas vector resize, I have to put this beside resizing!!
931 // ShapeFunctionsGradientsType
933  rResult.swap( temp );
934  }
935 
936  for ( IndexType i = 0; i < rResult.size(); i++ )
937  {
939  rResult[i].swap( temp );
940  }
941 
942  for ( unsigned int i = 0; i < this->PointsNumber(); i++ )
943  {
944  for ( unsigned int j = 0; j < 2; j++ )
945  {
946  rResult[i][j].resize( 2, 2, false );
947  noalias( rResult[i][j] ) = ZeroMatrix( 2, 2 );
948  }
949  }
950 
951 // double fx1 = 0.5 * ( rPoint[0] - 1 ) * rPoint[0];
952 
953 // double fx2 = 0.5 * ( rPoint[0] + 1 ) * rPoint[0];
954 // double fx3 = 1 - rPoint[0] * rPoint[0];
955 // double fy1 = 0.5 * ( rPoint[1] - 1 ) * rPoint[1];
956 // double fy2 = 0.5 * ( rPoint[1] + 1 ) * rPoint[1];
957 // double fy3 = 1 - rPoint[1] * rPoint[1];
958 
959  double gx1 = 0.5 * ( 2 * rPoint[0] - 1 );
960  double gx2 = 0.5 * ( 2 * rPoint[0] + 1 );
961  double gx3 = -2.0 * rPoint[0];
962  double gy1 = 0.5 * ( 2 * rPoint[1] - 1 );
963  double gy2 = 0.5 * ( 2 * rPoint[1] + 1 );
964  double gy3 = -2.0 * rPoint[1];
965 
966  double hx1 = 1.0;
967  double hx2 = 1.0;
968  double hx3 = -2.0;
969  double hy1 = 1.0;
970  double hy2 = 1.0;
971  double hy3 = -2.0;
972 
973  rResult[0][0]( 0, 0 ) = 0.0;
974  rResult[0][0]( 0, 1 ) = hx1 * gy1;
975  rResult[0][0]( 1, 0 ) = hx1 * gy1;
976  rResult[0][0]( 1, 1 ) = gx1 * hy1;
977  rResult[0][1]( 0, 0 ) = hx1 * gy1;
978  rResult[0][1]( 0, 1 ) = gx1 * hy1;
979  rResult[0][1]( 1, 0 ) = gx1 * hy1;
980  rResult[0][1]( 1, 1 ) = 0.0;
981 
982  rResult[1][0]( 0, 0 ) = 0.0;
983  rResult[1][0]( 0, 1 ) = hx2 * gy1;
984  rResult[1][0]( 1, 0 ) = hx2 * gy1;
985  rResult[1][0]( 1, 1 ) = gx2 * hy1;
986  rResult[1][1]( 0, 0 ) = hx2 * gy1;
987  rResult[1][1]( 0, 1 ) = gx2 * hy1;
988  rResult[1][1]( 1, 0 ) = gx2 * hy1;
989  rResult[1][1]( 1, 1 ) = 0.0;
990 
991  rResult[2][0]( 0, 0 ) = 0.0;
992  rResult[2][0]( 0, 1 ) = hx2 * gy2;
993  rResult[2][0]( 1, 0 ) = hx2 * gy2;
994  rResult[2][0]( 1, 1 ) = gx2 * hy2;
995  rResult[2][1]( 0, 0 ) = hx2 * gy2;
996  rResult[2][1]( 0, 1 ) = gx2 * hy2;
997  rResult[2][1]( 1, 0 ) = gx2 * hy2;
998  rResult[2][1]( 1, 1 ) = 0.0;
999 
1000  rResult[3][0]( 0, 0 ) = 0.0;
1001  rResult[3][0]( 0, 1 ) = hx1 * gy2;
1002  rResult[3][0]( 1, 0 ) = hx1 * gy2;
1003  rResult[3][0]( 1, 1 ) = gx1 * hy2;
1004  rResult[3][1]( 0, 0 ) = hx1 * gy2;
1005  rResult[3][1]( 0, 1 ) = gx1 * hy2;
1006  rResult[3][1]( 1, 0 ) = gx1 * hy2;
1007  rResult[3][1]( 1, 1 ) = 0.0;
1008 
1009  rResult[4][0]( 0, 0 ) = 0.0;
1010  rResult[4][0]( 0, 1 ) = hx3 * gy1;
1011  rResult[4][0]( 1, 0 ) = hx3 * gy1;
1012  rResult[4][0]( 1, 1 ) = gx3 * hy1;
1013  rResult[4][1]( 0, 0 ) = hx3 * gy1;
1014  rResult[4][1]( 0, 1 ) = gx3 * hy1;
1015  rResult[4][1]( 1, 0 ) = gx3 * hy1;
1016  rResult[4][1]( 1, 1 ) = 0.0;
1017 
1018  rResult[5][0]( 0, 0 ) = 0.0;
1019  rResult[5][0]( 0, 1 ) = hx2 * gy3;
1020  rResult[5][0]( 1, 0 ) = hx2 * gy3;
1021  rResult[5][0]( 1, 1 ) = gx2 * hy3;
1022  rResult[5][1]( 0, 0 ) = hx2 * gy3;
1023  rResult[5][1]( 0, 1 ) = gx2 * hy3;
1024  rResult[5][1]( 1, 0 ) = gx2 * hy3;
1025  rResult[5][1]( 1, 1 ) = 0.0;
1026 
1027  rResult[6][0]( 0, 0 ) = 0.0;
1028  rResult[6][0]( 0, 1 ) = hx3 * gy2;
1029  rResult[6][0]( 1, 0 ) = hx3 * gy2;
1030  rResult[6][0]( 1, 1 ) = gx3 * hy2;
1031  rResult[6][1]( 0, 0 ) = hx3 * gy2;
1032  rResult[6][1]( 0, 1 ) = gx3 * hy2;
1033  rResult[6][1]( 1, 0 ) = gx3 * hy2;
1034  rResult[6][1]( 1, 1 ) = 0.0;
1035 
1036  rResult[7][0]( 0, 0 ) = 0.0;
1037  rResult[7][0]( 0, 1 ) = hx1 * gy3;
1038  rResult[7][0]( 1, 0 ) = hx1 * gy3;
1039  rResult[7][0]( 1, 1 ) = gx1 * hy3;
1040  rResult[7][1]( 0, 0 ) = hx1 * gy3;
1041  rResult[7][1]( 0, 1 ) = gx1 * hy3;
1042  rResult[7][1]( 1, 0 ) = gx1 * hy3;
1043  rResult[7][1]( 1, 1 ) = 0.0;
1044 
1045  rResult[8][0]( 0, 0 ) = 0.0;
1046  rResult[8][0]( 0, 1 ) = hx3 * gy3;
1047  rResult[8][0]( 1, 0 ) = hx3 * gy3;
1048  rResult[8][0]( 1, 1 ) = gx3 * hy3;
1049  rResult[8][1]( 0, 0 ) = hx3 * gy3;
1050  rResult[8][1]( 0, 1 ) = gx3 * hy3;
1051  rResult[8][1]( 1, 0 ) = gx3 * hy3;
1052  rResult[8][1]( 1, 1 ) = 0.0;
1053 
1054 
1055  return rResult;
1056  }
1057 
1058 protected:
1059 
1060 
1066 private:
1067 
1072  static const GeometryData msGeometryData;
1073 
1074  static const GeometryDimension msGeometryDimension;
1075 
1079 
1080  friend class Serializer;
1081 
1082  void save( Serializer& rSerializer ) const override
1083  {
1085  }
1086 
1087  void load( Serializer& rSerializer ) override
1088  {
1090  }
1091 
1092  Quadrilateral2D9(): BaseType( PointsArrayType(), &msGeometryData ) {}
1093 
1112  static Matrix CalculateShapeFunctionsIntegrationPointsValues(
1113  typename BaseType::IntegrationMethod ThisMethod )
1114  {
1115  IntegrationPointsContainerType all_integration_points = AllIntegrationPoints();
1116  IntegrationPointsArrayType integration_points = all_integration_points[static_cast<int>(ThisMethod)];
1117  //number of integration points
1118  const int integration_points_number = integration_points.size();
1119  //number of nodes in current geometry
1120  const int points_number = 9;
1121  //setting up return matrix
1122  Matrix shape_function_values( integration_points_number, points_number );
1123  //loop over all integration points
1124 
1125  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
1126  {
1127  double fx1 = 0.5 * ( integration_points[pnt].X() - 1 ) * integration_points[pnt].X();
1128  double fx2 = 0.5 * ( integration_points[pnt].X() + 1 ) * integration_points[pnt].X();
1129  double fx3 = 1 - integration_points[pnt].X() * integration_points[pnt].X();
1130  double fy1 = 0.5 * ( integration_points[pnt].Y() - 1 ) * integration_points[pnt].Y();
1131  double fy2 = 0.5 * ( integration_points[pnt].Y() + 1 ) * integration_points[pnt].Y();
1132  double fy3 = 1 - integration_points[pnt].Y() * integration_points[pnt].Y();
1133 
1134  shape_function_values( pnt, 0 ) = ( fx1 * fy1 );
1135  shape_function_values( pnt, 1 ) = ( fx2 * fy1 );
1136  shape_function_values( pnt, 2 ) = ( fx2 * fy2 );
1137  shape_function_values( pnt, 3 ) = ( fx1 * fy2 );
1138  shape_function_values( pnt, 4 ) = ( fx3 * fy1 );
1139  shape_function_values( pnt, 5 ) = ( fx2 * fy3 );
1140  shape_function_values( pnt, 6 ) = ( fx3 * fy2 );
1141  shape_function_values( pnt, 7 ) = ( fx1 * fy3 );
1142  shape_function_values( pnt, 8 ) = ( fx3 * fy3 );
1143  }
1144 
1145  return shape_function_values;
1146  }
1147 
1160  static ShapeFunctionsGradientsType CalculateShapeFunctionsIntegrationPointsLocalGradients(
1161  typename BaseType::IntegrationMethod ThisMethod )
1162  {
1163  IntegrationPointsContainerType all_integration_points = AllIntegrationPoints();
1164  IntegrationPointsArrayType integration_points = all_integration_points[static_cast<int>(ThisMethod)];
1165  //number of integration points
1166  const int integration_points_number = integration_points.size();
1167  ShapeFunctionsGradientsType d_shape_f_values( integration_points_number );
1168  //initialising container
1169  //std::fill(d_shape_f_values.begin(), d_shape_f_values.end(), Matrix(4,2));
1170  //loop over all integration points
1171 
1172  for ( int pnt = 0; pnt < integration_points_number; pnt++ )
1173  {
1174  double fx1 = 0.5 * ( integration_points[pnt].X() - 1 ) * integration_points[pnt].X();
1175  double fx2 = 0.5 * ( integration_points[pnt].X() + 1 ) * integration_points[pnt].X();
1176  double fx3 = 1 - integration_points[pnt].X() * integration_points[pnt].X();
1177  double fy1 = 0.5 * ( integration_points[pnt].Y() - 1 ) * integration_points[pnt].Y();
1178  double fy2 = 0.5 * ( integration_points[pnt].Y() + 1 ) * integration_points[pnt].Y();
1179  double fy3 = 1 - integration_points[pnt].Y() * integration_points[pnt].Y();
1180 
1181  double gx1 = 0.5 * ( 2 * integration_points[pnt].X() - 1 );
1182  double gx2 = 0.5 * ( 2 * integration_points[pnt].X() + 1 );
1183  double gx3 = -2.0 * integration_points[pnt].X();
1184  double gy1 = 0.5 * ( 2 * integration_points[pnt].Y() - 1 );
1185  double gy2 = 0.5 * ( 2 * integration_points[pnt].Y() + 1 );
1186  double gy3 = -2.0 * integration_points[pnt].Y();
1187 
1188  Matrix result( 9, 2 );
1189  result( 0, 0 ) = gx1 * fy1;
1190  result( 0, 1 ) = fx1 * gy1;
1191  result( 1, 0 ) = gx2 * fy1;
1192  result( 1, 1 ) = fx2 * gy1;
1193  result( 2, 0 ) = gx2 * fy2;
1194  result( 2, 1 ) = fx2 * gy2;
1195  result( 3, 0 ) = gx1 * fy2;
1196  result( 3, 1 ) = fx1 * gy2;
1197  result( 4, 0 ) = gx3 * fy1;
1198  result( 4, 1 ) = fx3 * gy1;
1199  result( 5, 0 ) = gx2 * fy3;
1200  result( 5, 1 ) = fx2 * gy3;
1201  result( 6, 0 ) = gx3 * fy2;
1202  result( 6, 1 ) = fx3 * gy2;
1203  result( 7, 0 ) = gx1 * fy3;
1204  result( 7, 1 ) = fx1 * gy3;
1205  result( 8, 0 ) = gx3 * fy3;
1206  result( 8, 1 ) = fx3 * gy3;
1207 
1208  d_shape_f_values[pnt] = result;
1209  }
1210 
1211  return d_shape_f_values;
1212  }
1213 
1217  static const IntegrationPointsContainerType AllIntegrationPoints()
1218  {
1219  IntegrationPointsContainerType integration_points =
1220  {
1221  {
1222  Quadrature < QuadrilateralGaussLegendreIntegrationPoints1,
1223  2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1224  Quadrature < QuadrilateralGaussLegendreIntegrationPoints2,
1225  2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1226  Quadrature < QuadrilateralGaussLegendreIntegrationPoints3,
1227  2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1228  Quadrature < QuadrilateralGaussLegendreIntegrationPoints4,
1229  2, IntegrationPoint<3> >::GenerateIntegrationPoints(),
1230  Quadrature < QuadrilateralGaussLegendreIntegrationPoints5,
1231  2, IntegrationPoint<3> >::GenerateIntegrationPoints()
1232  }
1233  };
1234  return integration_points;
1235  }
1236 
1240  static const ShapeFunctionsValuesContainerType AllShapeFunctionsValues()
1241  {
1242  ShapeFunctionsValuesContainerType shape_functions_values =
1243  {
1244  {
1245  Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1247  Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1249  Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1251  Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1253  Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsValues(
1255  }
1256  };
1257  return shape_functions_values;
1258  }
1259 
1263  static const ShapeFunctionsLocalGradientsContainerType AllShapeFunctionsLocalGradients()
1264  {
1265  ShapeFunctionsLocalGradientsContainerType shape_functions_local_gradients =
1266  {
1267  {
1268  Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients
1270  Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients
1272  Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients
1274  Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients
1276  Quadrilateral2D9<TPointType>::CalculateShapeFunctionsIntegrationPointsLocalGradients
1278  }
1279  };
1280  return shape_functions_local_gradients;
1281  }
1282 
1287  template<class TOtherPointType> friend class Quadrilateral2D9;
1288 
1293 }; // Class Quadrilateral2D9
1294 
1300 template< class TPointType > inline std::istream& operator >> (
1301  std::istream& rIStream,
1303 
1305 template< class TPointType > inline std::ostream& operator << (
1306  std::ostream& rOStream,
1307  const Quadrilateral2D9<TPointType>& rThis )
1308 {
1309  rThis.PrintInfo( rOStream );
1310  rOStream << std::endl;
1311  rThis.PrintData( rOStream );
1312  return rOStream;
1313 }
1314 
1315 template<class TPointType>
1316 const GeometryData Quadrilateral2D9<TPointType>::msGeometryData(
1319  Quadrilateral2D9<TPointType>::AllIntegrationPoints(),
1320  Quadrilateral2D9<TPointType>::AllShapeFunctionsValues(),
1321  AllShapeFunctionsLocalGradients()
1322 );
1323 
1324 template<class TPointType>
1325 const GeometryDimension Quadrilateral2D9<TPointType>::msGeometryDimension(2, 2);
1326 
1327 } // namespace Kratos.
Definition: geometry_data.h:60
KratosGeometryType
Definition: geometry_data.h:110
SizeType IntegrationPointsNumber() const
Definition: geometry_data.h:430
IntegrationMethod
Definition: geometry_data.h:76
KratosGeometryFamily
Definition: geometry_data.h:91
IntegrationMethod DefaultIntegrationMethod() const
Definition: geometry_data.h:425
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
Short class definition.
Definition: integration_point.h:52
static double ComputeArea2DGeometry(const Geometry< TPointType > &rGeometry)
This method calculates and returns the volume of the geometry from a 3D geometry.
Definition: integration_utilities.h:107
Definition: amatrix_interface.h:41
void swap(Matrix &Other)
Definition: amatrix_interface.h:289
void resize(std::size_t NewSize1, std::size_t NewSize2, bool preserve=0)
Definition: amatrix_interface.h:224
An three node 2D line geometry with quadratic shape functions.
Definition: line_2d_3.h:63
This class defines the node.
Definition: node.h:65
double Y() const
Definition: point.h:187
double X() const
Definition: point.h:181
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 nine node 2D quadrilateral geometry with quadratic shape functions.
Definition: quadrilateral_2d_9.h:49
Quadrilateral2D9 & operator=(Quadrilateral2D9< TOtherPointType > const &rOther)
Definition: quadrilateral_2d_9.h:335
BaseType::CoordinatesArrayType CoordinatesArrayType
Definition: quadrilateral_2d_9.h:89
BaseType::JacobiansType JacobiansType
Definition: quadrilateral_2d_9.h:154
SizeType PointsNumberInDirection(IndexType LocalDirectionIndex) const override
Returns number of points per direction.
Definition: quadrilateral_2d_9.h:381
Quadrilateral2D9(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)
Definition: quadrilateral_2d_9.h:207
GeometryData::KratosGeometryFamily GetGeometryFamily() const override
Definition: quadrilateral_2d_9.h:291
void PrintData(std::ostream &rOStream) const override
Definition: quadrilateral_2d_9.h:644
Quadrilateral2D9 & operator=(const Quadrilateral2D9 &rOther)
Definition: quadrilateral_2d_9.h:316
BaseType::Pointer Create(const IndexType rNewGeometryId, PointsArrayType const &rThisPoints) const override
Creates a new geometry pointer.
Definition: quadrilateral_2d_9.h:352
Vector & ShapeFunctionsValues(Vector &rResult, const CoordinatesArrayType &rCoordinates) const override
Definition: quadrilateral_2d_9.h:583
bool IsInside(const CoordinatesArrayType &rPoint, CoordinatesArrayType &rResult, const double Tolerance=std::numeric_limits< double >::epsilon()) const override
Definition: quadrilateral_2d_9.h:468
double Volume() const override
This method calculates and returns the volume of this geometry.
Definition: quadrilateral_2d_9.h:435
ShapeFunctionsThirdDerivativesType & ShapeFunctionsThirdDerivatives(ShapeFunctionsThirdDerivativesType &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_9.h:924
BaseType::PointsArrayType PointsArrayType
Definition: quadrilateral_2d_9.h:111
std::string Info() const override
Definition: quadrilateral_2d_9.h:619
BaseType::ShapeFunctionsSecondDerivativesType ShapeFunctionsSecondDerivativesType
Definition: quadrilateral_2d_9.h:169
BaseType::SizeType SizeType
Definition: quadrilateral_2d_9.h:105
BaseType::ShapeFunctionsValuesContainerType ShapeFunctionsValuesContainerType
Definition: quadrilateral_2d_9.h:140
BaseType::NormalType NormalType
Definition: quadrilateral_2d_9.h:183
virtual ShapeFunctionsGradientsType ShapeFunctionsLocalGradients()
Definition: quadrilateral_2d_9.h:683
ShapeFunctionsSecondDerivativesType & ShapeFunctionsSecondDerivatives(ShapeFunctionsSecondDerivativesType &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_9.h:833
SizeType EdgesNumber() const override
This method gives you number of all edges of this geometry.
Definition: quadrilateral_2d_9.h:496
double DomainSize() const override
Definition: quadrilateral_2d_9.h:455
GeometriesArrayType GenerateEdges() const override
This method gives you all edges of this geometry.
Definition: quadrilateral_2d_9.h:509
virtual ShapeFunctionsGradientsType ShapeFunctionsLocalGradients(IntegrationMethod ThisMethod)
Definition: quadrilateral_2d_9.h:662
KRATOS_CLASS_POINTER_DEFINITION(Quadrilateral2D9)
BaseType::IndexType IndexType
Definition: quadrilateral_2d_9.h:98
Quadrilateral2D9(Quadrilateral2D9< TOtherPointType > const &rOther)
Definition: quadrilateral_2d_9.h:280
BaseType::Pointer Create(const IndexType NewGeometryId, const BaseType &rGeometry) const override
Creates a new geometry pointer.
Definition: quadrilateral_2d_9.h:366
Quadrilateral2D9(const std::string &rGeometryName, const PointsArrayType &rThisPoints)
Constructor with Geometry Name.
Definition: quadrilateral_2d_9.h:246
double Length() const override
Definition: quadrilateral_2d_9.h:406
Quadrilateral2D9(IndexType GeometryId, const PointsArrayType &rThisPoints)
Constructor with Geometry Id.
Definition: quadrilateral_2d_9.h:237
Matrix & ShapeFunctionsLocalGradients(Matrix &rResult, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_9.h:708
Quadrilateral2D9(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)
Definition: quadrilateral_2d_9.h:189
BaseType::ShapeFunctionsGradientsType ShapeFunctionsGradientsType
Definition: quadrilateral_2d_9.h:161
virtual Matrix & ShapeFunctionsGradients(Matrix &rResult, PointType &rPoint)
Definition: quadrilateral_2d_9.h:787
GeometryData::IntegrationMethod IntegrationMethod
Definition: quadrilateral_2d_9.h:73
Quadrilateral2D9(const PointsArrayType &ThisPoints)
Definition: quadrilateral_2d_9.h:229
Geometry< TPointType > BaseType
Definition: quadrilateral_2d_9.h:58
Quadrilateral2D9(Quadrilateral2D9 const &rOther)
Definition: quadrilateral_2d_9.h:263
TPointType PointType
Definition: quadrilateral_2d_9.h:84
Matrix & PointsLocalCoordinates(Matrix &rResult) const override
Definition: quadrilateral_2d_9.h:754
double Area() const override
Definition: quadrilateral_2d_9.h:423
BaseType::GeometriesArrayType GeometriesArrayType
Definition: quadrilateral_2d_9.h:79
~Quadrilateral2D9() override
Definition: quadrilateral_2d_9.h:289
void PrintInfo(std::ostream &rOStream) const override
Definition: quadrilateral_2d_9.h:630
friend class Quadrilateral2D9
Definition: quadrilateral_2d_9.h:1287
GeometryData::KratosGeometryType GetGeometryType() const override
Definition: quadrilateral_2d_9.h:296
double ShapeFunctionValue(IndexType ShapeFunctionIndex, const CoordinatesArrayType &rPoint) const override
Definition: quadrilateral_2d_9.h:535
BaseType::ShapeFunctionsThirdDerivativesType ShapeFunctionsThirdDerivativesType
Definition: quadrilateral_2d_9.h:177
BaseType::IntegrationPointsContainerType IntegrationPointsContainerType
Definition: quadrilateral_2d_9.h:133
Line2D3< TPointType > EdgeType
Definition: quadrilateral_2d_9.h:63
BaseType::IntegrationPointType IntegrationPointType
Definition: quadrilateral_2d_9.h:118
BaseType::IntegrationPointsArrayType IntegrationPointsArrayType
Definition: quadrilateral_2d_9.h:126
BaseType::ShapeFunctionsLocalGradientsContainerType ShapeFunctionsLocalGradientsContainerType
Definition: quadrilateral_2d_9.h:147
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
#define KRATOS_WARNING(label)
Definition: logger.h:265
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
std::istream & operator>>(std::istream &rIStream, LinearMasterSlaveConstraint &rThis)
input stream function
const GeometryData BrepCurve< TContainerPointType, TContainerPointEmbeddedType >::msGeometryData & msGeometryDimension
Definition: brep_curve.h:483
T & noalias(T &TheMatrix)
Definition: amatrix_interface.h:484
std::ostream & operator<<(std::ostream &rOStream, const LinearMasterSlaveConstraint &rThis)
output stream function
Definition: linear_master_slave_constraint.h:432
const GeometryData Quadrilateral2D9< TPointType >::msGeometryData & msGeometryDimension(), Quadrilateral2D9< TPointType >::AllShapeFunctionsValues(), AllShapeFunctionsLocalGradients()
Definition: brep_curve.h:483
def load(f)
Definition: ode_solve.py:307
int j
Definition: quadrature.py:648
float temp
Definition: rotating_cone.py:85
integer i
Definition: TensorModule.f:17