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.
poro_element_utilities.hpp
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: Ignasi de Pouplana
11 //
12 
13 
14 #if !defined(KRATOS_PORO_ELEMENT_UTILITIES )
15 #define KRATOS_PORO_ELEMENT_UTILITIES
16 
17 // System includes
18 //#include <cmath>
19 
20 // Project includes
21 #include "utilities/math_utils.h"
22 #include "includes/element.h"
23 
24 // Application includes
25 #include "poromechanics_application_variables.h"
26 
27 namespace Kratos
28 {
29 
30 class PoroElementUtilities
31 {
32 
33 typedef std::size_t IndexType;
34 
35 public:
36 
37 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
38 
39  static inline void CalculateNuMatrix(BoundedMatrix<double,2,6>& rNu, const Matrix& Ncontainer, const unsigned int& GPoint)
40  {
41  //Triangle_2d_3
42  rNu(0,0) = Ncontainer(GPoint,0); rNu(0,2) = Ncontainer(GPoint,1); rNu(0,4) = Ncontainer(GPoint,2);
43  rNu(1,1) = Ncontainer(GPoint,0); rNu(1,3) = Ncontainer(GPoint,1); rNu(1,5) = Ncontainer(GPoint,2);
44  }
45 
46  //----------------------------------------------------------------------------------------
47 
48  static inline void CalculateNuMatrix(BoundedMatrix<double,2,8>& rNu, const Matrix& Ncontainer, const unsigned int& GPoint)
49  {
50  //Quadrilateral_2d_4
51  rNu(0,0) = Ncontainer(GPoint,0); rNu(0,2) = Ncontainer(GPoint,1); rNu(0,4) = Ncontainer(GPoint,2); rNu(0,6) = Ncontainer(GPoint,3);
52  rNu(1,1) = Ncontainer(GPoint,0); rNu(1,3) = Ncontainer(GPoint,1); rNu(1,5) = Ncontainer(GPoint,2); rNu(1,7) = Ncontainer(GPoint,3);
53  }
54 
55  //----------------------------------------------------------------------------------------
56 
57  static inline void CalculateNuMatrix(BoundedMatrix<double,3,12>& rNu, const Matrix& Ncontainer, const unsigned int& GPoint)
58  {
59  //Tetrahedra_3d_4
60  rNu(0,0) = Ncontainer(GPoint,0); rNu(0,3) = Ncontainer(GPoint,1); rNu(0,6) = Ncontainer(GPoint,2); rNu(0,9) = Ncontainer(GPoint,3);
61  rNu(1,1) = Ncontainer(GPoint,0); rNu(1,4) = Ncontainer(GPoint,1); rNu(1,7) = Ncontainer(GPoint,2); rNu(1,10) = Ncontainer(GPoint,3);
62  rNu(2,2) = Ncontainer(GPoint,0); rNu(2,5) = Ncontainer(GPoint,1); rNu(2,8) = Ncontainer(GPoint,2); rNu(2,11) = Ncontainer(GPoint,3);
63  }
64 
65  //----------------------------------------------------------------------------------------
66 
67  static inline void CalculateNuMatrix(BoundedMatrix<double,3,18>& rNu, const Matrix& Ncontainer, const unsigned int& GPoint)
68  {
69  //Prism_3d_6
70  rNu(0,0) = Ncontainer(GPoint,0); rNu(0,3) = Ncontainer(GPoint,1); rNu(0,6) = Ncontainer(GPoint,2);
71  rNu(1,1) = Ncontainer(GPoint,0); rNu(1,4) = Ncontainer(GPoint,1); rNu(1,7) = Ncontainer(GPoint,2);
72  rNu(2,2) = Ncontainer(GPoint,0); rNu(2,5) = Ncontainer(GPoint,1); rNu(2,8) = Ncontainer(GPoint,2);
73 
74  rNu(0,9) = Ncontainer(GPoint,3); rNu(0,12) = Ncontainer(GPoint,4); rNu(0,15) = Ncontainer(GPoint,5);
75  rNu(1,10) = Ncontainer(GPoint,3); rNu(1,13) = Ncontainer(GPoint,4); rNu(1,16) = Ncontainer(GPoint,5);
76  rNu(2,11) = Ncontainer(GPoint,3); rNu(2,14) = Ncontainer(GPoint,4); rNu(2,17) = Ncontainer(GPoint,5);
77  }
78 
79  //----------------------------------------------------------------------------------------
80 
81  static inline void CalculateNuMatrix(BoundedMatrix<double,3,24>& rNu, const Matrix& Ncontainer, const unsigned int& GPoint)
82  {
83  //Hexahedron_3d_8
84  rNu(0,0) = Ncontainer(GPoint,0); rNu(0,3) = Ncontainer(GPoint,1); rNu(0,6) = Ncontainer(GPoint,2); rNu(0,9) = Ncontainer(GPoint,3);
85  rNu(1,1) = Ncontainer(GPoint,0); rNu(1,4) = Ncontainer(GPoint,1); rNu(1,7) = Ncontainer(GPoint,2); rNu(1,10) = Ncontainer(GPoint,3);
86  rNu(2,2) = Ncontainer(GPoint,0); rNu(2,5) = Ncontainer(GPoint,1); rNu(2,8) = Ncontainer(GPoint,2); rNu(2,11) = Ncontainer(GPoint,3);
87 
88  rNu(0,12) = Ncontainer(GPoint,4); rNu(0,15) = Ncontainer(GPoint,5); rNu(0,18) = Ncontainer(GPoint,6); rNu(0,21) = Ncontainer(GPoint,7);
89  rNu(1,13) = Ncontainer(GPoint,4); rNu(1,16) = Ncontainer(GPoint,5); rNu(1,19) = Ncontainer(GPoint,6); rNu(1,22) = Ncontainer(GPoint,7);
90  rNu(2,14) = Ncontainer(GPoint,4); rNu(2,17) = Ncontainer(GPoint,5); rNu(2,20) = Ncontainer(GPoint,6); rNu(2,23) = Ncontainer(GPoint,7);
91  }
92 
93 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
94 
95  static inline void CalculateNuElementMatrix(BoundedMatrix<double,3,9>& rNut, const Matrix& Ncontainer, const unsigned int& GPoint)
96  {
97  //Triangle_2d_3
98  rNut(0,0) = Ncontainer(GPoint,0); rNut(0,3) = Ncontainer(GPoint,1); rNut(0,6) = Ncontainer(GPoint,2);
99  rNut(1,1) = Ncontainer(GPoint,0); rNut(1,4) = Ncontainer(GPoint,1); rNut(1,7) = Ncontainer(GPoint,2);
100  }
101 
102  //----------------------------------------------------------------------------------------
103 
104  static inline void CalculateNuElementMatrix(BoundedMatrix<double,3,12>& rNut, const Matrix& Ncontainer, const unsigned int& GPoint)
105  {
106  //Quadrilateral_2d_4
107  rNut(0,0) = Ncontainer(GPoint,0); rNut(0,3) = Ncontainer(GPoint,1); rNut(0,6) = Ncontainer(GPoint,2); rNut(0,9) = Ncontainer(GPoint,3);
108  rNut(1,1) = Ncontainer(GPoint,0); rNut(1,4) = Ncontainer(GPoint,1); rNut(1,7) = Ncontainer(GPoint,2); rNut(1,10) = Ncontainer(GPoint,3);
109  }
110 
111  //----------------------------------------------------------------------------------------
112 
113  static inline void CalculateNuElementMatrix(BoundedMatrix<double,4,16>& rNut, const Matrix& Ncontainer, const unsigned int& GPoint)
114  {
115  //Tetrahedra_3d_4
116  rNut(0,0) = Ncontainer(GPoint,0); rNut(0,4) = Ncontainer(GPoint,1); rNut(0,8) = Ncontainer(GPoint,2); rNut(0,12) = Ncontainer(GPoint,3);
117  rNut(1,1) = Ncontainer(GPoint,0); rNut(1,5) = Ncontainer(GPoint,1); rNut(1,9) = Ncontainer(GPoint,2); rNut(1,13) = Ncontainer(GPoint,3);
118  rNut(2,2) = Ncontainer(GPoint,0); rNut(2,6) = Ncontainer(GPoint,1); rNut(2,10) = Ncontainer(GPoint,2); rNut(2,14) = Ncontainer(GPoint,3);
119  }
120 
121  //----------------------------------------------------------------------------------------
122 
123  static inline void CalculateNuElementMatrix(BoundedMatrix<double,4,24>& rNut, const Matrix& Ncontainer, const unsigned int& GPoint)
124  {
125  //Prism_3d_6
126  rNut(0,0) = Ncontainer(GPoint,0); rNut(0,4) = Ncontainer(GPoint,1); rNut(0,8) = Ncontainer(GPoint,2);
127  rNut(1,1) = Ncontainer(GPoint,0); rNut(1,5) = Ncontainer(GPoint,1); rNut(1,9) = Ncontainer(GPoint,2);
128  rNut(2,2) = Ncontainer(GPoint,0); rNut(2,6) = Ncontainer(GPoint,1); rNut(2,10) = Ncontainer(GPoint,2);
129 
130  rNut(0,12) = Ncontainer(GPoint,3); rNut(0,16) = Ncontainer(GPoint,4); rNut(0,20) = Ncontainer(GPoint,5);
131  rNut(1,13) = Ncontainer(GPoint,3); rNut(1,17) = Ncontainer(GPoint,4); rNut(1,21) = Ncontainer(GPoint,5);
132  rNut(2,14) = Ncontainer(GPoint,3); rNut(2,18) = Ncontainer(GPoint,4); rNut(2,22) = Ncontainer(GPoint,5);
133  }
134 
135  //----------------------------------------------------------------------------------------
136 
137  static inline void CalculateNuElementMatrix(BoundedMatrix<double,4,32>& rNut, const Matrix& Ncontainer, const unsigned int& GPoint)
138  {
139  //Hexahedron_3d_8
140  rNut(0,0) = Ncontainer(GPoint,0); rNut(0,4) = Ncontainer(GPoint,1); rNut(0,8) = Ncontainer(GPoint,2); rNut(0,12) = Ncontainer(GPoint,3);
141  rNut(1,1) = Ncontainer(GPoint,0); rNut(1,5) = Ncontainer(GPoint,1); rNut(1,9) = Ncontainer(GPoint,2); rNut(1,13) = Ncontainer(GPoint,3);
142  rNut(2,2) = Ncontainer(GPoint,0); rNut(2,6) = Ncontainer(GPoint,1); rNut(2,10) = Ncontainer(GPoint,2); rNut(2,14) = Ncontainer(GPoint,3);
143 
144  rNut(0,16) = Ncontainer(GPoint,4); rNut(0,20) = Ncontainer(GPoint,5); rNut(0,24) = Ncontainer(GPoint,6); rNut(0,28) = Ncontainer(GPoint,7);
145  rNut(1,17) = Ncontainer(GPoint,4); rNut(1,21) = Ncontainer(GPoint,5); rNut(1,25) = Ncontainer(GPoint,6); rNut(1,29) = Ncontainer(GPoint,7);
146  rNut(2,18) = Ncontainer(GPoint,4); rNut(2,22) = Ncontainer(GPoint,5); rNut(2,26) = Ncontainer(GPoint,6); rNut(2,30) = Ncontainer(GPoint,7);
147  }
148 
149 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
150 
151  static inline void InterpolateVariableWithComponents(array_1d<double,2>& rVector,const Matrix& Ncontainer,
152  const array_1d<double,6>& VariableWithComponents,const unsigned int& GPoint)
153  {
154  //Triangle_2d_3
155  noalias(rVector) = ZeroVector(2);
156 
157  unsigned int index = 0;
158  for(unsigned int i=0; i<3; i++)
159  {
160  rVector[0] += Ncontainer(GPoint,i)*VariableWithComponents[index++];
161  rVector[1] += Ncontainer(GPoint,i)*VariableWithComponents[index++];
162  }
163  }
164 
165  //----------------------------------------------------------------------------------------
166 
167  static inline void InterpolateVariableWithComponents(array_1d<double,2>& rVector,const Matrix& Ncontainer,
168  const array_1d<double,8>& VariableWithComponents,const unsigned int& GPoint)
169  {
170  //Quadrilateral_2d_4
171  noalias(rVector) = ZeroVector(2);
172 
173  unsigned int index = 0;
174  for(unsigned int i=0; i<4; i++)
175  {
176  rVector[0] += Ncontainer(GPoint,i)*VariableWithComponents[index++];
177  rVector[1] += Ncontainer(GPoint,i)*VariableWithComponents[index++];
178  }
179  }
180 
181  //----------------------------------------------------------------------------------------
182 
183  static inline void InterpolateVariableWithComponents(array_1d<double,3>& rVector,const Matrix& Ncontainer,
184  const array_1d<double,12>& VariableWithComponents,const unsigned int& GPoint)
185  {
186  //Tetrahedra_3d_4
187  noalias(rVector) = ZeroVector(3);
188 
189  unsigned int index = 0;
190  for(unsigned int i=0; i<4; i++)
191  {
192  rVector[0] += Ncontainer(GPoint,i)*VariableWithComponents[index++];
193  rVector[1] += Ncontainer(GPoint,i)*VariableWithComponents[index++];
194  rVector[2] += Ncontainer(GPoint,i)*VariableWithComponents[index++];
195  }
196  }
197 
198  //----------------------------------------------------------------------------------------
199 
200  static inline void InterpolateVariableWithComponents(array_1d<double,3>& rVector,const Matrix& Ncontainer,
201  const array_1d<double,18>& VariableWithComponents,const unsigned int& GPoint)
202  {
203  //Prism_3d_6
204  noalias(rVector) = ZeroVector(3);
205 
206  unsigned int index = 0;
207  for(unsigned int i=0; i<6; i++)
208  {
209  rVector[0] += Ncontainer(GPoint,i)*VariableWithComponents[index++];
210  rVector[1] += Ncontainer(GPoint,i)*VariableWithComponents[index++];
211  rVector[2] += Ncontainer(GPoint,i)*VariableWithComponents[index++];
212  }
213  }
214 
215  //----------------------------------------------------------------------------------------
216 
217  static inline void InterpolateVariableWithComponents(array_1d<double,3>& rVector,const Matrix& Ncontainer,
218  const array_1d<double,24>& VariableWithComponents,const unsigned int& GPoint)
219  {
220  //Hexahedra_3d_8
221  noalias(rVector) = ZeroVector(3);
222 
223  unsigned int index = 0;
224  for(unsigned int i=0; i<8; i++)
225  {
226  rVector[0] += Ncontainer(GPoint,i)*VariableWithComponents[index++];
227  rVector[1] += Ncontainer(GPoint,i)*VariableWithComponents[index++];
228  rVector[2] += Ncontainer(GPoint,i)*VariableWithComponents[index++];
229  }
230  }
231 
232 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
233 
234  static inline void FillArray1dOutput(array_1d<double,3>& rOutputValue, const array_1d<double,2>& ComputedValue)
235  {
236  rOutputValue[0] = ComputedValue[0];
237  rOutputValue[1] = ComputedValue[1];
238  rOutputValue[2] = 0.0;
239  }
240 
241  //----------------------------------------------------------------------------------------
242 
243  static inline void FillArray1dOutput(array_1d<double,3>& rOutputValue, const array_1d<double,3>& ComputedValue)
244  {
245  rOutputValue[0] = ComputedValue[0];
246  rOutputValue[1] = ComputedValue[1];
247  rOutputValue[2] = ComputedValue[2];
248  }
249 
250 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
251 
252  static inline void GetNodalVariableVector(array_1d<double,6>& rNodalVariableVector, const Element::GeometryType& Geom,
253  const Variable<array_1d<double,3>>& Variable, IndexType SolutionStepIndex = 0)
254  {
255  //Triangle_2d_3
256  array_1d<double,3> NodalVariableAux;
257  unsigned int index = 0;
258  for(unsigned int i=0; i<3; i++)
259  {
260  noalias(NodalVariableAux) = Geom[i].FastGetSolutionStepValue(Variable,SolutionStepIndex);
261  rNodalVariableVector[index++] = NodalVariableAux[0];
262  rNodalVariableVector[index++] = NodalVariableAux[1];
263  }
264  }
265 
266  //----------------------------------------------------------------------------------------
267 
268  static inline void GetNodalVariableVector(array_1d<double,8>& rNodalVariableVector, const Element::GeometryType& Geom,
269  const Variable<array_1d<double,3>>& Variable, IndexType SolutionStepIndex = 0)
270  {
271  //Quadrilateral_2d_4
272  array_1d<double,3> NodalVariableAux;
273  unsigned int index = 0;
274  for(unsigned int i=0; i<4; i++)
275  {
276  noalias(NodalVariableAux) = Geom[i].FastGetSolutionStepValue(Variable,SolutionStepIndex);
277  rNodalVariableVector[index++] = NodalVariableAux[0];
278  rNodalVariableVector[index++] = NodalVariableAux[1];
279  }
280  }
281 
282  //----------------------------------------------------------------------------------------
283 
284  static inline void GetNodalVariableVector(array_1d<double,12>& rNodalVariableVector, const Element::GeometryType& Geom,
285  const Variable<array_1d<double,3>>& Variable, IndexType SolutionStepIndex = 0)
286  {
287  //Tetrahedra_3d_4
288  array_1d<double,3> NodalVariableAux;
289  unsigned int index = 0;
290  for(unsigned int i=0; i<4; i++)
291  {
292  noalias(NodalVariableAux) = Geom[i].FastGetSolutionStepValue(Variable,SolutionStepIndex);
293  rNodalVariableVector[index++] = NodalVariableAux[0];
294  rNodalVariableVector[index++] = NodalVariableAux[1];
295  rNodalVariableVector[index++] = NodalVariableAux[2];
296  }
297  }
298 
299  //----------------------------------------------------------------------------------------
300 
301  static inline void GetNodalVariableVector(array_1d<double,18>& rNodalVariableVector, const Element::GeometryType& Geom,
302  const Variable<array_1d<double,3>>& Variable, IndexType SolutionStepIndex = 0)
303  {
304  //Prism_3d_6
305  array_1d<double,3> NodalVariableAux;
306  unsigned int index = 0;
307  for(unsigned int i=0; i<6; i++)
308  {
309  noalias(NodalVariableAux) = Geom[i].FastGetSolutionStepValue(Variable,SolutionStepIndex);
310  rNodalVariableVector[index++] = NodalVariableAux[0];
311  rNodalVariableVector[index++] = NodalVariableAux[1];
312  rNodalVariableVector[index++] = NodalVariableAux[2];
313  }
314  }
315 
316  //----------------------------------------------------------------------------------------
317 
318  static inline void GetNodalVariableVector(array_1d<double,24>& rNodalVariableVector, const Element::GeometryType& Geom,
319  const Variable<array_1d<double,3>>& Variable, IndexType SolutionStepIndex = 0)
320  {
321  //Hexahedron_3d_8
322  array_1d<double,3> NodalVariableAux;
323  unsigned int index = 0;
324  for(unsigned int i=0; i<8; i++)
325  {
326  noalias(NodalVariableAux) = Geom[i].FastGetSolutionStepValue(Variable,SolutionStepIndex);
327  rNodalVariableVector[index++] = NodalVariableAux[0];
328  rNodalVariableVector[index++] = NodalVariableAux[1];
329  rNodalVariableVector[index++] = NodalVariableAux[2];
330  }
331  }
332 
333 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
334 
335  static inline void CalculatePermeabilityMatrix(BoundedMatrix<double,2,2>& rPermeabilityMatrix,
336  const Element::PropertiesType& Prop)
337  {
338  //2D
339  rPermeabilityMatrix(0,0) = Prop[PERMEABILITY_XX];
340  rPermeabilityMatrix(1,1) = Prop[PERMEABILITY_YY];
341 
342  rPermeabilityMatrix(0,1) = Prop[PERMEABILITY_XY];
343  rPermeabilityMatrix(1,0) = rPermeabilityMatrix(0,1);
344  }
345 
346  //----------------------------------------------------------------------------------------
347 
348  static inline void CalculatePermeabilityMatrix(BoundedMatrix<double,3,3>& rPermeabilityMatrix,
349  const Element::PropertiesType& Prop)
350  {
351  //3D
352  rPermeabilityMatrix(0,0) = Prop[PERMEABILITY_XX];
353  rPermeabilityMatrix(1,1) = Prop[PERMEABILITY_YY];
354  rPermeabilityMatrix(2,2) = Prop[PERMEABILITY_ZZ];
355 
356  rPermeabilityMatrix(0,1) = Prop[PERMEABILITY_XY];
357  rPermeabilityMatrix(1,0) = rPermeabilityMatrix(0,1);
358 
359  rPermeabilityMatrix(1,2) = Prop[PERMEABILITY_YZ];
360  rPermeabilityMatrix(2,1) = rPermeabilityMatrix(1,2);
361 
362  rPermeabilityMatrix(2,0) = Prop[PERMEABILITY_ZX];
363  rPermeabilityMatrix(0,2) = rPermeabilityMatrix(2,0);
364  }
365 
366  //----------------------------------------------------------------------------------------
367 
368  static inline void CalculatePermeabilityMatrix(Matrix& rPermeabilityMatrix,
369  const Element::PropertiesType& Prop,
370  const unsigned int& Dim)
371  {
372  if ( rPermeabilityMatrix.size1() != Dim )
373  rPermeabilityMatrix.resize( Dim, Dim, false );
374 
375  rPermeabilityMatrix(0,0) = Prop[PERMEABILITY_XX];
376  rPermeabilityMatrix(1,1) = Prop[PERMEABILITY_YY];
377  rPermeabilityMatrix(0,1) = Prop[PERMEABILITY_XY];
378  rPermeabilityMatrix(1,0) = rPermeabilityMatrix(0,1);
379  if(Dim==3)
380  {
381  rPermeabilityMatrix(2,2) = Prop[PERMEABILITY_ZZ];
382  rPermeabilityMatrix(2,0) = Prop[PERMEABILITY_ZX];
383  rPermeabilityMatrix(1,2) = Prop[PERMEABILITY_YZ];
384  rPermeabilityMatrix(0,2) = rPermeabilityMatrix(2,0);
385  rPermeabilityMatrix(2,1) = rPermeabilityMatrix(1,2);
386  }
387  }
388 
389 
390 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
391 
392  static inline void InvertMatrix2(BoundedMatrix<double,2,2>& rInvertedMatrix,
393  const BoundedMatrix<double,2,2>& InputMatrix)
394  {
395  double InputMatrixDet = InputMatrix(0,0)*InputMatrix(1,1)-InputMatrix(0,1)*InputMatrix(1,0);
396 
397  rInvertedMatrix(0,0) = InputMatrix(1,1)/InputMatrixDet;
398  rInvertedMatrix(0,1) = -InputMatrix(0,1)/InputMatrixDet;
399  rInvertedMatrix(1,0) = -InputMatrix(1,0)/InputMatrixDet;
400  rInvertedMatrix(1,1) = InputMatrix(0,0)/InputMatrixDet;
401  }
402 
403  //----------------------------------------------------------------------------------------
404 /*
405  template<class T>
406  bool InvertMatrix(const T& input, T& inverse)
407  {
408  typedef permutation_matrix<std::size_t> pmatrix;
409 
410  // create a working copy of the input
411  T A(input);
412 
413  // create a permutation matrix for the LU-factorization
414  pmatrix pm(A.size1());
415 
416  // perform LU-factorization
417  int res = lu_factorize(A, pm);
418  if (res != 0)
419  return false;
420 
421  // create identity matrix of "inverse"
422  inverse.assign(identity_matrix<double> (A.size1()));
423 
424  // backsubstitute to get the inverse
425  lu_substitute(A, pm, inverse);
426 
427  return true;
428  }
429 */
430 
431 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
432 
433  static inline void AssembleUBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,6,6>& UBlockMatrix)
434  {
435  //Triangle_2d_3
436  unsigned int Global_i, Global_j, Local_i, Local_j;
437 
438  for(unsigned int i = 0; i < 3; i++)
439  {
440  Global_i = i * (2 + 1);
441  Local_i = i * 2;
442 
443  for(unsigned int j = 0; j < 3; j++)
444  {
445  Global_j = j * (2 + 1);
446  Local_j = j * 2;
447 
448  rLeftHandSideMatrix(Global_i,Global_j) += UBlockMatrix(Local_i,Local_j);
449  rLeftHandSideMatrix(Global_i,Global_j+1) += UBlockMatrix(Local_i,Local_j+1);
450  rLeftHandSideMatrix(Global_i+1,Global_j) += UBlockMatrix(Local_i+1,Local_j);
451  rLeftHandSideMatrix(Global_i+1,Global_j+1) += UBlockMatrix(Local_i+1,Local_j+1);
452  }
453  }
454  }
455 
456  //----------------------------------------------------------------------------------------
457 
458  static inline void AssembleUBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,8,8>& UBlockMatrix)
459  {
460  //Quadrilateral_2d_4
461  unsigned int Global_i, Global_j, Local_i, Local_j;
462 
463  for(unsigned int i = 0; i < 4; i++)
464  {
465  Global_i = i * (2 + 1);
466  Local_i = i * 2;
467 
468  for(unsigned int j = 0; j < 4; j++)
469  {
470  Global_j = j * (2 + 1);
471  Local_j = j * 2;
472 
473  rLeftHandSideMatrix(Global_i,Global_j) += UBlockMatrix(Local_i,Local_j);
474  rLeftHandSideMatrix(Global_i,Global_j+1) += UBlockMatrix(Local_i,Local_j+1);
475  rLeftHandSideMatrix(Global_i+1,Global_j) += UBlockMatrix(Local_i+1,Local_j);
476  rLeftHandSideMatrix(Global_i+1,Global_j+1) += UBlockMatrix(Local_i+1,Local_j+1);
477  }
478  }
479  }
480 
481  //----------------------------------------------------------------------------------------
482 
483  static inline void AssembleUBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,12,12>& UBlockMatrix)
484  {
485  //Tetrahedra_3d_4
486  unsigned int Global_i, Global_j, Local_i, Local_j;
487 
488  for(unsigned int i = 0; i < 4; i++)
489  {
490  Global_i = i * (3 + 1);
491  Local_i = i * 3;
492 
493  for(unsigned int j = 0; j < 4; j++)
494  {
495  Global_j = j * (3 + 1);
496  Local_j = j * 3;
497 
498  rLeftHandSideMatrix(Global_i,Global_j) += UBlockMatrix(Local_i,Local_j);
499  rLeftHandSideMatrix(Global_i,Global_j+1) += UBlockMatrix(Local_i,Local_j+1);
500  rLeftHandSideMatrix(Global_i+1,Global_j) += UBlockMatrix(Local_i+1,Local_j);
501  rLeftHandSideMatrix(Global_i+1,Global_j+1) += UBlockMatrix(Local_i+1,Local_j+1);
502 
503  rLeftHandSideMatrix(Global_i,Global_j+2) += UBlockMatrix(Local_i,Local_j+2);
504  rLeftHandSideMatrix(Global_i+1,Global_j+2) += UBlockMatrix(Local_i+1,Local_j+2);
505  rLeftHandSideMatrix(Global_i+2,Global_j+1) += UBlockMatrix(Local_i+2,Local_j+1);
506  rLeftHandSideMatrix(Global_i+2,Global_j) += UBlockMatrix(Local_i+2,Local_j);
507  rLeftHandSideMatrix(Global_i+2,Global_j+2) += UBlockMatrix(Local_i+2,Local_j+2);
508  }
509  }
510  }
511 
512  //----------------------------------------------------------------------------------------
513 
514  static inline void AssembleUBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,18,18>& UBlockMatrix)
515  {
516  //Prism_3d_6
517  unsigned int Global_i, Global_j, Local_i, Local_j;
518 
519  for(unsigned int i = 0; i < 6; i++)
520  {
521  Global_i = i * (3 + 1);
522  Local_i = i * 3;
523 
524  for(unsigned int j = 0; j < 6; j++)
525  {
526  Global_j = j * (3 + 1);
527  Local_j = j * 3;
528 
529  rLeftHandSideMatrix(Global_i,Global_j) += UBlockMatrix(Local_i,Local_j);
530  rLeftHandSideMatrix(Global_i,Global_j+1) += UBlockMatrix(Local_i,Local_j+1);
531  rLeftHandSideMatrix(Global_i+1,Global_j) += UBlockMatrix(Local_i+1,Local_j);
532  rLeftHandSideMatrix(Global_i+1,Global_j+1) += UBlockMatrix(Local_i+1,Local_j+1);
533 
534  rLeftHandSideMatrix(Global_i,Global_j+2) += UBlockMatrix(Local_i,Local_j+2);
535  rLeftHandSideMatrix(Global_i+1,Global_j+2) += UBlockMatrix(Local_i+1,Local_j+2);
536  rLeftHandSideMatrix(Global_i+2,Global_j+1) += UBlockMatrix(Local_i+2,Local_j+1);
537  rLeftHandSideMatrix(Global_i+2,Global_j) += UBlockMatrix(Local_i+2,Local_j);
538  rLeftHandSideMatrix(Global_i+2,Global_j+2) += UBlockMatrix(Local_i+2,Local_j+2);
539  }
540  }
541  }
542 
543  //----------------------------------------------------------------------------------------
544 
545  static inline void AssembleUBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,24,24>& UBlockMatrix)
546  {
547  //Hexahedra_3d_8
548  unsigned int Global_i, Global_j, Local_i, Local_j;
549 
550  for(unsigned int i = 0; i < 8; i++)
551  {
552  Global_i = i * (3 + 1);
553  Local_i = i * 3;
554 
555  for(unsigned int j = 0; j < 8; j++)
556  {
557  Global_j = j * (3 + 1);
558  Local_j = j * 3;
559 
560  rLeftHandSideMatrix(Global_i,Global_j) += UBlockMatrix(Local_i,Local_j);
561  rLeftHandSideMatrix(Global_i,Global_j+1) += UBlockMatrix(Local_i,Local_j+1);
562  rLeftHandSideMatrix(Global_i+1,Global_j) += UBlockMatrix(Local_i+1,Local_j);
563  rLeftHandSideMatrix(Global_i+1,Global_j+1) += UBlockMatrix(Local_i+1,Local_j+1);
564 
565  rLeftHandSideMatrix(Global_i,Global_j+2) += UBlockMatrix(Local_i,Local_j+2);
566  rLeftHandSideMatrix(Global_i+1,Global_j+2) += UBlockMatrix(Local_i+1,Local_j+2);
567  rLeftHandSideMatrix(Global_i+2,Global_j+1) += UBlockMatrix(Local_i+2,Local_j+1);
568  rLeftHandSideMatrix(Global_i+2,Global_j) += UBlockMatrix(Local_i+2,Local_j);
569  rLeftHandSideMatrix(Global_i+2,Global_j+2) += UBlockMatrix(Local_i+2,Local_j+2);
570  }
571  }
572  }
573 
574  //----------------------------------------------------------------------------------------
575 
576  static inline void AssembleUPBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,6,3>& UPBlockMatrix)
577  {
578  //Triangle_2d_3
579  unsigned int Global_i, Global_j, Local_i;
580 
581  for(unsigned int i = 0; i < 3; i++)
582  {
583  Global_i = i * (2 + 1);
584  Local_i = i * 2;
585 
586  for(unsigned int j = 0; j < 3; j++)
587  {
588  Global_j = j * (2 + 1) + 2;
589 
590  rLeftHandSideMatrix(Global_i,Global_j) += UPBlockMatrix(Local_i,j);
591  rLeftHandSideMatrix(Global_i+1,Global_j) += UPBlockMatrix(Local_i+1,j);
592  }
593  }
594  }
595 
596  //----------------------------------------------------------------------------------------
597 
598  static inline void AssembleUPBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,8,4>& UPBlockMatrix)
599  {
600  //Quadrilateral_2d_4
601  unsigned int Global_i, Global_j, Local_i;
602 
603  for(unsigned int i = 0; i < 4; i++)
604  {
605  Global_i = i * (2 + 1);
606  Local_i = i * 2;
607 
608  for(unsigned int j = 0; j < 4; j++)
609  {
610  Global_j = j * (2 + 1) + 2;
611 
612  rLeftHandSideMatrix(Global_i,Global_j) += UPBlockMatrix(Local_i,j);
613  rLeftHandSideMatrix(Global_i+1,Global_j) += UPBlockMatrix(Local_i+1,j);
614  }
615  }
616  }
617 
618  //----------------------------------------------------------------------------------------
619 
620  static inline void AssembleUPBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,12,4>& UPBlockMatrix)
621  {
622  //Tetrahedra_3d_4
623  unsigned int Global_i, Global_j, Local_i;
624 
625  for(unsigned int i = 0; i < 4; i++)
626  {
627  Global_i = i * (3 + 1);
628  Local_i = i * 3;
629 
630  for(unsigned int j = 0; j < 4; j++)
631  {
632  Global_j = j * (3 + 1) + 3;
633 
634  rLeftHandSideMatrix(Global_i,Global_j) += UPBlockMatrix(Local_i,j);
635  rLeftHandSideMatrix(Global_i+1,Global_j) += UPBlockMatrix(Local_i+1,j);
636  rLeftHandSideMatrix(Global_i+2,Global_j) += UPBlockMatrix(Local_i+2,j);
637  }
638  }
639  }
640 
641  //----------------------------------------------------------------------------------------
642 
643  static inline void AssembleUPBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,18,6>& UPBlockMatrix)
644  {
645  //Prism_3d_6
646  unsigned int Global_i, Global_j, Local_i;
647 
648  for(unsigned int i = 0; i < 6; i++)
649  {
650  Global_i = i * (3 + 1);
651  Local_i = i * 3;
652 
653  for(unsigned int j = 0; j < 6; j++)
654  {
655  Global_j = j * (3 + 1) + 3;
656 
657  rLeftHandSideMatrix(Global_i,Global_j) += UPBlockMatrix(Local_i,j);
658  rLeftHandSideMatrix(Global_i+1,Global_j) += UPBlockMatrix(Local_i+1,j);
659  rLeftHandSideMatrix(Global_i+2,Global_j) += UPBlockMatrix(Local_i+2,j);
660  }
661  }
662  }
663 
664  //----------------------------------------------------------------------------------------
665 
666  static inline void AssembleUPBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,24,8>& UPBlockMatrix)
667  {
668  //Hexahedra_3d_8
669  unsigned int Global_i, Global_j, Local_i;
670 
671  for(unsigned int i = 0; i < 8; i++)
672  {
673  Global_i = i * (3 + 1);
674  Local_i = i * 3;
675 
676  for(unsigned int j = 0; j < 8; j++)
677  {
678  Global_j = j * (3 + 1) + 3;
679 
680  rLeftHandSideMatrix(Global_i,Global_j) += UPBlockMatrix(Local_i,j);
681  rLeftHandSideMatrix(Global_i+1,Global_j) += UPBlockMatrix(Local_i+1,j);
682  rLeftHandSideMatrix(Global_i+2,Global_j) += UPBlockMatrix(Local_i+2,j);
683  }
684  }
685  }
686 
687  //----------------------------------------------------------------------------------------
688 
689  static inline void AssemblePUBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,3,6>& PUBlockMatrix)
690  {
691  //Triangle_2d_3
692  unsigned int Global_i, Global_j, Local_j;
693 
694  for(unsigned int i = 0; i < 3; i++)
695  {
696  Global_i = i * (2 + 1) + 2;
697 
698  for(unsigned int j = 0; j < 3; j++)
699  {
700  Global_j = j * (2 + 1);
701  Local_j = j * 2;
702 
703  rLeftHandSideMatrix(Global_i,Global_j) += PUBlockMatrix(i,Local_j);
704  rLeftHandSideMatrix(Global_i,Global_j+1) += PUBlockMatrix(i,Local_j+1);
705  }
706  }
707  }
708 
709  //----------------------------------------------------------------------------------------
710 
711  static inline void AssemblePUBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,4,8>& PUBlockMatrix)
712  {
713  //Quadrilateral_2d_4
714  unsigned int Global_i, Global_j, Local_j;
715 
716  for(unsigned int i = 0; i < 4; i++)
717  {
718  Global_i = i * (2 + 1) + 2;
719 
720  for(unsigned int j = 0; j < 4; j++)
721  {
722  Global_j = j * (2 + 1);
723  Local_j = j * 2;
724 
725  rLeftHandSideMatrix(Global_i,Global_j) += PUBlockMatrix(i,Local_j);
726  rLeftHandSideMatrix(Global_i,Global_j+1) += PUBlockMatrix(i,Local_j+1);
727  }
728  }
729  }
730 
731  //----------------------------------------------------------------------------------------
732 
733  static inline void AssemblePUBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,4,12>& PUBlockMatrix)
734  {
735  //Tetrahedra_3d_4
736  unsigned int Global_i, Global_j, Local_j;
737 
738  for(unsigned int i = 0; i < 4; i++)
739  {
740  Global_i = i * (3 + 1) + 3;
741 
742  for(unsigned int j = 0; j < 4; j++)
743  {
744  Global_j = j * (3 + 1);
745  Local_j = j * 3;
746 
747  rLeftHandSideMatrix(Global_i,Global_j) += PUBlockMatrix(i,Local_j);
748  rLeftHandSideMatrix(Global_i,Global_j+1) += PUBlockMatrix(i,Local_j+1);
749  rLeftHandSideMatrix(Global_i,Global_j+2) += PUBlockMatrix(i,Local_j+2);
750  }
751  }
752  }
753 
754  //----------------------------------------------------------------------------------------
755 
756  static inline void AssemblePUBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,6,18>& PUBlockMatrix)
757  {
758  //Prism_3d_6
759  unsigned int Global_i, Global_j, Local_j;
760 
761  for(unsigned int i = 0; i < 6; i++)
762  {
763  Global_i = i * (3 + 1) + 3;
764 
765  for(unsigned int j = 0; j < 6; j++)
766  {
767  Global_j = j * (3 + 1);
768  Local_j = j * 3;
769 
770  rLeftHandSideMatrix(Global_i,Global_j) += PUBlockMatrix(i,Local_j);
771  rLeftHandSideMatrix(Global_i,Global_j+1) += PUBlockMatrix(i,Local_j+1);
772  rLeftHandSideMatrix(Global_i,Global_j+2) += PUBlockMatrix(i,Local_j+2);
773  }
774  }
775  }
776 
777  //----------------------------------------------------------------------------------------
778 
779  static inline void AssemblePUBlockMatrix(Matrix& rLeftHandSideMatrix, const BoundedMatrix<double,8,24>& PUBlockMatrix)
780  {
781  //Hexahedra_3d_8
782  unsigned int Global_i, Global_j, Local_j;
783 
784  for(unsigned int i = 0; i < 8; i++)
785  {
786  Global_i = i * (3 + 1) + 3;
787 
788  for(unsigned int j = 0; j < 8; j++)
789  {
790  Global_j = j * (3 + 1);
791  Local_j = j * 3;
792 
793  rLeftHandSideMatrix(Global_i,Global_j) += PUBlockMatrix(i,Local_j);
794  rLeftHandSideMatrix(Global_i,Global_j+1) += PUBlockMatrix(i,Local_j+1);
795  rLeftHandSideMatrix(Global_i,Global_j+2) += PUBlockMatrix(i,Local_j+2);
796  }
797  }
798  }
799 
800  //----------------------------------------------------------------------------------------
801 
802  template< class TMatrixType >
803  static inline void AssemblePBlockMatrix(Matrix& rLeftHandSideMatrix,const TMatrixType& PBlockMatrix, const unsigned int& Dim, const unsigned int& NumNodes)
804  {
805  unsigned int Global_i, Global_j;
806 
807  for(unsigned int i = 0; i < NumNodes; i++)
808  {
809  Global_i = i * (Dim + 1) + Dim;
810 
811  for(unsigned int j = 0; j < NumNodes; j++)
812  {
813  Global_j = j * (Dim + 1) + Dim;
814 
815  rLeftHandSideMatrix(Global_i,Global_j) += PBlockMatrix(i,j);
816  }
817  }
818  }
819 
820  //----------------------------------------------------------------------------------------
821 
822  static inline void AssembleUBlockVector(Vector& rRightHandSideVector, const array_1d<double,6>& UBlockVector)
823  {
824  //Triangle_2d_3
825  unsigned int Global_i, Local_i;
826 
827  for(unsigned int i = 0; i < 3; i++)
828  {
829  Global_i = i * (2 + 1);
830  Local_i = i * 2;
831 
832  rRightHandSideVector[Global_i] += UBlockVector[Local_i];
833  rRightHandSideVector[Global_i+1] += UBlockVector[Local_i+1];
834  }
835  }
836 
837  //----------------------------------------------------------------------------------------
838 
839  static inline void AssembleUBlockVector(Vector& rRightHandSideVector, const array_1d<double,8>& UBlockVector)
840  {
841  //Quadrilateral_2d_4
842  unsigned int Global_i, Local_i;
843 
844  for(unsigned int i = 0; i < 4; i++)
845  {
846  Global_i = i * (2 + 1);
847  Local_i = i * 2;
848 
849  rRightHandSideVector[Global_i] += UBlockVector[Local_i];
850  rRightHandSideVector[Global_i+1] += UBlockVector[Local_i+1];
851  }
852  }
853 
854  //----------------------------------------------------------------------------------------
855 
856  static inline void AssembleUBlockVector(Vector& rRightHandSideVector, const array_1d<double,12>& UBlockVector)
857  {
858  //Tetrahedra_3d_4
859  unsigned int Global_i, Local_i;
860 
861  for(unsigned int i = 0; i < 4; i++)
862  {
863  Global_i = i * (3 + 1);
864  Local_i = i * 3;
865 
866  rRightHandSideVector[Global_i] += UBlockVector[Local_i];
867  rRightHandSideVector[Global_i+1] += UBlockVector[Local_i+1];
868  rRightHandSideVector[Global_i+2] += UBlockVector[Local_i+2];
869  }
870  }
871 
872  //----------------------------------------------------------------------------------------
873 
874  static inline void AssembleUBlockVector(Vector& rRightHandSideVector, const array_1d<double,18>& UBlockVector)
875  {
876  //Prism_3d_6
877  unsigned int Global_i, Local_i;
878 
879  for(unsigned int i = 0; i < 6; i++)
880  {
881  Global_i = i * (3 + 1);
882  Local_i = i * 3;
883 
884  rRightHandSideVector[Global_i] += UBlockVector[Local_i];
885  rRightHandSideVector[Global_i+1] += UBlockVector[Local_i+1];
886  rRightHandSideVector[Global_i+2] += UBlockVector[Local_i+2];
887  }
888  }
889 
890  //----------------------------------------------------------------------------------------
891 
892  static inline void AssembleUBlockVector(Vector& rRightHandSideVector, const array_1d<double,24>& UBlockVector)
893  {
894  //Hexahedra_3d_8
895  unsigned int Global_i, Local_i;
896 
897  for(unsigned int i = 0; i < 8; i++)
898  {
899  Global_i = i * (3 + 1);
900  Local_i = i * 3;
901 
902  rRightHandSideVector[Global_i] += UBlockVector[Local_i];
903  rRightHandSideVector[Global_i+1] += UBlockVector[Local_i+1];
904  rRightHandSideVector[Global_i+2] += UBlockVector[Local_i+2];
905  }
906  }
907 
908  //----------------------------------------------------------------------------------------
909 
910  template< class TVectorType >
911  static inline void AssemblePBlockVector(Vector& rRightHandSideVector,const TVectorType& PBlockVector, const unsigned int& Dim, const unsigned int& NumNodes)
912  {
913  unsigned int Global_i;
914 
915  for(unsigned int i = 0; i < NumNodes; i++)
916  {
917  Global_i = i * (Dim + 1) + Dim;
918 
919  rRightHandSideVector[Global_i] += PBlockVector[i];
920  }
921  }
922 
923 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
924 
929 
930  static inline void Calculate2DExtrapolationMatrix(BoundedMatrix<double,3,3>& rExtrapolationMatrix)
931  {
932  //Triangle_2d_3
933  //GI_GAUSS_2
934 
935  rExtrapolationMatrix(0,0) = 1.6666666666666666666; rExtrapolationMatrix(0,1) = -0.33333333333333333333; rExtrapolationMatrix(0,2) = -0.33333333333333333333;
936  rExtrapolationMatrix(1,0) = -0.33333333333333333333; rExtrapolationMatrix(1,1) = 1.6666666666666666666; rExtrapolationMatrix(1,2) = -0.33333333333333333333;
937  rExtrapolationMatrix(2,0) = -0.33333333333333333333; rExtrapolationMatrix(2,1) = -0.33333333333333333333; rExtrapolationMatrix(2,2) = 1.6666666666666666666;
938  }
939 
940  //----------------------------------------------------------------------------------------
941 
942  static inline void Calculate2DExtrapolationMatrix(BoundedMatrix<double,4,4>& rExtrapolationMatrix)
943  {
944  //Quadrilateral_2d_4
945  //GI_GAUSS_2
946 
947  rExtrapolationMatrix(0,0) = 1.8660254037844386; rExtrapolationMatrix(0,1) = -0.5; rExtrapolationMatrix(0,2) = 0.13397459621556132; rExtrapolationMatrix(0,3) = -0.5;
948  rExtrapolationMatrix(1,0) = -0.5; rExtrapolationMatrix(1,1) = 1.8660254037844386; rExtrapolationMatrix(1,2) = -0.5; rExtrapolationMatrix(1,3) = 0.13397459621556132;
949  rExtrapolationMatrix(2,0) = 0.13397459621556132; rExtrapolationMatrix(2,1) = -0.5; rExtrapolationMatrix(2,2) = 1.8660254037844386; rExtrapolationMatrix(2,3) = -0.5;
950  rExtrapolationMatrix(3,0) = -0.5; rExtrapolationMatrix(3,1) = 0.13397459621556132; rExtrapolationMatrix(3,2) = -0.5; rExtrapolationMatrix(3,3) = 1.8660254037844386;
951  }
952 
953  //----------------------------------------------------------------------------------------
954 
955  static inline void Calculate3DExtrapolationMatrix(BoundedMatrix<double,4,4>& rExtrapolationMatrix)
956  {
957  //Tetrahedra_3d_4
958  //GI_GAUSS_2
959 
960  rExtrapolationMatrix(0,0) = -0.309016988749894905; rExtrapolationMatrix(0,1) = -0.3090169887498949046; rExtrapolationMatrix(0,2) = -0.309016988749894905; rExtrapolationMatrix(0,3) = 1.9270509662496847144;
961  rExtrapolationMatrix(1,0) = 1.9270509662496847144; rExtrapolationMatrix(1,1) = -0.30901698874989490481; rExtrapolationMatrix(1,2) = -0.3090169887498949049; rExtrapolationMatrix(1,3) = -0.30901698874989490481;
962  rExtrapolationMatrix(2,0) = -0.30901698874989490473; rExtrapolationMatrix(2,1) = 1.9270509662496847143; rExtrapolationMatrix(2,2) = -0.3090169887498949049; rExtrapolationMatrix(2,3) = -0.30901698874989490481;
963  rExtrapolationMatrix(3,0) = -0.3090169887498949048; rExtrapolationMatrix(3,1) = -0.30901698874989490471; rExtrapolationMatrix(3,2) = 1.9270509662496847143; rExtrapolationMatrix(3,3) = -0.30901698874989490481;
964  }
965 
966  //----------------------------------------------------------------------------------------
967 
968  static inline void Calculate3DExtrapolationMatrix(BoundedMatrix<double,8,8>& rExtrapolationMatrix)
969  {
970  //Hexahedra_3d_8
971  //GI_GAUSS_2
972 
973  rExtrapolationMatrix(0,0) = 2.549038105676658; rExtrapolationMatrix(0,1) = -0.6830127018922192; rExtrapolationMatrix(0,2) = 0.18301270189221927; rExtrapolationMatrix(0,3) = -0.6830127018922192;
974  rExtrapolationMatrix(0,4) = -0.6830127018922192; rExtrapolationMatrix(0,5) = 0.18301270189221927; rExtrapolationMatrix(0,6) = -0.04903810567665795; rExtrapolationMatrix(0,7) = 0.18301270189221927;
975 
976  rExtrapolationMatrix(1,0) = -0.6830127018922192; rExtrapolationMatrix(1,1) = 2.549038105676658; rExtrapolationMatrix(1,2) = -0.6830127018922192; rExtrapolationMatrix(1,3) = 0.18301270189221927;
977  rExtrapolationMatrix(1,4) = 0.18301270189221927; rExtrapolationMatrix(1,5) = -0.6830127018922192; rExtrapolationMatrix(1,6) = 0.18301270189221927; rExtrapolationMatrix(1,7) = -0.04903810567665795;
978 
979  rExtrapolationMatrix(2,0) = 0.18301270189221927; rExtrapolationMatrix(2,1) = -0.6830127018922192; rExtrapolationMatrix(2,2) = 2.549038105676658; rExtrapolationMatrix(2,3) = -0.6830127018922192;
980  rExtrapolationMatrix(2,4) = -0.04903810567665795; rExtrapolationMatrix(2,5) = 0.18301270189221927; rExtrapolationMatrix(2,6) = -0.6830127018922192; rExtrapolationMatrix(2,7) = 0.18301270189221927;
981 
982  rExtrapolationMatrix(3,0) = -0.6830127018922192; rExtrapolationMatrix(3,1) = 0.18301270189221927; rExtrapolationMatrix(3,2) = -0.6830127018922192; rExtrapolationMatrix(3,3) = 2.549038105676658;
983  rExtrapolationMatrix(3,4) = 0.18301270189221927; rExtrapolationMatrix(3,5) = -0.04903810567665795; rExtrapolationMatrix(3,6) = 0.18301270189221927; rExtrapolationMatrix(3,7) = -0.6830127018922192;
984 
985  rExtrapolationMatrix(4,0) = -0.6830127018922192; rExtrapolationMatrix(4,1) = 0.18301270189221927; rExtrapolationMatrix(4,2) = -0.04903810567665795; rExtrapolationMatrix(4,3) = 0.18301270189221927;
986  rExtrapolationMatrix(4,4) = 2.549038105676658; rExtrapolationMatrix(4,5) = -0.6830127018922192; rExtrapolationMatrix(4,6) = 0.18301270189221927; rExtrapolationMatrix(4,7) = -0.6830127018922192;
987 
988  rExtrapolationMatrix(5,0) = 0.18301270189221927; rExtrapolationMatrix(5,1) = -0.6830127018922192; rExtrapolationMatrix(5,2) = 0.18301270189221927; rExtrapolationMatrix(5,3) = -0.04903810567665795;
989  rExtrapolationMatrix(5,4) = -0.6830127018922192; rExtrapolationMatrix(5,5) = 2.549038105676658; rExtrapolationMatrix(5,6) = -0.6830127018922192; rExtrapolationMatrix(5,7) = 0.18301270189221927;
990 
991  rExtrapolationMatrix(6,0) = -0.04903810567665795; rExtrapolationMatrix(6,1) = 0.18301270189221927; rExtrapolationMatrix(6,2) = -0.6830127018922192; rExtrapolationMatrix(6,3) = 0.18301270189221927;
992  rExtrapolationMatrix(6,4) = 0.18301270189221927; rExtrapolationMatrix(6,5) = -0.6830127018922192; rExtrapolationMatrix(6,6) = 2.549038105676658; rExtrapolationMatrix(6,7) = -0.6830127018922192;
993 
994  rExtrapolationMatrix(7,0) = 0.18301270189221927; rExtrapolationMatrix(7,1) = -0.04903810567665795; rExtrapolationMatrix(7,2) = 0.18301270189221927; rExtrapolationMatrix(7,3) = -0.6830127018922192;
995  rExtrapolationMatrix(7,4) = -0.6830127018922192; rExtrapolationMatrix(7,5) = 0.18301270189221927; rExtrapolationMatrix(7,6) = -0.6830127018922192; rExtrapolationMatrix(7,7) = 2.549038105676658;
996  }
997 
998 }; /* Class PoroElementUtilities*/
999 } /* namespace Kratos.*/
1000 
1001 #endif /* KRATOS_PORO_ELEMENT_UTILITIES defined */
Kratos::Geometry
Geometry base class.
Definition: geometry.h:71
Kratos::Internals::Matrix
Definition: amatrix_interface.h:41
Kratos::Internals::Matrix::resize
void resize(std::size_t NewSize1, std::size_t NewSize2, bool preserve=0)
Definition: amatrix_interface.h:224
Kratos::PoroElementUtilities
Definition: poro_element_utilities.hpp:31
Kratos::PoroElementUtilities::AssemblePUBlockMatrix
static void AssemblePUBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 6, 18 > &PUBlockMatrix)
Definition: poro_element_utilities.hpp:756
Kratos::PoroElementUtilities::Calculate3DExtrapolationMatrix
static void Calculate3DExtrapolationMatrix(BoundedMatrix< double, 4, 4 > &rExtrapolationMatrix)
Definition: poro_element_utilities.hpp:955
Kratos::PoroElementUtilities::InterpolateVariableWithComponents
static void InterpolateVariableWithComponents(array_1d< double, 3 > &rVector, const Matrix &Ncontainer, const array_1d< double, 18 > &VariableWithComponents, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:200
Kratos::PoroElementUtilities::FillArray1dOutput
static void FillArray1dOutput(array_1d< double, 3 > &rOutputValue, const array_1d< double, 2 > &ComputedValue)
Definition: poro_element_utilities.hpp:234
Kratos::PoroElementUtilities::GetNodalVariableVector
static void GetNodalVariableVector(array_1d< double, 6 > &rNodalVariableVector, const Element::GeometryType &Geom, const Variable< array_1d< double, 3 >> &Variable, IndexType SolutionStepIndex=0)
Definition: poro_element_utilities.hpp:252
Kratos::PoroElementUtilities::GetNodalVariableVector
static void GetNodalVariableVector(array_1d< double, 12 > &rNodalVariableVector, const Element::GeometryType &Geom, const Variable< array_1d< double, 3 >> &Variable, IndexType SolutionStepIndex=0)
Definition: poro_element_utilities.hpp:284
Kratos::PoroElementUtilities::CalculateNuElementMatrix
static void CalculateNuElementMatrix(BoundedMatrix< double, 3, 9 > &rNut, const Matrix &Ncontainer, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:95
Kratos::PoroElementUtilities::AssembleUBlockMatrix
static void AssembleUBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 24, 24 > &UBlockMatrix)
Definition: poro_element_utilities.hpp:545
Kratos::PoroElementUtilities::CalculateNuMatrix
static void CalculateNuMatrix(BoundedMatrix< double, 3, 24 > &rNu, const Matrix &Ncontainer, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:81
Kratos::PoroElementUtilities::InterpolateVariableWithComponents
static void InterpolateVariableWithComponents(array_1d< double, 2 > &rVector, const Matrix &Ncontainer, const array_1d< double, 8 > &VariableWithComponents, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:167
Kratos::PoroElementUtilities::AssemblePBlockVector
static void AssemblePBlockVector(Vector &rRightHandSideVector, const TVectorType &PBlockVector, const unsigned int &Dim, const unsigned int &NumNodes)
Definition: poro_element_utilities.hpp:911
Kratos::PoroElementUtilities::CalculatePermeabilityMatrix
static void CalculatePermeabilityMatrix(BoundedMatrix< double, 3, 3 > &rPermeabilityMatrix, const Element::PropertiesType &Prop)
Definition: poro_element_utilities.hpp:348
Kratos::PoroElementUtilities::InvertMatrix2
static void InvertMatrix2(BoundedMatrix< double, 2, 2 > &rInvertedMatrix, const BoundedMatrix< double, 2, 2 > &InputMatrix)
Definition: poro_element_utilities.hpp:392
Kratos::PoroElementUtilities::AssembleUPBlockMatrix
static void AssembleUPBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 12, 4 > &UPBlockMatrix)
Definition: poro_element_utilities.hpp:620
Kratos::PoroElementUtilities::CalculateNuElementMatrix
static void CalculateNuElementMatrix(BoundedMatrix< double, 4, 16 > &rNut, const Matrix &Ncontainer, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:113
Kratos::PoroElementUtilities::AssembleUBlockMatrix
static void AssembleUBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 12, 12 > &UBlockMatrix)
Definition: poro_element_utilities.hpp:483
Kratos::PoroElementUtilities::AssembleUBlockMatrix
static void AssembleUBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 8, 8 > &UBlockMatrix)
Definition: poro_element_utilities.hpp:458
Kratos::PoroElementUtilities::CalculateNuMatrix
static void CalculateNuMatrix(BoundedMatrix< double, 3, 18 > &rNu, const Matrix &Ncontainer, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:67
Kratos::PoroElementUtilities::AssemblePUBlockMatrix
static void AssemblePUBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 4, 12 > &PUBlockMatrix)
Definition: poro_element_utilities.hpp:733
Kratos::PoroElementUtilities::CalculateNuElementMatrix
static void CalculateNuElementMatrix(BoundedMatrix< double, 3, 12 > &rNut, const Matrix &Ncontainer, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:104
Kratos::PoroElementUtilities::AssemblePUBlockMatrix
static void AssemblePUBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 8, 24 > &PUBlockMatrix)
Definition: poro_element_utilities.hpp:779
Kratos::PoroElementUtilities::AssembleUBlockMatrix
static void AssembleUBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 6, 6 > &UBlockMatrix)
Definition: poro_element_utilities.hpp:433
Kratos::PoroElementUtilities::CalculateNuElementMatrix
static void CalculateNuElementMatrix(BoundedMatrix< double, 4, 32 > &rNut, const Matrix &Ncontainer, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:137
Kratos::PoroElementUtilities::AssembleUPBlockMatrix
static void AssembleUPBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 24, 8 > &UPBlockMatrix)
Definition: poro_element_utilities.hpp:666
Kratos::PoroElementUtilities::AssembleUBlockMatrix
static void AssembleUBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 18, 18 > &UBlockMatrix)
Definition: poro_element_utilities.hpp:514
Kratos::PoroElementUtilities::GetNodalVariableVector
static void GetNodalVariableVector(array_1d< double, 24 > &rNodalVariableVector, const Element::GeometryType &Geom, const Variable< array_1d< double, 3 >> &Variable, IndexType SolutionStepIndex=0)
Definition: poro_element_utilities.hpp:318
Kratos::PoroElementUtilities::AssembleUBlockVector
static void AssembleUBlockVector(Vector &rRightHandSideVector, const array_1d< double, 8 > &UBlockVector)
Definition: poro_element_utilities.hpp:839
Kratos::PoroElementUtilities::InterpolateVariableWithComponents
static void InterpolateVariableWithComponents(array_1d< double, 3 > &rVector, const Matrix &Ncontainer, const array_1d< double, 24 > &VariableWithComponents, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:217
Kratos::PoroElementUtilities::AssembleUPBlockMatrix
static void AssembleUPBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 6, 3 > &UPBlockMatrix)
Definition: poro_element_utilities.hpp:576
Kratos::PoroElementUtilities::FillArray1dOutput
static void FillArray1dOutput(array_1d< double, 3 > &rOutputValue, const array_1d< double, 3 > &ComputedValue)
Definition: poro_element_utilities.hpp:243
Kratos::PoroElementUtilities::AssemblePBlockMatrix
static void AssemblePBlockMatrix(Matrix &rLeftHandSideMatrix, const TMatrixType &PBlockMatrix, const unsigned int &Dim, const unsigned int &NumNodes)
Definition: poro_element_utilities.hpp:803
Kratos::PoroElementUtilities::CalculatePermeabilityMatrix
static void CalculatePermeabilityMatrix(BoundedMatrix< double, 2, 2 > &rPermeabilityMatrix, const Element::PropertiesType &Prop)
Definition: poro_element_utilities.hpp:335
Kratos::PoroElementUtilities::CalculatePermeabilityMatrix
static void CalculatePermeabilityMatrix(Matrix &rPermeabilityMatrix, const Element::PropertiesType &Prop, const unsigned int &Dim)
Definition: poro_element_utilities.hpp:368
Kratos::PoroElementUtilities::AssembleUBlockVector
static void AssembleUBlockVector(Vector &rRightHandSideVector, const array_1d< double, 6 > &UBlockVector)
Definition: poro_element_utilities.hpp:822
Kratos::PoroElementUtilities::GetNodalVariableVector
static void GetNodalVariableVector(array_1d< double, 18 > &rNodalVariableVector, const Element::GeometryType &Geom, const Variable< array_1d< double, 3 >> &Variable, IndexType SolutionStepIndex=0)
Definition: poro_element_utilities.hpp:301
Kratos::PoroElementUtilities::CalculateNuMatrix
static void CalculateNuMatrix(BoundedMatrix< double, 2, 8 > &rNu, const Matrix &Ncontainer, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:48
Kratos::PoroElementUtilities::CalculateNuElementMatrix
static void CalculateNuElementMatrix(BoundedMatrix< double, 4, 24 > &rNut, const Matrix &Ncontainer, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:123
Kratos::PoroElementUtilities::Calculate2DExtrapolationMatrix
static void Calculate2DExtrapolationMatrix(BoundedMatrix< double, 4, 4 > &rExtrapolationMatrix)
Definition: poro_element_utilities.hpp:942
Kratos::PoroElementUtilities::AssemblePUBlockMatrix
static void AssemblePUBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 4, 8 > &PUBlockMatrix)
Definition: poro_element_utilities.hpp:711
Kratos::PoroElementUtilities::InterpolateVariableWithComponents
static void InterpolateVariableWithComponents(array_1d< double, 2 > &rVector, const Matrix &Ncontainer, const array_1d< double, 6 > &VariableWithComponents, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:151
Kratos::PoroElementUtilities::GetNodalVariableVector
static void GetNodalVariableVector(array_1d< double, 8 > &rNodalVariableVector, const Element::GeometryType &Geom, const Variable< array_1d< double, 3 >> &Variable, IndexType SolutionStepIndex=0)
Definition: poro_element_utilities.hpp:268
Kratos::PoroElementUtilities::AssembleUBlockVector
static void AssembleUBlockVector(Vector &rRightHandSideVector, const array_1d< double, 18 > &UBlockVector)
Definition: poro_element_utilities.hpp:874
Kratos::PoroElementUtilities::AssembleUBlockVector
static void AssembleUBlockVector(Vector &rRightHandSideVector, const array_1d< double, 12 > &UBlockVector)
Definition: poro_element_utilities.hpp:856
Kratos::PoroElementUtilities::CalculateNuMatrix
static void CalculateNuMatrix(BoundedMatrix< double, 3, 12 > &rNu, const Matrix &Ncontainer, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:57
Kratos::PoroElementUtilities::AssemblePUBlockMatrix
static void AssemblePUBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 3, 6 > &PUBlockMatrix)
Definition: poro_element_utilities.hpp:689
Kratos::PoroElementUtilities::AssembleUPBlockMatrix
static void AssembleUPBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 8, 4 > &UPBlockMatrix)
Definition: poro_element_utilities.hpp:598
Kratos::PoroElementUtilities::AssembleUPBlockMatrix
static void AssembleUPBlockMatrix(Matrix &rLeftHandSideMatrix, const BoundedMatrix< double, 18, 6 > &UPBlockMatrix)
Definition: poro_element_utilities.hpp:643
Kratos::PoroElementUtilities::InterpolateVariableWithComponents
static void InterpolateVariableWithComponents(array_1d< double, 3 > &rVector, const Matrix &Ncontainer, const array_1d< double, 12 > &VariableWithComponents, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:183
Kratos::PoroElementUtilities::Calculate3DExtrapolationMatrix
static void Calculate3DExtrapolationMatrix(BoundedMatrix< double, 8, 8 > &rExtrapolationMatrix)
Definition: poro_element_utilities.hpp:968
Kratos::PoroElementUtilities::Calculate2DExtrapolationMatrix
static void Calculate2DExtrapolationMatrix(BoundedMatrix< double, 3, 3 > &rExtrapolationMatrix)
Definition: poro_element_utilities.hpp:930
Kratos::PoroElementUtilities::AssembleUBlockVector
static void AssembleUBlockVector(Vector &rRightHandSideVector, const array_1d< double, 24 > &UBlockVector)
Definition: poro_element_utilities.hpp:892
Kratos::PoroElementUtilities::CalculateNuMatrix
static void CalculateNuMatrix(BoundedMatrix< double, 2, 6 > &rNu, const Matrix &Ncontainer, const unsigned int &GPoint)
Definition: poro_element_utilities.hpp:39
Kratos::Properties
Properties encapsulates data shared by different Elements or Conditions. It can store any type of dat...
Definition: properties.h:69
Kratos::Variable
Variable class contains all information needed to store and retrive data from a data container.
Definition: variable.h:63
Kratos::array_1d< double, 2 >
math_utils.h
Kratos
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
Kratos::ZeroVector
KratosZeroVector< double > ZeroVector
Definition: amatrix_interface.h:561
Kratos::noalias
T & noalias(T &TheMatrix)
Definition: amatrix_interface.h:484
quadrature.j
int j
Definition: quadrature.py:648
tensors::i
integer i
Definition: TensorModule.f:17
poromechanics_application_variables.h
Generated on Sat Feb 3 2024 02:11:12 for KratosMultiphysics by doxygen 1.9.1